From bf935a6ab11d600f87b8e5285831ff2c8e448c8a Mon Sep 17 00:00:00 2001 From: ebanks Date: Mon, 23 Nov 2009 01:31:48 +0000 Subject: [PATCH] 1. Fixed bug in PrimaryBaseSecondaryBaseSymmetry code (not checking for null before trying to access object's methods) which was causing Integration Tests to fail. 2. Retired allele frequency range from UG, which wasn't very useful. git-svn-id: file:///humgen/gsa-scr1/gsa-engineering/svn_contents/trunk@2113 348d0f76-0448-11de-a6fe-93d51630548a --- .../PrimaryBaseSecondaryBaseSymmetry.java | 6 +++- .../genotyper/GenotypeCalculationModel.java | 4 --- ...JointEstimateGenotypeCalculationModel.java | 34 ------------------ .../genotyper/UnifiedArgumentCollection.java | 3 -- .../utils/genotype/AlleleFrequencyBacked.java | 36 ------------------- .../genotype/vcf/VCFGenotypeLocusData.java | 20 +---------- .../vcf/VCFGenotypeWriterAdapter.java | 3 -- .../UnifiedGenotyperIntegrationTest.java | 4 +-- 8 files changed, 8 insertions(+), 102 deletions(-) diff --git a/java/src/org/broadinstitute/sting/gatk/walkers/annotator/PrimaryBaseSecondaryBaseSymmetry.java b/java/src/org/broadinstitute/sting/gatk/walkers/annotator/PrimaryBaseSecondaryBaseSymmetry.java index d6381332c..2f653eec6 100755 --- a/java/src/org/broadinstitute/sting/gatk/walkers/annotator/PrimaryBaseSecondaryBaseSymmetry.java +++ b/java/src/org/broadinstitute/sting/gatk/walkers/annotator/PrimaryBaseSecondaryBaseSymmetry.java @@ -71,7 +71,11 @@ public class PrimaryBaseSecondaryBaseSymmetry implements VariantAnnotation{ return null; } - char[] secondaryPileup = p.getSecondaryBasePileup().toCharArray(); + String secondaryPileupStr = p.getSecondaryBasePileup(); + if ( secondaryPileupStr == null ) + return null; + + char[] secondaryPileup = secondaryPileupStr.toCharArray(); int depth = p.size(); int support = 0; diff --git a/java/src/org/broadinstitute/sting/gatk/walkers/genotyper/GenotypeCalculationModel.java b/java/src/org/broadinstitute/sting/gatk/walkers/genotyper/GenotypeCalculationModel.java index f84e68351..36c1f3298 100755 --- a/java/src/org/broadinstitute/sting/gatk/walkers/genotyper/GenotypeCalculationModel.java +++ b/java/src/org/broadinstitute/sting/gatk/walkers/genotyper/GenotypeCalculationModel.java @@ -38,7 +38,6 @@ public abstract class GenotypeCalculationModel implements Cloneable { protected int POOL_SIZE; protected double CONFIDENCE_THRESHOLD; protected double MINIMUM_ALLELE_FREQUENCY; - protected double ALLELE_FREQUENCY_RANGE; protected boolean REPORT_SLOD; protected int maxDeletionsInPileup; protected String assumedSingleSample; @@ -71,9 +70,6 @@ public abstract class GenotypeCalculationModel implements Cloneable { POOL_SIZE = UAC.POOLSIZE; CONFIDENCE_THRESHOLD = UAC.CONFIDENCE_THRESHOLD; MINIMUM_ALLELE_FREQUENCY = UAC.MINIMUM_ALLELE_FREQUENCY; - ALLELE_FREQUENCY_RANGE = UAC.ALLELE_FREQUENCY_RANGE; - if ( ALLELE_FREQUENCY_RANGE < 0.0 || ALLELE_FREQUENCY_RANGE > 1.0 ) - throw new StingException("Allele frequency fraction must be a value between 0.0 and 1.0"); maxDeletionsInPileup = UAC.MAX_DELETIONS; assumedSingleSample = UAC.ASSUME_SINGLE_SAMPLE; if ( UAC.VERBOSE != null ) { diff --git a/java/src/org/broadinstitute/sting/gatk/walkers/genotyper/JointEstimateGenotypeCalculationModel.java b/java/src/org/broadinstitute/sting/gatk/walkers/genotyper/JointEstimateGenotypeCalculationModel.java index c578eb0ec..b262e5ede 100644 --- a/java/src/org/broadinstitute/sting/gatk/walkers/genotyper/JointEstimateGenotypeCalculationModel.java +++ b/java/src/org/broadinstitute/sting/gatk/walkers/genotyper/JointEstimateGenotypeCalculationModel.java @@ -261,11 +261,6 @@ public abstract class JointEstimateGenotypeCalculationModel extends GenotypeCalc } if ( locusdata instanceof AlleleFrequencyBacked ) { ((AlleleFrequencyBacked)locusdata).setAlleleFrequency((double)bestAFguess / (double)(frequencyEstimationPoints-1)); - // frequenc range doesn't make sense for single samples - if ( getNSamples(contexts) > 1 ) { - AlleleFrequencyBacked.AlleleFrequencyRange range = computeAFrange(alleleFrequencyPosteriors[indexOfMax], frequencyEstimationPoints-1, bestAFguess, ALLELE_FREQUENCY_RANGE); - ((AlleleFrequencyBacked)locusdata).setAlleleFrequencyRange(range); - } } if ( locusdata instanceof IDBacked ) { rodDbSNP dbsnp = getDbSNP(tracker); @@ -308,33 +303,4 @@ public abstract class JointEstimateGenotypeCalculationModel extends GenotypeCalc return new Pair, GenotypeLocusData>(calls, locusdata); } - - // computes the range of allele frequencies making up the given fraction of the total probability - private static AlleleFrequencyBacked.AlleleFrequencyRange computeAFrange(double[] alleleFrequencyProbs, int N, int bestAFguess, double fraction) { - double totalProb = alleleFrequencyProbs[bestAFguess]; - int lowIndex = bestAFguess; - int highIndex = bestAFguess; - - // it's possible that AF=0 contains more probability than 1-fraction - if ( alleleFrequencyProbs[0] >= (1.0 - fraction) ) { - // in this case, the range is all possible AFs - lowIndex = 1; - highIndex = N; - } - // otherwise, find the range moving out from the best AF guess - else { - while ( totalProb < fraction ) { - if ( lowIndex > 1 ) { - lowIndex--; - totalProb += alleleFrequencyProbs[lowIndex]; - } - if ( highIndex < N ) { - highIndex++; - totalProb += alleleFrequencyProbs[highIndex]; - } - } - } - - return new AlleleFrequencyBacked.AlleleFrequencyRange((double)lowIndex / (double)N, (double)highIndex / (double)N, fraction); - } } \ No newline at end of file diff --git a/java/src/org/broadinstitute/sting/gatk/walkers/genotyper/UnifiedArgumentCollection.java b/java/src/org/broadinstitute/sting/gatk/walkers/genotyper/UnifiedArgumentCollection.java index 902caffc5..7b399d953 100755 --- a/java/src/org/broadinstitute/sting/gatk/walkers/genotyper/UnifiedArgumentCollection.java +++ b/java/src/org/broadinstitute/sting/gatk/walkers/genotyper/UnifiedArgumentCollection.java @@ -84,7 +84,4 @@ public class UnifiedArgumentCollection { @Argument(fullName = "min_allele_frequency", shortName = "min_freq", doc = "The minimum possible allele frequency in a population (advanced)", required = false) public double MINIMUM_ALLELE_FREQUENCY = 1e-8; - - @Argument(fullName = "allele_frequency_range", shortName = "freq_range", doc = "The range/fraction to emit of the total probability over all frequencies (in JOINT_ESTIMATION model only)", required = false) - public double ALLELE_FREQUENCY_RANGE = 0.95; } diff --git a/java/src/org/broadinstitute/sting/utils/genotype/AlleleFrequencyBacked.java b/java/src/org/broadinstitute/sting/utils/genotype/AlleleFrequencyBacked.java index 0c413f0f8..043ef8896 100755 --- a/java/src/org/broadinstitute/sting/utils/genotype/AlleleFrequencyBacked.java +++ b/java/src/org/broadinstitute/sting/utils/genotype/AlleleFrequencyBacked.java @@ -20,40 +20,4 @@ public interface AlleleFrequencyBacked { * @param frequency the allele frequency for this genotype */ public void setAlleleFrequency(double frequency); - - /** - * - * @return returns the allele frequency for this genotype - */ - public AlleleFrequencyRange getAlleleFrequencyRange(); - - /** - * - * @param range the allele frequency range for this genotype - */ - public void setAlleleFrequencyRange(AlleleFrequencyRange range); - - - /** - * A class representing a range of allele frequencies that make up a given - * fraction of the total probability over all frequencies. - */ - public class AlleleFrequencyRange { - - private double mLow, mHigh, mFraction; - - public AlleleFrequencyRange(double low, double high, double fraction) { - mLow = low; - mHigh = high; - mFraction = fraction; - } - - public double getLowEnd() { return mLow; } - public double getHighEnd() { return mHigh; } - public double getPercentageOfProbability() { return mFraction; } - - public String toString() { - return String.format("%.2f-%.2f,%.0f%%", mLow, mHigh, 100*mFraction); - } - } } \ No newline at end of file diff --git a/java/src/org/broadinstitute/sting/utils/genotype/vcf/VCFGenotypeLocusData.java b/java/src/org/broadinstitute/sting/utils/genotype/vcf/VCFGenotypeLocusData.java index 674b7ecdd..83c80ea07 100755 --- a/java/src/org/broadinstitute/sting/utils/genotype/vcf/VCFGenotypeLocusData.java +++ b/java/src/org/broadinstitute/sting/utils/genotype/vcf/VCFGenotypeLocusData.java @@ -21,9 +21,8 @@ public class VCFGenotypeLocusData implements GenotypeLocusData, ConfidenceBacked // the strand score lod private double mSLOD = 0.0; - // the allele frequency fields + // the allele frequency field private double mAlleleFrequency = 0.0; - private AlleleFrequencyRange mAFrange = null; // the location private GenomeLoc mLoc; @@ -133,23 +132,6 @@ public class VCFGenotypeLocusData implements GenotypeLocusData, ConfidenceBacked mAlleleFrequency = frequency; } - /** - * get the allele frequency range - * - * @return the allele frequency range - */ - public AlleleFrequencyRange getAlleleFrequencyRange() { - return mAFrange; - } - - /** - * - * @param range the allele frequency range for this genotype - */ - public void setAlleleFrequencyRange(AlleleFrequencyRange range) { - mAFrange = range; - } - /** * @return returns the dbsnp id for this genotype */ diff --git a/java/src/org/broadinstitute/sting/utils/genotype/vcf/VCFGenotypeWriterAdapter.java b/java/src/org/broadinstitute/sting/utils/genotype/vcf/VCFGenotypeWriterAdapter.java index de5fe6198..afb978e06 100644 --- a/java/src/org/broadinstitute/sting/utils/genotype/vcf/VCFGenotypeWriterAdapter.java +++ b/java/src/org/broadinstitute/sting/utils/genotype/vcf/VCFGenotypeWriterAdapter.java @@ -189,9 +189,6 @@ public class VCFGenotypeWriterAdapter implements GenotypeWriter { if ( locusdata != null ) { infoFields.put("SB", String.format("%.2f", locusdata.getSLOD())); infoFields.put("AF", String.format("%.2f", locusdata.getAlleleFrequency())); - AlleleFrequencyBacked.AlleleFrequencyRange range = locusdata.getAlleleFrequencyRange(); - if ( range != null ) - infoFields.put("AFrange", range.toString()); Map otherFields = locusdata.getFields(); if ( otherFields != null ) { infoFields.putAll(otherFields); diff --git a/java/test/org/broadinstitute/sting/gatk/walkers/genotyper/UnifiedGenotyperIntegrationTest.java b/java/test/org/broadinstitute/sting/gatk/walkers/genotyper/UnifiedGenotyperIntegrationTest.java index 4023fdb09..910d960eb 100755 --- a/java/test/org/broadinstitute/sting/gatk/walkers/genotyper/UnifiedGenotyperIntegrationTest.java +++ b/java/test/org/broadinstitute/sting/gatk/walkers/genotyper/UnifiedGenotyperIntegrationTest.java @@ -68,7 +68,7 @@ public class UnifiedGenotyperIntegrationTest extends WalkerTest { public void testMultiSamplePilot1Joint() { WalkerTest.WalkerTestSpec spec = new WalkerTest.WalkerTestSpec( "-T UnifiedGenotyper -R /broad/1KG/reference/human_b36_both.fasta -I /humgen/gsa-scr1/GATK_Data/Validation_Data/low_coverage_CEU.chr1.10k-11k.bam -varout %s -L 1:10,023,400-10,024,000 -bm empirical -gm JOINT_ESTIMATE -confidence 30", 1, - Arrays.asList("6f2753dbbf6b4ff98a75510e7ff2d5ae")); + Arrays.asList("a80abdc63117564ccd65a13ee9ca46be")); executeTest("testMultiSamplePilot1 - Joint Estimate", spec); } @@ -76,7 +76,7 @@ public class UnifiedGenotyperIntegrationTest extends WalkerTest { public void testMultiSamplePilot2Joint() { WalkerTest.WalkerTestSpec spec = new WalkerTest.WalkerTestSpec( "-T UnifiedGenotyper -R /broad/1KG/reference/human_b36_both.fasta -I /humgen/gsa-scr1/GATK_Data/Validation_Data/pilot2_daughters.chr20.10k-11k.bam -varout %s -L 20:10,000,000-10,010,000 -bm empirical -gm JOINT_ESTIMATE -confidence 30", 1, - Arrays.asList("979c1d88f2464548e83027daf685c2ad")); + Arrays.asList("47926784f90123ab058d36e12922a1ee")); executeTest("testMultiSamplePilot2 - Joint Estimate", spec); }