How is QD calculate for multiple heterozygous samples?
Morning,
I have a question related to the calculation of QD in mutiple sample VCFs.
In GATK documentation QD is defined as QUAL/AD. Then, I wonder how QD is calculated when there are multiple heterozygous samples at a given locus. I mean if a single value is given per locus for QD in my VCF files that contain 6 samples is QD averaged among the samples or is QD denominator the average of AD samples?
I am using GATK v4.3.0.0.
Many thanks for your answer!
Best,
Vinciane
-
When QD is calculated only depth from samples with non homozygous reference calls is taken into accounts.
Below is the code snippet that calculates the DP for each site. We don't take averages of DP values when calculating QD for multiple samples. The sum of DPs calculated by the below function is used to divide the QUAL value to find QD.
public static int getDepth(final GenotypesContext genotypes, final AlleleLikelihoods<GATKRead, Allele> likelihoods) {
int depth = 0;
int ADrestrictedDepth = 0;
for ( final Genotype genotype : genotypes ) {
// we care only about variant calls with likelihoods
if ( !genotype.isHet() && !genotype.isHomVar() ) {
continue;
}
// if we have the AD values for this sample, let's make sure that the variant depth is greater than 1!
if ( genotype.hasAD() ) {
final int[] AD = genotype.getAD();
final int totalADdepth = (int) MathUtils.sum(AD);
if ( totalADdepth != 0 ) {
if (totalADdepth - AD[0] > 1) {
ADrestrictedDepth += totalADdepth;
}
depth += totalADdepth;
continue;
}
}
// if there is no AD value or it is a dummy value, we want to look to other means to get the depth
if (likelihoods != null) {
depth += likelihoods.sampleEvidenceCount(likelihoods.indexOfSample(genotype.getSampleName()));
} else if ( genotype.hasDP() ) {
depth += genotype.getDP();
}
}
// if the AD-restricted depth is a usable value (i.e. not zero), then we should use that one going forward
if ( ADrestrictedDepth > 0 ) {
depth = ADrestrictedDepth;
}
return depth;
}I hope this helps.
Please sign in to leave a comment.
1 comment