From 9648399630797d914d2b12f94c3fe2f495d2c894 Mon Sep 17 00:00:00 2001 From: delangel Date: Mon, 10 Jan 2011 14:56:28 +0000 Subject: [PATCH] Boneheaded silly bug in indel caller - posterior probability computation was using priors gotten from SNP heterozygosity, not indel heterozygosity. Added then indel het. argument to command line and hook it up (not a radical change in calls though, just a few dubious calls around the edges fall off) git-svn-id: file:///humgen/gsa-scr1/gsa-engineering/svn_contents/trunk@4967 348d0f76-0448-11de-a6fe-93d51630548a --- .../walkers/genotyper/UnifiedArgumentCollection.java | 6 ++++++ .../gatk/walkers/genotyper/UnifiedGenotyperEngine.java | 9 ++++++++- 2 files changed, 14 insertions(+), 1 deletion(-) 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 50e326346..89f7d3377 100755 --- a/java/src/org/broadinstitute/sting/gatk/walkers/genotyper/UnifiedArgumentCollection.java +++ b/java/src/org/broadinstitute/sting/gatk/walkers/genotyper/UnifiedArgumentCollection.java @@ -92,6 +92,9 @@ public class UnifiedArgumentCollection { @Argument(fullName = "min_indel_count_for_genotyping", shortName = "minIndelCnt", doc = "Minimum number of consensus indels required to trigger genotyping run", required = false) public int MIN_INDEL_COUNT_FOR_GENOTYPING = 5; + @Argument(fullName = "indel_heterozygosity", shortName = "indelHeterozygosity", doc = "Heterozygosity for indel calling", required = false) + public double INDEL_HETEROZYGOSITY = 1.0/8000; + public UnifiedArgumentCollection clone() { UnifiedArgumentCollection uac = new UnifiedArgumentCollection(); @@ -110,6 +113,9 @@ public class UnifiedArgumentCollection { uac.MAX_MISMATCHES = MAX_MISMATCHES; uac.USE_BADLY_MATED_READS = USE_BADLY_MATED_READS; uac.MAX_DELETION_FRACTION = MAX_DELETION_FRACTION; + uac.GET_ALLELES_FROM_VCF = GET_ALLELES_FROM_VCF; + uac.MIN_INDEL_COUNT_FOR_GENOTYPING = MIN_INDEL_COUNT_FOR_GENOTYPING; + uac.INDEL_HETEROZYGOSITY = INDEL_HETEROZYGOSITY; return uac; } diff --git a/java/src/org/broadinstitute/sting/gatk/walkers/genotyper/UnifiedGenotyperEngine.java b/java/src/org/broadinstitute/sting/gatk/walkers/genotyper/UnifiedGenotyperEngine.java index b8a4b09de..7e9574a4f 100755 --- a/java/src/org/broadinstitute/sting/gatk/walkers/genotyper/UnifiedGenotyperEngine.java +++ b/java/src/org/broadinstitute/sting/gatk/walkers/genotyper/UnifiedGenotyperEngine.java @@ -629,8 +629,15 @@ public class UnifiedGenotyperEngine { protected void computeAlleleFrequencyPriors(int N) { // calculate the allele frequency priors for 1-N double sum = 0.0; + double heterozygosity; + + if (UAC.GLmodel == GenotypeLikelihoodsCalculationModel.Model.DINDEL) + heterozygosity = UAC.INDEL_HETEROZYGOSITY; + else + heterozygosity = UAC.heterozygosity; + for (int i = 1; i <= N; i++) { - double value = UAC.heterozygosity / (double)i; + double value = heterozygosity / (double)i; log10AlleleFrequencyPriors[i] = Math.log10(value); sum += value; }