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

Merge pull request #1358 from broadinstitute/sl_issue_1327 

Changed calls with GQ=0 to no-call for HaplotypeCaller in normal mode.
This commit is contained in:
samuelklee 2016-05-19 12:26:35 -04:00
parent edd9fcb7e5 bf4b1a5421
commit 7fdd3c2a0c
4 changed files with 41 additions and 6 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

29
protected/gatk-tools-protected/src/main/java/org/broadinstitute/gatk/tools/walkers/haplotypecaller/HaplotypeCallerGenotypingEngine.java
View File

@ -238,8 +238,6 @@ public class HaplotypeCallerGenotypingEngine extends GenotypingEngine<AssemblyBa
if( mergedVC == null )
continue;
final GenotypeLikelihoodsCalculationModel.Model calculationModel = mergedVC.isSNP()
? GenotypeLikelihoodsCalculationModel.Model.SNP : GenotypeLikelihoodsCalculationModel.Model.INDEL;
@ -286,6 +284,11 @@ public class HaplotypeCallerGenotypingEngine extends GenotypingEngine<AssemblyBa
calledHaplotypes.addAll(haplotypeList);
}
if ( !emitReferenceConfidence ) {
// set GTs to no-call when GQ is 0 in normal mode
annotatedCall = clearUnconfidentGenotypeCalls(annotatedCall);
}
returnCalls.add( annotatedCall );
}
}
@ -809,4 +812,26 @@ public class HaplotypeCallerGenotypingEngine extends GenotypingEngine<AssemblyBa
public GenotypingModel getGenotypingModel() {
return genotypingModel;
}
/**
* Cleans up genotype-level annotations that need to be updated
* (similar to {@link org.broadinstitute.gatk.tools.walkers.variantutils.GenotypeGVCFs#cleanupGenotypeAnnotations},
* but here we only check for unconfident calls (i.e., those with GQ = 0), which are set to no-call).
*/
private VariantContext clearUnconfidentGenotypeCalls(final VariantContext VC) {
final GenotypesContext oldGTs = VC.getGenotypes();
final List<Genotype> clearedGTs = new ArrayList<>(oldGTs.size());
for ( final Genotype oldGT : oldGTs ) {
// set GT to no-call when GQ is 0
if (oldGT.hasGQ() && oldGT.getGQ() == 0) {
final int ploidy = oldGT.getPloidy();
final List<Allele> noCallAlleles = GATKVariantContextUtils.noCallAlleles(ploidy);
final Genotype noCallGT = new GenotypeBuilder().alleles(noCallAlleles).make();
clearedGTs.add(noCallGT);
} else {
clearedGTs.add(oldGT);
}
}
return new VariantContextBuilder(VC).genotypes(clearedGTs).make();
}
}

2
protected/gatk-tools-protected/src/test/java/org/broadinstitute/gatk/tools/walkers/haplotypecaller/HaplotypeCallerComplexAndSymbolicVariantsIntegrationTest.java
View File

@ -72,7 +72,7 @@ public class HaplotypeCallerComplexAndSymbolicVariantsIntegrationTest extends Wa
@Test
public void testHaplotypeCallerMultiSampleComplex1() {
HCTestComplexVariants(privateTestDir + "AFR.complex.variants.bam", "", "32ce23b3830f5f2c693161b40de8b15e");
HCTestComplexVariants(privateTestDir + "AFR.complex.variants.bam", "", "88255eda0e29e4a6e128ddb7177a03ab");
}
private void HCTestSymbolicVariants(String bam, String args, String md5) {

6
protected/gatk-tools-protected/src/test/java/org/broadinstitute/gatk/tools/walkers/haplotypecaller/HaplotypeCallerGVCFIntegrationTest.java
View File

@ -126,7 +126,7 @@ public class HaplotypeCallerGVCFIntegrationTest extends WalkerTest {
final String WExIntervals = "-L 20:10,000,000-10,100,000 -isr INTERSECTION -L " + hg19Chr20Intervals;
// this functionality can be adapted to provide input data for whatever you might want in your data
tests.add(new Object[]{NA12878_PCRFREE, ReferenceConfidenceMode.NONE, PCRFreeIntervals, "0c87e26fdd7ab5629eb33f36833e3607"});
tests.add(new Object[]{NA12878_PCRFREE, ReferenceConfidenceMode.NONE, PCRFreeIntervals, "8bf132d73cf6b0851ae73c6799f19ba9"});
tests.add(new Object[]{NA12878_PCRFREE, ReferenceConfidenceMode.BP_RESOLUTION, PCRFreeIntervals, "90b25f3050435c9e67aa0ee325c24167"});
tests.add(new Object[]{NA12878_PCRFREE, ReferenceConfidenceMode.GVCF, PCRFreeIntervals, "5f329540dc5c4556ab029d0e2cfcabcb"});
tests.add(new Object[]{NA12878_WEx, ReferenceConfidenceMode.NONE, WExIntervals, "6ad7855dbf6dda2060aa93a3ee010b3e"});
@ -144,10 +144,10 @@ public class HaplotypeCallerGVCFIntegrationTest extends WalkerTest {
final String WExIntervals = "-L 20:10,000,000-10,100,000 -isr INTERSECTION -L " + hg19Chr20Intervals;
// this functionality can be adapted to provide input data for whatever you might want in your data
tests.add(new Object[]{NA12878_PCRFREE, ReferenceConfidenceMode.NONE, PCRFreeIntervals, "0820ae1d19ba0a2da25737ded8e2c96f"});
tests.add(new Object[]{NA12878_PCRFREE, ReferenceConfidenceMode.NONE, PCRFreeIntervals, "6662cfc41393257dfd6c39f1af1e3843"});
tests.add(new Object[]{NA12878_PCRFREE, ReferenceConfidenceMode.BP_RESOLUTION, PCRFreeIntervals, "70ee4e60d9f86b63aaab09075a71ddd3"});
tests.add(new Object[]{NA12878_PCRFREE, ReferenceConfidenceMode.GVCF, PCRFreeIntervals, "700d79df3b0b481444e81471204e242e"});
tests.add(new Object[]{NA12878_WEx, ReferenceConfidenceMode.NONE, WExIntervals, "09d1ae38586465b98dea0a0e432a7146"});
tests.add(new Object[]{NA12878_WEx, ReferenceConfidenceMode.NONE, WExIntervals, "af0fe243e3b96e59097187cd16ba1597"});
tests.add(new Object[]{NA12878_WEx, ReferenceConfidenceMode.BP_RESOLUTION, WExIntervals, "228e1d2ec2e729a5f79c37f3f2557708"});
tests.add(new Object[]{NA12878_WEx, ReferenceConfidenceMode.GVCF, WExIntervals, "2fc7020457dde4439b4133c098d9ab9b"});

10
protected/gatk-tools-protected/src/test/java/org/broadinstitute/gatk/tools/walkers/haplotypecaller/HaplotypeCallerIntegrationTest.java
View File

@ -490,5 +490,15 @@ public class HaplotypeCallerIntegrationTest extends WalkerTest {
public void testHBaseCountsBySample() throws IOException{
HCTest(NA12878_BAM, " -L 20:10001000-10010000 -A BaseCountsBySample", "f5ad4e03c0faaa806ee6ae536af8a479");
}
@Test
public void testSetZeroGQsToNoCall() throws IOException{
final File testBAM = new File(privateTestDir + "set.zero.GQs.no.call.bam");
final String md5 = "b90da12d97fce42f5127bcb6cad07b09";
final String base = String.format("-T HaplotypeCaller -R %s -I %s -L 8:17312375-17312975 ", REF, testBAM) +
" --no_cmdline_in_header -o %s";
final WalkerTestSpec spec = new WalkerTestSpec(base, Arrays.asList(md5));
executeTest("testSetZeroGQsToNoCall", spec);
}
}