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

Updated VQSR integration tests because of the new Omni file. Fixed overflow condition in FisherStrand when the depth is too high. 

git-svn-id: file:///humgen/gsa-scr1/gsa-engineering/svn_contents/trunk@5910 348d0f76-0448-11de-a6fe-93d51630548a
This commit is contained in:
rpoplin 2011-06-01 14:20:37 +00:00
parent 6cc84c3ce2
commit 2b5683909e
3 changed files with 8 additions and 20 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

20
java/src/org/broadinstitute/sting/gatk/walkers/annotator/FisherStrand.java
View File

@ -49,42 +49,30 @@ public class FisherStrand implements InfoFieldAnnotation, WorkInProgressAnnotati
private static final String ALTFWD = "ALTFWD"; private static final String ALTFWD = "ALTFWD";
private static final String ALTREV = "ALTREV"; private static final String ALTREV = "ALTREV";
private static final String FS = "FS"; private static final String FS = "FS";
private static final double MIN_PVALUE = 1E-320;
public Map<String, Object> annotate(RefMetaDataTracker tracker, ReferenceContext ref, Map<String, AlignmentContext> stratifiedContexts, VariantContext vc) { public Map<String, Object> annotate(RefMetaDataTracker tracker, ReferenceContext ref, Map<String, AlignmentContext> stratifiedContexts, VariantContext vc) {
if ( ! vc.isVariant() || vc.isFiltered() || ! vc.isBiallelic() || ! vc.isSNP() ) if ( ! vc.isVariant() || vc.isFiltered() || ! vc.isBiallelic() || ! vc.isSNP() )
return null; return null;
int[][] table = getContingencyTable(stratifiedContexts, vc.getReference(), vc.getAlternateAllele(0)); int[][] table = getContingencyTable(stratifiedContexts, vc.getReference(), vc.getAlternateAllele(0));
Double pvalue = pValueForContingencyTable(table); Double pvalue = Math.max(pValueForContingencyTable(table), MIN_PVALUE);
if ( pvalue == null ) if ( pvalue == null )
return null; return null;
// use Math.abs to prevent -0's // use Math.abs to prevent -0's
Map<String, Object> map = new HashMap<String, Object>(); Map<String, Object> map = new HashMap<String, Object>();
int phredPValue = (int)Math.round(QualityUtils.phredScaleErrorRate(pvalue)); map.put(FS, String.format("%.3f", QualityUtils.phredScaleErrorRate(pvalue)));
addInt(map, REFFWD, table[0][0]);
addInt(map, REFREV, table[0][1]);
addInt(map, ALTFWD, table[1][0]);
addInt(map, ALTREV, table[1][1]);
addInt(map, FS, phredPValue);
return map; return map;
} }
private static void addInt(Map<String, Object> map, String key, int value) {
map.put(key, String.format("%d", value));
}
public List<String> getKeyNames() { public List<String> getKeyNames() {
return Arrays.asList(REFFWD,REFREV,ALTFWD,ALTREV,FS); return Arrays.asList(REFFWD,REFREV,ALTFWD,ALTREV,FS);
} }
public List<VCFInfoHeaderLine> getDescriptions() { public List<VCFInfoHeaderLine> getDescriptions() {
return Arrays.asList( return Arrays.asList(
new VCFInfoHeaderLine(REFFWD, 1, VCFHeaderLineType.Integer, "Count of bases with REF allele, forward strand"), new VCFInfoHeaderLine(FS, 1, VCFHeaderLineType.Float, "Phred-scaled p-value using Fisher's exact test to detect strand bias"));
new VCFInfoHeaderLine(REFREV, 1, VCFHeaderLineType.Integer, "Count of bases with REF allele, reverse strand"),
new VCFInfoHeaderLine(ALTFWD, 1, VCFHeaderLineType.Integer, "Count of bases with ALT allele, forward strand"),
new VCFInfoHeaderLine(ALTREV, 1, VCFHeaderLineType.Integer, "Count of bases with ALT allele, reverse strand"),
new VCFInfoHeaderLine(FS, 1, VCFHeaderLineType.Integer, "Integer-rounded Phred-scaled p-value using Fisher's exact test to detect strand bias"));
} }
private Double pValueForContingencyTable(int[][] originalTable) { private Double pValueForContingencyTable(int[][] originalTable) {

2
java/src/org/broadinstitute/sting/gatk/walkers/variantrecalibration/VariantDataManager.java
View File

@ -214,7 +214,7 @@ public class VariantDataManager {
} else { } else {
value = Double.parseDouble( (String)vc.getAttribute( annotationKey ) ); value = Double.parseDouble( (String)vc.getAttribute( annotationKey ) );
if( Double.isInfinite(value) ) { value = Double.NaN; } if( Double.isInfinite(value) ) { value = Double.NaN; }
if( jitter && ( annotationKey.equalsIgnoreCase("HRUN") || annotationKey.equalsIgnoreCase("FS") ) ) { // Integer valued annotations must be jittered a bit to work in this GMM if( jitter && ( annotationKey.equalsIgnoreCase("HRUN") ) ) { // Integer valued annotations must be jittered a bit to work in this GMM
value += -0.25 + 0.5 * GenomeAnalysisEngine.getRandomGenerator().nextDouble(); value += -0.25 + 0.5 * GenomeAnalysisEngine.getRandomGenerator().nextDouble();
} }
if( annotationKey.equals("HaplotypeScore") && MathUtils.compareDoubles(value, 0.0, 0.0001) == 0 ) { value = -0.2 + 0.4*GenomeAnalysisEngine.getRandomGenerator().nextDouble(); } if( annotationKey.equals("HaplotypeScore") && MathUtils.compareDoubles(value, 0.0, 0.0001) == 0 ) { value = -0.2 + 0.4*GenomeAnalysisEngine.getRandomGenerator().nextDouble(); }

6
java/test/org/broadinstitute/sting/gatk/walkers/variantrecalibration/VariantRecalibrationWalkersIntegrationTest.java
View File

@ -25,9 +25,9 @@ public class VariantRecalibrationWalkersIntegrationTest extends WalkerTest {
} }
VRTest lowPass = new VRTest("phase1.projectConsensus.chr20.raw.snps.vcf", VRTest lowPass = new VRTest("phase1.projectConsensus.chr20.raw.snps.vcf",
"920b12d7765eb4f6f4a1bab045679b31", // tranches "d33212a84368e821cbedecd4f59756d6", // tranches
"41bbc5f07c8a9573d5bb638f01808bba", // recal file "a35cd067f378442eee8cd5edeea92be0", // recal file
"d9aa7a0f8fb886df4394f1636605adca"); // cut VCF "7259b7daefe57b11ae9e537e38569160"); // cut VCF
@DataProvider(name = "VRTest") @DataProvider(name = "VRTest")
public Object[][] createData1() { public Object[][] createData1() {