How could I set any threshold to control FILTER column which PASS or RefCall?

[Rochsen]

Hi,

I have two problem about this:

First question, this picture is VCF line, this position's FILTER is RefCall called by deepvariant:1.6.0 (docker's image)

And GATK4-HaplotypeCaller（The Genome Analysis Toolkit (GATK) v4.2.6.1,HTSJDK Version: 2.24.1,Picard Version: 2.27.1) Call the same position Like follow:

chr5 147499874 . G GA 52.60 . AC=1;AF=0.500;AN=2;BaseQRankSum=-2.210;DP=61;ExcessHet=0.0000;FS=0.000;MLEAC=1;MLEAF=0.500;MQ=60.00;MQRankSum=0.000;QD=0.86;ReadPosRankSum=2.170;SOR=0.846 GT:AD:DP:GQ:PL 0/1:52,9:61:60:60,0,1625

Why not support INFO column in deepvariant for user to filter snp or indel by them, how could I make FILTER=PASS ?

Last question, Why deepvariant call this position's QUAL is zero, and GATK call this position's QUAL is 52.60? which is more accurate?

[Rochsen]

Note：this issue's code is the same as #768 in closed issue

[AndrewCarroll]

Hi @Luosanmu,

For your filtering question. We don't directly calculate many of the statistics in the INFO field for variant calling. We observe that the most effective way to adjust filters is to use the Genotype Quality (GQ) property or the Phred Likelihood (PL) fields (GQ is mathematically derived from PL). Do you want to increase sensitivity? The best way is likely to post-process the VCF to extract the low-confidence REF calls using this field.

For the question - why does DeepVariant make a call that differs from GATK. For any single call, it's difficult to say the exact reasons. Sometimes, looking at the reads and the reference in the region can give clues about why a call would be made. If you have an IGV screenshot showing the region it might be informative. DeepVariant does seem very confident that there isn't a variant here.

It would also be helpful to know something about the sequencing and prep. Is this Illumina data? PacBio data? Is this a PCR-free prep, or does it include PCR? things like that.

Thank you, Andrew

[Rochsen]

Hi, @AndrewCarroll

chr5-147499874-G-GA with IGV pictures are here

And this sample is Illumina NGS data.

Thanks, Luosanmu

[AndrewCarroll]

Thank you @Luosanmu

I see. From your image, I think I understand why DeepVariant would make a REF call here. The variant in question is a 1-bp extension of a homopolymer (10A -> 11A). Homopolymers are generally difficult to sequence through. The number of reference-supporting reads are 47 and alternate-supporting reads are 10 (~16%), which is far from the typically-expected 50% if the position is heterozygous.

DeepVariant's model has to weigh which probability is more likely: that this is a real HET event and the random sampling of the alleles causes the observations to be skewed as far as 16%, or is there a sufficiently recurring 1bp insertion error during sequencing that explains these insertions at this ratio.

Presumably, over the bulk of DeepVariant's training, when it has seen similar situations, in more cases these are insertion errors. Now, whether that is what is truly going on in your sample, it's difficult for me as a human to say.

[pichuan]

Hi @Luosanmu , Due to inactivity on this issue, I'll close it. Please feel free to follow up if you have more questions.