From a19389528ddafb7d1f3e269ad704529cff19a6b5 Mon Sep 17 00:00:00 2001 From: delangel Date: Tue, 3 May 2011 22:38:33 +0000 Subject: [PATCH] Bring back from the dead the old likelihoods model for indels, which has worse performance but is about 4x faster. Enabled with argument -GSA_PRODUCTION_ONLY in UG git-svn-id: file:///humgen/gsa-scr1/gsa-engineering/svn_contents/trunk@5748 348d0f76-0448-11de-a6fe-93d51630548a --- ...elGenotypeLikelihoodsCalculationModel.java | 35 +++++++++++++++++-- .../genotyper/UnifiedArgumentCollection.java | 4 +++ 2 files changed, 36 insertions(+), 3 deletions(-) diff --git a/java/src/org/broadinstitute/sting/gatk/walkers/genotyper/IndelGenotypeLikelihoodsCalculationModel.java b/java/src/org/broadinstitute/sting/gatk/walkers/genotyper/IndelGenotypeLikelihoodsCalculationModel.java index 12f55e60f..2a5ebf027 100755 --- a/java/src/org/broadinstitute/sting/gatk/walkers/genotyper/IndelGenotypeLikelihoodsCalculationModel.java +++ b/java/src/org/broadinstitute/sting/gatk/walkers/genotyper/IndelGenotypeLikelihoodsCalculationModel.java @@ -56,12 +56,33 @@ public class IndelGenotypeLikelihoodsCalculationModel extends GenotypeLikelihood private boolean DEBUG = false; private PairHMMIndelErrorModel pairModel; - +// gdebug removeme +// todo -cleanup +private HaplotypeIndelErrorModel model; + private boolean useOldWrongHorribleHackedUpLikelihoodModel = false; +// private GenomeLoc lastSiteVisited; private List alleleList; protected IndelGenotypeLikelihoodsCalculationModel(UnifiedArgumentCollection UAC, Logger logger) { super(UAC, logger); + if (UAC.GSA_PRODUCTION_ONLY == false) { + pairModel = new PairHMMIndelErrorModel(UAC.INDEL_GAP_OPEN_PENALTY,UAC.INDEL_GAP_CONTINUATION_PENALTY, + UAC.OUTPUT_DEBUG_INDEL_INFO, UAC.DO_CONTEXT_DEPENDENT_PENALTIES, UAC.dovit, UAC.GET_GAP_PENALTIES_FROM_DATA, UAC.INDEL_RECAL_FILE); + useOldWrongHorribleHackedUpLikelihoodModel = false; + } + else { + useOldWrongHorribleHackedUpLikelihoodModel = true; + double INSERTION_START_PROBABILITY = 1e-3; + + double INSERTION_END_PROBABILITY = 0.5; + + double ALPHA_DELETION_PROBABILITY = 1e-3; + + + model = new HaplotypeIndelErrorModel(3, INSERTION_START_PROBABILITY, + INSERTION_END_PROBABILITY,ALPHA_DELETION_PROBABILITY,UAC.INDEL_HAPLOTYPE_SIZE, false, UAC.OUTPUT_DEBUG_INDEL_INFO); + } pairModel = new PairHMMIndelErrorModel(UAC.INDEL_GAP_OPEN_PENALTY,UAC.INDEL_GAP_CONTINUATION_PENALTY, UAC.OUTPUT_DEBUG_INDEL_INFO, UAC.DO_CONTEXT_DEPENDENT_PENALTIES, UAC.dovit, UAC.GET_GAP_PENALTIES_FROM_DATA, UAC.INDEL_RECAL_FILE); alleleList = new ArrayList(); @@ -322,7 +343,12 @@ public class IndelGenotypeLikelihoodsCalculationModel extends GenotypeLikelihood List haplotypesInVC; int hsize = (int)ref.getWindow().size()-Math.abs(eventLength)-1; - int numPrefBases = ref.getLocus().getStart()-ref.getWindow().getStart()+1; + int numPrefBases= ref.getLocus().getStart()-ref.getWindow().getStart()+1; + + if (useOldWrongHorribleHackedUpLikelihoodModel) { + numPrefBases = 20; + hsize=80; + } if (DEBUG) System.out.format("hsize: %d eventLength: %d refSize: %d, locStart: %d numpr: %d\n",hsize,eventLength, (int)ref.getWindow().size(), loc.getStart(), numPrefBases); @@ -348,8 +374,11 @@ public class IndelGenotypeLikelihoodsCalculationModel extends GenotypeLikelihood pileup = context.getBasePileup(); if (pileup != null ) { + if (useOldWrongHorribleHackedUpLikelihoodModel) + haplotypeLikehoodMatrix = model.computeReadHaplotypeLikelihoods( pileup, haplotypesInVC); + else + haplotypeLikehoodMatrix = pairModel.computeReadHaplotypeLikelihoods( pileup, haplotypesInVC, ref, HAPLOTYPE_SIZE, eventLength); - haplotypeLikehoodMatrix = pairModel.computeReadHaplotypeLikelihoods( pileup, haplotypesInVC, ref, HAPLOTYPE_SIZE, eventLength); double[] genotypeLikelihoods = HaplotypeIndelErrorModel.getHaplotypeLikelihoods( haplotypeLikehoodMatrix); 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 67be9cfea..7e38b5222 100755 --- a/java/src/org/broadinstitute/sting/gatk/walkers/genotyper/UnifiedArgumentCollection.java +++ b/java/src/org/broadinstitute/sting/gatk/walkers/genotyper/UnifiedArgumentCollection.java @@ -117,6 +117,9 @@ public class UnifiedArgumentCollection { @Argument(fullName = "dovit", shortName = "dovit", doc = "Output indel debug info", required = false) public boolean dovit = false; @Hidden + @Argument(fullName = "GSA_PRODUCTION_ONLY", shortName = "GSA_PRODUCTION_ONLY", doc = "don't ever use me", required = false) + public boolean GSA_PRODUCTION_ONLY = false; + @Hidden @Argument(fullName = "exactCalculation", shortName = "exactCalculation", doc = "expt", required = false) public ExactAFCalculationModel.ExactCalculation EXACT_CALCULATION_TYPE = ExactAFCalculationModel.ExactCalculation.LINEAR_EXPERIMENTAL; @@ -161,6 +164,7 @@ public class UnifiedArgumentCollection { // todo- arguments to remove uac.COVERAGE_AT_WHICH_TO_ABORT = COVERAGE_AT_WHICH_TO_ABORT; uac.dovit = dovit; + uac.GSA_PRODUCTION_ONLY = GSA_PRODUCTION_ONLY; return uac; }