From dce6b91f7d94ab9265b56e65b691c502c84f64d1 Mon Sep 17 00:00:00 2001 From: Eric Banks Date: Fri, 16 Mar 2012 11:14:37 -0400 Subject: [PATCH] Add a conversion from the deprecated PL ordering to the new one. We need this for the DiploidSNPGenotypeLikelihoods which still use the old ordering. My intention is for this to be a temporary patch, but changing the ordering in DiploidSNPGenotypeLikelihoods is not appriopriate for committing to stable as it will break all of the external tools (e.g. MuTec) that are built on top of the class. We will have to talk to e.g. Kristian to see how disruptive this will be. Added unit tests to the GL conversions and indexing. --- .../walkers/genotyper/DiploidGenotype.java | 1 - ...NPGenotypeLikelihoodsCalculationModel.java | 2 +- .../variantcontext/GenotypeLikelihoods.java | 23 +++++++++++++++++++ .../GenotypeLikelihoodsUnitTest.java | 5 +--- 4 files changed, 25 insertions(+), 6 deletions(-) diff --git a/public/java/src/org/broadinstitute/sting/gatk/walkers/genotyper/DiploidGenotype.java b/public/java/src/org/broadinstitute/sting/gatk/walkers/genotyper/DiploidGenotype.java index 09936c112..4aa580052 100755 --- a/public/java/src/org/broadinstitute/sting/gatk/walkers/genotyper/DiploidGenotype.java +++ b/public/java/src/org/broadinstitute/sting/gatk/walkers/genotyper/DiploidGenotype.java @@ -64,7 +64,6 @@ public enum DiploidGenotype { return r != base2; else return base2 == r; - //return MathUtils.countOccurrences(r, this.toString()) == 1; } public boolean isHom() { diff --git a/public/java/src/org/broadinstitute/sting/gatk/walkers/genotyper/SNPGenotypeLikelihoodsCalculationModel.java b/public/java/src/org/broadinstitute/sting/gatk/walkers/genotyper/SNPGenotypeLikelihoodsCalculationModel.java index 20b1b0b14..23d4e4905 100755 --- a/public/java/src/org/broadinstitute/sting/gatk/walkers/genotyper/SNPGenotypeLikelihoodsCalculationModel.java +++ b/public/java/src/org/broadinstitute/sting/gatk/walkers/genotyper/SNPGenotypeLikelihoodsCalculationModel.java @@ -176,7 +176,7 @@ public class SNPGenotypeLikelihoodsCalculationModel extends GenotypeLikelihoodsC final double[] likelihoods = sampleData.GL.getLikelihoods(); final int PLindexOfBestGL = MathUtils.maxElementIndex(likelihoods); if ( PLindexOfBestGL != PLindexOfRef ) { - GenotypeLikelihoods.GenotypeLikelihoodsAllelePair alleles = GenotypeLikelihoods.getAllelePair(3, PLindexOfBestGL); + GenotypeLikelihoods.GenotypeLikelihoodsAllelePair alleles = GenotypeLikelihoods.getAllelePairUsingDeprecatedOrdering(PLindexOfBestGL); if ( alleles.alleleIndex1 != baseIndexOfRef ) likelihoodSums[alleles.alleleIndex1] += likelihoods[PLindexOfBestGL] - likelihoods[PLindexOfRef]; // don't double-count it diff --git a/public/java/src/org/broadinstitute/sting/utils/variantcontext/GenotypeLikelihoods.java b/public/java/src/org/broadinstitute/sting/utils/variantcontext/GenotypeLikelihoods.java index 54b4cf230..cc0d9788b 100755 --- a/public/java/src/org/broadinstitute/sting/utils/variantcontext/GenotypeLikelihoods.java +++ b/public/java/src/org/broadinstitute/sting/utils/variantcontext/GenotypeLikelihoods.java @@ -272,6 +272,29 @@ public class GenotypeLikelihoods { return PLIndexToAlleleIndex[PLindex]; } + // An index conversion from the deprecated PL ordering to the new VCF-based ordering for up to 3 alternate alleles + protected static int[] PLindexConversion = new int[]{0, 1, 3, 6, 2, 4, 7, 5, 8, 9}; + + /** + * get the allele index pair for the given PL using the deprecated PL ordering: + * AA,AB,AC,AD,BB,BC,BD,CC,CD,DD instead of AA,AB,BB,AC,BC,CC,AD,BD,CD,DD. + * Although it's painful to keep this conversion around, our DiploidSNPGenotypeLikelihoods class uses the deprecated + * ordering and I know with certainty that external users have built code on top of it; changing it now would + * cause a whole lot of heartache for our collaborators, so for now at least there's a standard conversion method. + * This method assumes at most 3 alternate alleles. + * TODO -- address this issue at the source by updating DiploidSNPGenotypeLikelihoods. + * + * @param PLindex the PL index + * @return the allele index pair + */ + public static GenotypeLikelihoodsAllelePair getAllelePairUsingDeprecatedOrdering(final int PLindex) { + // make sure that we've cached enough data + if ( PLindex >= PLIndexToAlleleIndex.length ) + calculatePLcache(PLindex); + + return PLIndexToAlleleIndex[PLindexConversion[PLindex]]; + } + /** * get the PL indexes (AA, AB, BB) for the given allele pair; assumes allele1Index <= allele2Index. * diff --git a/public/java/test/org/broadinstitute/sting/utils/variantcontext/GenotypeLikelihoodsUnitTest.java b/public/java/test/org/broadinstitute/sting/utils/variantcontext/GenotypeLikelihoodsUnitTest.java index 6320a140d..638fd2531 100755 --- a/public/java/test/org/broadinstitute/sting/utils/variantcontext/GenotypeLikelihoodsUnitTest.java +++ b/public/java/test/org/broadinstitute/sting/utils/variantcontext/GenotypeLikelihoodsUnitTest.java @@ -137,13 +137,10 @@ public class GenotypeLikelihoodsUnitTest { @Test public void testCalculatePLindex(){ - // AA,AB,BB,AC,BC,CC,AD,BD,CD,DD called in the order of AA,AB,AC,AD,BB,BC,BD,CC,CD,DD - int[] indexes = new int[]{0, 1, 3, 6, 2, 4, 7, 5, 8, 9}; - int counter = 0; for ( int i = 0; i <= 3; i++ ) { for ( int j = i; j <= 3; j++ ) { - Assert.assertEquals(GenotypeLikelihoods.calculatePLindex(i, j), indexes[counter++], "PL index of alleles " + i + "," + j + " was not calculated correctly"); + Assert.assertEquals(GenotypeLikelihoods.calculatePLindex(i, j), GenotypeLikelihoods.PLindexConversion[counter++], "PL index of alleles " + i + "," + j + " was not calculated correctly"); } } }