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

Merged bug fix from Stable into Unstable. Updating md5s 

Conflicts:
	protected/java/test/org/broadinstitute/sting/gatk/walkers/genotyper/UnifiedGenotyperIntegrationTest.java
This commit is contained in:
Mark DePristo 2012-11-28 12:54:36 -05:00
parent fb1e525e0f a1d6461121
commit c676853731
4 changed files with 69 additions and 12 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

6
protected/java/test/org/broadinstitute/sting/gatk/walkers/genotyper/UnifiedGenotyperIntegrationTest.java
View File

@ -62,7 +62,7 @@ public class UnifiedGenotyperIntegrationTest extends WalkerTest {
public void testMultipleSNPAlleles() {
WalkerTest.WalkerTestSpec spec = new WalkerTest.WalkerTestSpec(
"-T UnifiedGenotyper -R " + b37KGReference + " --no_cmdline_in_header -glm BOTH --dbsnp " + b37dbSNP129 + " -I " + privateTestDir + "multiallelic.snps.bam -o %s -L " + privateTestDir + "multiallelic.snps.intervals", 1,
Arrays.asList("704888987baacff8c7b273b8ab9938d0"));
Arrays.asList("d20c7a143b899f0239bf64b652ad3edb"));
executeTest("test Multiple SNP alleles", spec);
}
@ -197,7 +197,7 @@ public class UnifiedGenotyperIntegrationTest extends WalkerTest {
@Test
public void testOutputParameterAllSites() {
testOutputParameters("--output_mode EMIT_ALL_SITES", "81fff490c0f59890f1e75dc290833434");
testOutputParameters("--output_mode EMIT_ALL_SITES", "8b26088a035e579c4afd3b46737291e4");
}
private void testOutputParameters(final String args, final String md5) {
@ -345,7 +345,7 @@ public class UnifiedGenotyperIntegrationTest extends WalkerTest {
public void testMultiSampleIndels1() {
WalkerTest.WalkerTestSpec spec1 = new WalkerTest.WalkerTestSpec(
baseCommandIndels + " -I " + validationDataLocation + "low_coverage_CEU.chr1.10k-11k.bam -o %s -L 1:10450700-10551000", 1,
Arrays.asList("a4761d7f25e7a62f34494801c98a0da7"));
Arrays.asList("69df7a00f800204564ca3726e1871132"));
List<File> result = executeTest("test MultiSample Pilot1 CEU indels", spec1).getFirst();
WalkerTest.WalkerTestSpec spec2 = new WalkerTest.WalkerTestSpec(

56
protected/java/test/org/broadinstitute/sting/gatk/walkers/genotyper/afcalc/AFCalcResultUnitTest.java
View File

@ -2,16 +2,14 @@ package org.broadinstitute.sting.gatk.walkers.genotyper.afcalc;
import org.broadinstitute.sting.BaseTest;
import org.broadinstitute.sting.utils.MathUtils;
import org.broadinstitute.sting.utils.QualityUtils;
import org.broadinstitute.sting.utils.Utils;
import org.broadinstitute.sting.utils.variantcontext.Allele;
import org.testng.Assert;
import org.testng.annotations.DataProvider;
import org.testng.annotations.Test;
import java.util.ArrayList;
import java.util.Arrays;
import java.util.Collections;
import java.util.List;
import java.util.*;
public class AFCalcResultUnitTest extends BaseTest {
private static class MyTest {
@ -79,4 +77,54 @@ public class AFCalcResultUnitTest extends BaseTest {
final double[] actualPosteriors = new double[]{result.getLog10PosteriorOfAFEq0(), result.getLog10PosteriorOfAFGT0()};
Assert.assertEquals(MathUtils.sumLog10(actualPosteriors), 1.0, 1e-3, "Posteriors don't sum to 1 with 1e-3 precision");
}
@DataProvider(name = "TestIsPolymorphic")
public Object[][] makeTestIsPolymorphic() {
List<Object[]> tests = new ArrayList<Object[]>();
final List<Double> pValues = new LinkedList<Double>();
for ( final double p : Arrays.asList(0.01, 0.1, 0.9, 0.99, 0.999) )
for ( final double espilon : Arrays.asList(-1e-5, 0.0, 1e-5) )
pValues.add(p + espilon);
for ( final double pNonRef : pValues ) {
for ( final double pThreshold : pValues ) {
final boolean shouldBePoly = pNonRef >= pThreshold;
if ( pNonRef != pThreshold)
// let's not deal with numerical instability
tests.add(new Object[]{ pNonRef, pThreshold, shouldBePoly });
}
}
return tests.toArray(new Object[][]{});
}
private AFCalcResult makePolymorphicTestData(final double pNonRef) {
return new AFCalcResult(
new int[]{0},
1,
alleles,
MathUtils.normalizeFromLog10(new double[]{1 - pNonRef, pNonRef}, true, false),
log10Even,
Collections.singletonMap(C, Math.log10(pNonRef)));
}
@Test(enabled = true, dataProvider = "TestIsPolymorphic")
private void testIsPolymorphic(final double pNonRef, final double pThreshold, final boolean shouldBePoly) {
final AFCalcResult result = makePolymorphicTestData(pNonRef);
final boolean actualIsPoly = result.isPolymorphic(C, Math.log10(pThreshold));
Assert.assertEquals(actualIsPoly, shouldBePoly,
"isPolymorphic with pNonRef " + pNonRef + " and threshold " + pThreshold + " returned "
+ actualIsPoly + " but the expected result is " + shouldBePoly);
}
@Test(enabled = true, dataProvider = "TestIsPolymorphic")
private void testIsPolymorphicQual(final double pNonRef, final double pThreshold, final boolean shouldBePoly) {
final AFCalcResult result = makePolymorphicTestData(pNonRef);
final double qual = QualityUtils.phredScaleCorrectRate(pThreshold);
final boolean actualIsPoly = result.isPolymorphicPhredScaledQual(C, qual);
Assert.assertEquals(actualIsPoly, shouldBePoly,
"isPolymorphic with pNonRef " + pNonRef + " and threshold " + pThreshold + " returned "
+ actualIsPoly + " but the expected result is " + shouldBePoly);
}
}

8
public/java/src/org/broadinstitute/sting/gatk/walkers/genotyper/UnifiedGenotyperEngine.java
View File

@ -382,11 +382,9 @@ public class UnifiedGenotyperEngine {
if ( alternateAllele.isReference() )
continue;
// we are non-ref if the probability of being non-ref > the emit confidence.
// the emit confidence is phred-scaled, say 30 => 10^-3.
// the posterior AF > 0 is log10: -5 => 10^-5
// we are non-ref if 10^-5 < 10^-3 => -5 < -3
final boolean isNonRef = AFresult.isPolymorphic(alternateAllele, UAC.STANDARD_CONFIDENCE_FOR_EMITTING / -10.0);
// Compute if the site is considered polymorphic with sufficient confidence relative to our
// phred-scaled emission QUAL
final boolean isNonRef = AFresult.isPolymorphicPhredScaledQual(alternateAllele, UAC.STANDARD_CONFIDENCE_FOR_EMITTING);
// if the most likely AC is not 0, then this is a good alternate allele to use
if ( isNonRef ) {

11
public/java/src/org/broadinstitute/sting/gatk/walkers/genotyper/afcalc/AFCalcResult.java
View File

@ -28,6 +28,7 @@ package org.broadinstitute.sting.gatk.walkers.genotyper.afcalc;
import com.google.java.contract.Ensures;
import com.google.java.contract.Requires;
import org.broadinstitute.sting.utils.MathUtils;
import org.broadinstitute.sting.utils.QualityUtils;
import org.broadinstitute.sting.utils.Utils;
import org.broadinstitute.sting.utils.variantcontext.Allele;
@ -234,10 +235,20 @@ public class AFCalcResult {
*
* @return true if there's enough confidence (relative to log10minPNonRef) to reject AF == 0
*/
@Requires("MathUtils.goodLog10Probability(log10minPNonRef)")
public boolean isPolymorphic(final Allele allele, final double log10minPNonRef) {
return getLog10PosteriorOfAFGt0ForAllele(allele) >= log10minPNonRef;
}
/**
* Same as #isPolymorphic but takes a phred-scaled quality score as input
*/
public boolean isPolymorphicPhredScaledQual(final Allele allele, final double minPNonRefPhredScaledQual) {
if ( minPNonRefPhredScaledQual < 0 ) throw new IllegalArgumentException("phredScaledQual " + minPNonRefPhredScaledQual + " < 0 ");
final double log10Threshold = Math.log10(QualityUtils.qualToProb(minPNonRefPhredScaledQual));
return isPolymorphic(allele, log10Threshold);
}
/**
* Are any of the alleles polymorphic w.r.t. #isPolymorphic?
*