For Guillermo: only decide that something is a clear reference call if it is at least 10 times as likely as the next best genotype

git-svn-id: file:///humgen/gsa-scr1/gsa-engineering/svn_contents/trunk@4441 348d0f76-0448-11de-a6fe-93d51630548a
This commit is contained in:
ebanks 2010-10-06 15:16:41 +00:00
parent effcd26977
commit 3c5dc675ab
2 changed files with 6 additions and 4 deletions

View File

@ -52,6 +52,8 @@ public abstract class AlleleFrequencyCalculationModel implements Cloneable {
GRID_SEARCH
}
private static final double LOG10_REFERENCE_CALL_EPSILON = 1.0;
protected int N;
protected AlleleFrequencyMatrix AFMatrix;
protected Set<BiallelicGenotypeLikelihoods> refCalls;
@ -179,8 +181,7 @@ public abstract class AlleleFrequencyCalculationModel implements Cloneable {
AFMatrix.clear();
for ( BiallelicGenotypeLikelihoods GL : GLs.values() ) {
// todo - gdebug workaround for crash
if (false) { //isClearRefCall(GL) ) {
if ( isClearRefCall(GL) ) {
refCalls.add(GL);
} else {
AFMatrix.setLikelihoods(GL.getPosteriors(), GL.getSample());
@ -193,7 +194,8 @@ public abstract class AlleleFrequencyCalculationModel implements Cloneable {
return false;
double[] likelihoods = GL.getLikelihoods();
return ( likelihoods[0] > likelihoods[1] && likelihoods[0] > likelihoods[2]);
double refLikelihoodMinusEpsilon = likelihoods[0] - LOG10_REFERENCE_CALL_EPSILON;
return ( refLikelihoodMinusEpsilon > likelihoods[1] && refLikelihoodMinusEpsilon > likelihoods[2]);
}
protected class CalculatedAlleleFrequency {

View File

@ -351,7 +351,7 @@ public class UnifiedGenotyperEngine {
AFline.append(String.format("%.2f\t", ((float)i)/N));
AFline.append(String.format("%.8f\t", log10AlleleFrequencyPriors[i]));
if ( log10AlleleFrequencyPosteriors.get()[i] == AlleleFrequencyCalculationModel.VALUE_NOT_CALCULATED)
AFline.append("0.0\t");
AFline.append("0.00000000\t");
else
AFline.append(String.format("%.8f\t", log10AlleleFrequencyPosteriors.get()[i]));
AFline.append(String.format("%.8f\t", normalizedPosteriors[i]));