zzh
/
gatk-3.8
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

subtract largest posterior value from all values; this hopefully solves any precision issues 

git-svn-id: file:///humgen/gsa-scr1/gsa-engineering/svn_contents/trunk@1870 348d0f76-0448-11de-a6fe-93d51630548a
This commit is contained in:
ebanks 2009-10-18 05:20:15 +00:00
parent b9e8867287
commit 51f9ec0a5c
1 changed files with 10 additions and 1 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

11
java/src/org/broadinstitute/sting/gatk/walkers/genotyper/GenotypeLikelihoods.java
View File

@ -6,6 +6,7 @@ import org.broadinstitute.sting.utils.genotype.DiploidGenotype;
import static java.lang.Math.log10; import static java.lang.Math.log10;
import static java.lang.Math.pow; import static java.lang.Math.pow;
import java.util.Arrays;
/** /**
* Stable, error checking version of the Bayesian genotyper. Useful for calculating the likelihoods, priors, * Stable, error checking version of the Bayesian genotyper. Useful for calculating the likelihoods, priors,
@ -158,9 +159,17 @@ public abstract class GenotypeLikelihoods implements Cloneable {
double[] normalized = new double[posteriors.length]; double[] normalized = new double[posteriors.length];
double sum = 0.0; double sum = 0.0;
// for precision purposes, we need to add (or really subtract, since everything is negative)
// the largest posterior value from all entries so that numbers don't get too small
double maxValue = posteriors[0];
for (int i = 1; i < posteriors.length; i++) {
if ( maxValue < posteriors[i] )
maxValue = posteriors[i];
}
// collect the posteriors // collect the posteriors
for ( DiploidGenotype g : DiploidGenotype.values() ) { for ( DiploidGenotype g : DiploidGenotype.values() ) {
double posterior = Math.pow(10, getPosterior(g)); double posterior = Math.pow(10, getPosterior(g) - maxValue);
normalized[g.ordinal()] = posterior; normalized[g.ordinal()] = posterior;
sum += posterior; sum += posterior;
} }