diff --git a/protected/gatk-tools-protected/src/main/java/org/broadinstitute/gatk/tools/walkers/variantutils/CalculateGenotypePosteriors.java b/protected/gatk-tools-protected/src/main/java/org/broadinstitute/gatk/tools/walkers/variantutils/CalculateGenotypePosteriors.java index 224caa97f..41df5ddd1 100644 --- a/protected/gatk-tools-protected/src/main/java/org/broadinstitute/gatk/tools/walkers/variantutils/CalculateGenotypePosteriors.java +++ b/protected/gatk-tools-protected/src/main/java/org/broadinstitute/gatk/tools/walkers/variantutils/CalculateGenotypePosteriors.java @@ -185,6 +185,8 @@ import java.util.*; * -o output.withPosteriors.vcf * * + *

Caveat

+ *

If applying family priors, only diploid family genotypes are supported

*/ @DocumentedGATKFeature( groupName = HelpConstants.DOCS_CAT_VARDISC, extraDocs = {CommandLineGATK.class} ) public class CalculateGenotypePosteriors extends RodWalker { diff --git a/protected/gatk-tools-protected/src/main/java/org/broadinstitute/gatk/tools/walkers/variantutils/FamilyLikelihoodsUtils.java b/protected/gatk-tools-protected/src/main/java/org/broadinstitute/gatk/tools/walkers/variantutils/FamilyLikelihoodsUtils.java index 9049e5c79..2807e72c8 100644 --- a/protected/gatk-tools-protected/src/main/java/org/broadinstitute/gatk/tools/walkers/variantutils/FamilyLikelihoodsUtils.java +++ b/protected/gatk-tools-protected/src/main/java/org/broadinstitute/gatk/tools/walkers/variantutils/FamilyLikelihoodsUtils.java @@ -51,6 +51,7 @@ package org.broadinstitute.gatk.tools.walkers.variantutils; +import htsjdk.variant.vcf.VCFConstants; import org.apache.log4j.Logger; import org.broadinstitute.gatk.engine.samples.Sample; import org.broadinstitute.gatk.utils.MathUtils; @@ -263,7 +264,7 @@ public class FamilyLikelihoodsUtils { updateFamilyGenotypes(vc, mother, father, child, trioGenotypes); //replace uses sample names to match genotypes, so order doesn't matter - if (trioGenotypes.size() > 0) { + if (!trioGenotypes.isEmpty()) { genotypesContext.replace(trioGenotypes.get(0)); genotypesContext.replace(trioGenotypes.get(1)); genotypesContext.replace(trioGenotypes.get(2)); @@ -440,7 +441,7 @@ public class FamilyLikelihoodsUtils { //In case of null, unavailable or no call, all likelihoods are log10(1/3) else if(genotype == null || !hasCalledGT(genotype.getType()) || genotype.getLikelihoods() == null){ - likelihoods = new double[3]; + likelihoods = new double[NUM_CALLED_GENOTYPETYPES]; likelihoods[0] = LOG10_OF_ONE_THIRD; likelihoods[1] = LOG10_OF_ONE_THIRD; likelihoods[2] = LOG10_OF_ONE_THIRD; @@ -450,6 +451,13 @@ public class FamilyLikelihoodsUtils { else likelihoods = GeneralUtils.normalizeFromLog10(genotype.getLikelihoods().getAsVector(),true,true); + if (likelihoods.length != NUM_CALLED_GENOTYPETYPES) { + final String key = genotype.hasExtendedAttribute(GATKVCFConstants.PHRED_SCALED_POSTERIORS_KEY) ? + GATKVCFConstants.PHRED_SCALED_POSTERIORS_KEY : VCFConstants.GENOTYPE_PL_KEY; + throw new UserException(genotype + " has " + likelihoods.length + " " + key + " values, should be " + NUM_CALLED_GENOTYPETYPES + + " since only the diploid case is supported when applying family priors."); + } + likelihoodsMap.put(GenotypeType.HOM_REF,likelihoods[genotypeTypeToValue(GenotypeType.HOM_REF)]); likelihoodsMap.put(GenotypeType.HET,likelihoods[genotypeTypeToValue(GenotypeType.HET)]); likelihoodsMap.put(GenotypeType.HOM_VAR, likelihoods[genotypeTypeToValue(GenotypeType.HOM_VAR)]); diff --git a/protected/gatk-tools-protected/src/test/java/org/broadinstitute/gatk/tools/walkers/variantutils/CalculateGenotypePosteriorsIntegrationTest.java b/protected/gatk-tools-protected/src/test/java/org/broadinstitute/gatk/tools/walkers/variantutils/CalculateGenotypePosteriorsIntegrationTest.java index 0432795cf..92bdae781 100644 --- a/protected/gatk-tools-protected/src/test/java/org/broadinstitute/gatk/tools/walkers/variantutils/CalculateGenotypePosteriorsIntegrationTest.java +++ b/protected/gatk-tools-protected/src/test/java/org/broadinstitute/gatk/tools/walkers/variantutils/CalculateGenotypePosteriorsIntegrationTest.java @@ -52,6 +52,7 @@ package org.broadinstitute.gatk.tools.walkers.variantutils; import org.broadinstitute.gatk.engine.walkers.WalkerTest; +import org.broadinstitute.gatk.utils.exceptions.UserException; import org.testng.annotations.Test; import java.util.Arrays; @@ -61,6 +62,7 @@ public class CalculateGenotypePosteriorsIntegrationTest extends WalkerTest { private static String CEUtrioFamilyFile = privateTestDir + "CEUtrio.ped"; private static String CEUtrioTest = privateTestDir + "CEUtrioTest.vcf"; private static String CEUtrioPopPriorsTest = privateTestDir + "CEUtrioPopPriorsTest.vcf"; + private static String CEUtrioMixedPloidyTest = privateTestDir + "CEUtrioMixedPloidy.vcf"; private static String threeMemberNonTrioFamilyFile = privateTestDir + "threeMemberNonTrio.ped"; private static String getThreeMemberNonTrioTest = privateTestDir + "threeMemberNonTrioTest.vcf"; @@ -135,5 +137,16 @@ public class CalculateGenotypePosteriorsIntegrationTest extends WalkerTest { executeTest("testFamilyPriors", spec); } - + @Test(enabled = true) + public void testFamilyPriorsMixedPloidy() { + WalkerTestSpec spec = new WalkerTestSpec( + "-T CalculateGenotypePosteriors --no_cmdline_in_header " + + " -o %s" + + " -R " + b37KGReference + + " -ped " + CEUtrioFamilyFile + + " -V " + CEUtrioMixedPloidyTest, + 1, + UserException.class); + executeTest("testFamilyPriorsMixedPloidy", spec); + } }