Filtering VCF files based on VAF giving incorrect results
2
0
Entering edit mode
7 months ago
Arton ▴ 20

I'm trying to filter variants that have VAF<0.2. I used picard FilterVcf, bcftools filter, and tried to use GATK VariantFiltration, and neither of these worked. picard and bcftools only filter a portion of variants but I still get a long list of variants with VAF<0.2 labeled as "PASS". I haven't figured out how to use GATK VariantFiltration for VAF. Any help to figure out the problem would be greatly appreciated. Here are my codes:

java -jar picard.jar FilterVcf  \
    I=input.vcf
    MIN_AB=0.2 \
    O=output.vcf   

bcftools filter --include '(FMT/AD[0:1])/(FMT/AD[0:0]+FMT/AD[0:1]) >=0.25 && (FMT/DP[0:0] >=10)' \
        input.vcf \
        --output output.vcf \
        --soft-filter Filter1 \
        --mode +
GATK picard VCF bcftools • 772 views
ADD COMMENT
1
Entering edit mode
7 months ago

using vcfilterjdk: https://jvarkit.readthedocs.io/en/latest/VcfFilterJdk/

 java -jar ${JVARKIT_DIST}/jvarkit.jar vcffilterjdk -e 'return !variant.getGenotypes().stream().filter(G->G.isHet() && G.hasAD() && G.getAD().length==2).anyMatch(G->{double n1=G.getAD()[0];double n2=G.getAD()[1];if(n1+n2<=0) return false; double f=(n2/(n1+n2));return (f<0.2 || f>0.8);});'  --filter Filter1 input.vcf
ADD COMMENT
0
Entering edit mode

Thank you. I can't try a completely new tool but I will keep your command for future use. Do you think my commands were correct? For example what does mean that MIN_AB is double?

ADD REPLY
0
Entering edit mode

Do you think my commands were correct?

I'm not familiar enough with this bcftools syntax.

but what if there is no AD, what if both are '0', what if there is a NO-CALL at this position ? what if the variant is multiallelic ?

ADD REPLY

Login before adding your answer.

Traffic: 3701 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