This is a test version of Biostars. For the public version, visit https://www.biostars.org.
Bcftools: filtering by AD as a percentage of DP

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
• 1,772 views • link
written by Timotheus • Timotheus

1 answer

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 ;});'
• 0 views • link
written by Pierre Lindenbaum Pierre Lindenbaum

Thank you, Pierre: is there no bcftools solution?

• 0 views • link
written by Timotheus • Timotheus

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?

• 0 views • link
written by Timotheus • Timotheus

Log in to answer this question.

More posts like this