From e920badcc41a1407b4a790e217818f0ab262abd9 Mon Sep 17 00:00:00 2001 From: delangel Date: Fri, 1 Oct 2010 17:43:43 +0000 Subject: [PATCH] Temporary fix for case where genotype likelihoods are exactly (1,0,0) or (0,1,0) etc. at a site with new indel genotyper: this would make us blow up when converting to log space and try to assign genotypes at a site. A more robust solution is in the works. git-svn-id: file:///humgen/gsa-scr1/gsa-engineering/svn_contents/trunk@4401 348d0f76-0448-11de-a6fe-93d51630548a --- .../DindelGenotypeLikelihoodsCalculationModel.java | 9 ++++++++- 1 file changed, 8 insertions(+), 1 deletion(-) diff --git a/java/src/org/broadinstitute/sting/playground/gatk/walkers/genotyper/DindelGenotypeLikelihoodsCalculationModel.java b/java/src/org/broadinstitute/sting/playground/gatk/walkers/genotyper/DindelGenotypeLikelihoodsCalculationModel.java index 5a73a54a2..1da0f8221 100755 --- a/java/src/org/broadinstitute/sting/playground/gatk/walkers/genotyper/DindelGenotypeLikelihoodsCalculationModel.java +++ b/java/src/org/broadinstitute/sting/playground/gatk/walkers/genotyper/DindelGenotypeLikelihoodsCalculationModel.java @@ -47,6 +47,7 @@ public class DindelGenotypeLikelihoodsCalculationModel extends GenotypeLikelihoo private final double insertionEndProbability = 0.5; private final double alphaDeletionProbability = 1e-3; private final int HAPLOTYPE_SIZE = 80; + private static final double MINUS_INFINITY = -1e200; // todo - the following need to be exposed for command line argument control private final double indelHeterozygosity = 1.0/8000; @@ -125,9 +126,15 @@ public class DindelGenotypeLikelihoodsCalculationModel extends GenotypeLikelihoo double[] genotypeLikelihoods = HaplotypeIndelErrorModel.getPosteriorProbabilitesFromHaplotypeLikelihoods( haplotypeLikehoodMatrix); + // todo- cleaner solution for case where probability is of form (1,0,0) or similar + for (int k=0; k < 3; k++) { + genotypeLikelihoods[k] = Math.log10(genotypeLikelihoods[k]); + if (Double.isInfinite(genotypeLikelihoods[k])) + genotypeLikelihoods[k] = -MINUS_INFINITY; + } GLs.put(sample.getKey(), new BiallelicGenotypeLikelihoods(sample.getKey(),vc.getReference(), vc.getAlternateAllele(0), - Math.log10(genotypeLikelihoods[0]),Math.log10(genotypeLikelihoods[1]), Math.log10(genotypeLikelihoods[2]))); + genotypeLikelihoods[0],genotypeLikelihoods[1], genotypeLikelihoods[2])); } }