Bcftools: filtering by AD as a percentage of DP
1
0
Entering edit mode
20 months ago
Timotheus ▴ 40

Hello,

I'm filtering my two-sample diploid VCF by, among other filters, coverage (DP) as follows (example for sample with index 1):

bcftools filter -i 'FORMAT/DP[1]>90 && FORMAT/DP[1]<140' my.vcf -o my_filtered.vcf

I nonetheless discovered that allelic depth (AD) is markedly off at some het sites, i.e., the proportions of the two alleles are far from 50:50. How do I put bounds on these proportions, e.g., so that they are at most 40:60?

bcftools • 854 views
ADD COMMENT
1
Entering edit mode
20 months ago

using vcffilterjdk http://lindenb.github.io/jvarkit/VcfFilterJdk.html

java -jar dist/vcffilterjdk.jar -e 'return variant.getGenotypes().stream().filter(G->G.isHet() && G.hasAD()).map(G->G.getAD()).filter(A->A.length==2).allMatch(A->{double r=A[0],a=A[1];if(r+a==0) return false; double f=a/(r+a); return f < 0.6 && f>0.4 ;});'
ADD COMMENT
0
Entering edit mode

Thank you, Pierre: is there no bcftools solution?

ADD REPLY
0
Entering edit mode

But your script works very well, thanks! A pity perhaps it uses a Java expression... how should I apply the filter just to one sample, for example, with index 1?

ADD REPLY

Login before adding your answer.

Traffic: 1752 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6