Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Fix reported bug where no-call genotypes with no reads get genotype #5667

Merged
merged 3 commits into from
Feb 14, 2019
Merged
Show file tree
Hide file tree
Changes from 2 commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
Original file line number Diff line number Diff line change
Expand Up @@ -187,7 +187,26 @@ public static double[] parsePosteriorsIntoProbSpace(Genotype genotype) {

// return the double[] of likelihoods if available, otherwise null
private static double[] getLikelihoodsVector(Genotype genotype) {
return genotype.hasLikelihoods() ? genotype.getLikelihoods().getAsVector() : null;
return hasRealLikelihoods(genotype) ? genotype.getLikelihoods().getAsVector() : null;
}

/**
* return true if the genotype contains real data for likelihoods, not just zero placeholder
*/
private static boolean hasRealLikelihoods(Genotype genotype) {
if (!genotype.hasLikelihoods()) {
return false;
}
if (genotype.hasDP() && genotype.getDP() == 0) {
int[] pls = genotype.getPL();
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

How about return MathUtils.max(genotype.getPL()) > 0?

for (int i = 0; i < pls.length; i++) {
if (pls[i] > 0) {
return true;
}
}
return false; //if there's no depth and the PLs are all zero, then there's really no data here; don't rely on no-call
}
return true;
}

/**
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -357,7 +357,9 @@ public void testGVCFModeGenotypePosteriors() throws Exception {

for (final VariantContext vc : results.getRight()) {
final Genotype g = vc.getGenotype(0);
Assert.assertTrue(g.hasExtendedAttribute(GATKVCFConstants.PHRED_SCALED_POSTERIORS_KEY));
if (g.hasDP() && g.getDP() > 0 && g.hasGQ() && g.getGQ() > 0) {
Assert.assertTrue(g.hasExtendedAttribute(GATKVCFConstants.PHRED_SCALED_POSTERIORS_KEY));
}
if (isGVCFReferenceBlock(vc) ) {
Assert.assertTrue(!vc.hasAttribute(GATKVCFConstants.GENOTYPE_PRIOR_KEY));
}
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -57,6 +57,7 @@ public void testNumRefWithPanel() throws IOException {

@Test
//use the first 20 variants to save time; they have a nice range of AC from 4 to over 4000
//three variants have GTs with zero depth and PLs=[1,0,0], which should get PPs
public void testUsingDiscoveredAF() throws IOException {
final IntegrationTestSpec spec = new IntegrationTestSpec(
" -O %s" +
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -391,6 +391,18 @@ public void testDifferentRefAllelesAndNonRefCounts() {
Assert.assertEquals(arraysEq(expectedPPs, _mleparse((List<Integer>)test1result.getGenotype(0).getAnyAttribute(GATKVCFConstants.PHRED_SCALED_POSTERIORS_KEY))), "");
}

@Test
public void testGenotypesWithNoDataDoNotChange() {
final GenotypeBuilder gb = new GenotypeBuilder("s1");
gb.alleles(Arrays.asList(Allele.NO_CALL, Allele.NO_CALL)).DP(0).AD(new int[]{0,0}).GQ(0).PL(new int[]{0,0,0});
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

How about chaining .alleles etc with the previous line?

final VariantContext vc = makeVC("noCall", Arrays.asList(Allele.create("A",true), Allele.create("T")), gb.make());
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Allele.create("A", true) can be replaced with the pre-defined Allele.REF_A.

final List<VariantContext> resources = new ArrayList<>();
resources.add(new VariantContextBuilder(makeVC("2", Arrays.asList(Aref,T))).attribute(GATKVCFConstants.MLE_ALLELE_COUNT_KEY,500).attribute(VCFConstants.ALLELE_NUMBER_KEY,1000).make());
final VariantContext result = PosteriorProbabilitiesUtils.calculatePosteriorProbs(vc, resources, 0, defaultOptions);
Assert.assertTrue(result.getGenotype("s1").isNoCall());
Assert.assertFalse(result.getGenotype("s1").hasExtendedAttribute(GATKVCFConstants.PHRED_SCALED_POSTERIORS_KEY));
}

private double[] pl2gl(final int[] pl) {
final double[] gl = new double[pl.length];
for ( int idx = 0; idx < gl.length; idx++ ) {
Expand Down