Fixed bug where k=2N case wasn't properly being computed. Added optimization for BB genotype case not in old model. At this point, integration tests pass except for 1 case where QUALs differ by 0.01 (this is okay because I occasionally need to compute extra cells in the matrix which affects the approximations) and 2 cases where multi-allelic indels are being genotyped (some work still needs to be done to support them).

This commit is contained in:
Eric Banks 2011-12-03 23:12:04 -05:00
parent 71f793b71b
commit 29662be3d7
2 changed files with 8 additions and 9 deletions

View File

@ -317,11 +317,11 @@ public class ExactAFCalculationModel extends AlleleFrequencyCalculationModel {
}
}
public int linearExactMultiAllelic(GenotypesContext GLs,
int numAlternateAlleles,
double[] log10AlleleFrequencyPriors,
double[] log10AlleleFrequencyPosteriors,
boolean preserveData) {
static public int linearExactMultiAllelic(GenotypesContext GLs,
int numAlternateAlleles,
double[] log10AlleleFrequencyPriors,
double[] log10AlleleFrequencyPosteriors,
boolean preserveData) {
final ArrayList<double[]> genotypeLikelihoods = getGLs(GLs);
final int numSamples = genotypeLikelihoods.size()-1;
@ -350,8 +350,7 @@ public class ExactAFCalculationModel extends AlleleFrequencyCalculationModel {
maxLog10L = Math.max(maxLog10L, log10LofKs);
}
// TODO -- finish me
// TODO -- why do we need to return anything here?
return 0;
}
@ -471,7 +470,7 @@ public class ExactAFCalculationModel extends AlleleFrequencyCalculationModel {
log10ConformationLikelihoods[0] = MathUtils.log10Cache[2*j-totalK] + MathUtils.log10Cache[2*j-totalK-1] + set.log10Likelihoods[j-1] + gl[HOM_REF_INDEX];
// deal with the other possible conformations now
if ( totalK < 2*j ) {
if ( totalK <= 2*j ) { // skip impossible conformations
int conformationIndex = 1;
for ( Map.Entry<Integer, Integer> mapping : set.ACsetIndexToPLIndex.entrySet() )
log10ConformationLikelihoods[conformationIndex++] =

View File

@ -29,7 +29,7 @@ public class UnifiedGenotyperIntegrationTest extends WalkerTest {
public void testMultiSamplePilot1() {
WalkerTest.WalkerTestSpec spec = new WalkerTest.WalkerTestSpec(
baseCommand + " -I " + validationDataLocation + "low_coverage_CEU.chr1.10k-11k.bam -o %s -L 1:10,022,000-10,025,000", 1,
Arrays.asList("286f0de92e4ce57986ba861390c6019d"));
Arrays.asList("b70732a2f63f8409b61e41fa53eaae3e"));
executeTest("test MultiSample Pilot1", spec);
}