diff --git a/java/src/org/broadinstitute/sting/gatk/walkers/genotyper/UnifiedGenotyperEngine.java b/java/src/org/broadinstitute/sting/gatk/walkers/genotyper/UnifiedGenotyperEngine.java index 4d2326aeb..37c09dac7 100755 --- a/java/src/org/broadinstitute/sting/gatk/walkers/genotyper/UnifiedGenotyperEngine.java +++ b/java/src/org/broadinstitute/sting/gatk/walkers/genotyper/UnifiedGenotyperEngine.java @@ -158,8 +158,11 @@ public class UnifiedGenotyperEngine { */ public VariantCallContext calculateLikelihoodsAndGenotypes(RefMetaDataTracker tracker, ReferenceContext refContext, AlignmentContext rawContext) { Map stratifiedContexts = getFilteredAndStratifiedContexts(UAC, refContext, rawContext); + if ( stratifiedContexts == null ) + return (UAC.OutputMode != OUTPUT_MODE.EMIT_ALL_SITES ? null : new VariantCallContext(generateEmptyContext(tracker, refContext, stratifiedContexts, rawContext), refContext.getBase(), false)); + VariantContext vc = calculateLikelihoods(tracker, refContext, stratifiedContexts, StratifiedAlignmentContext.StratifiedContextType.COMPLETE, null); - if ( vc == null || !vc.hasGenotypes() ) + if ( vc == null ) return null; VariantCallContext vcc = calculateGenotypes(tracker, refContext, rawContext, stratifiedContexts, vc); @@ -178,13 +181,13 @@ public class UnifiedGenotyperEngine { */ public VariantContext calculateLikelihoods(RefMetaDataTracker tracker, ReferenceContext refContext, AlignmentContext rawContext, Allele alternateAlleleToUse) { Map stratifiedContexts = getFilteredAndStratifiedContexts(UAC, refContext, rawContext); + if ( stratifiedContexts == null ) + return null; VariantContext vc = calculateLikelihoods(tracker, refContext, stratifiedContexts, StratifiedAlignmentContext.StratifiedContextType.COMPLETE, alternateAlleleToUse); return GLsToPLs(vc); } private VariantContext calculateLikelihoods(RefMetaDataTracker tracker, ReferenceContext refContext, Map stratifiedContexts, StratifiedAlignmentContext.StratifiedContextType type, Allele alternateAlleleToUse) { - if ( stratifiedContexts == null ) - return null; // initialize the data for this thread if that hasn't been done yet if ( glcm.get() == null ) { @@ -201,6 +204,32 @@ public class UnifiedGenotyperEngine { return null; } + private VariantContext generateEmptyContext(RefMetaDataTracker tracker, ReferenceContext ref, Map stratifiedContexts, AlignmentContext rawContext) { + VariantContext vc; + if ( UAC.GenotypingMode == GenotypeLikelihoodsCalculationModel.GENOTYPING_MODE.GENOTYPE_GIVEN_ALLELES ) { + final VariantContext vcInput = tracker.getVariantContext(ref, "alleles", null, ref.getLocus(), true); + vc = new VariantContext("UG_call", vcInput.getChr(), vcInput.getStart(), vcInput.getEnd(), vcInput.getAlleles()); + } else { + Set alleles = new HashSet(); + alleles.add(Allele.create(ref.getBase(), true)); + vc = new VariantContext("UG_call", ref.getLocus().getContig(), ref.getLocus().getStart(), ref.getLocus().getStart(), alleles); + } + + if ( annotationEngine != null ) { + // we want to use the *unfiltered* context for the annotations + ReadBackedPileup pileup = null; + if (rawContext.hasExtendedEventPileup()) + pileup = rawContext.getExtendedEventPileup(); + else if (rawContext.hasBasePileup()) + pileup = rawContext.getBasePileup(); + stratifiedContexts = StratifiedAlignmentContext.splitContextBySampleName(pileup, UAC.ASSUME_SINGLE_SAMPLE); + + vc = annotationEngine.annotateContext(tracker, ref, stratifiedContexts, vc).iterator().next(); + } + + return vc; + } + private VariantContext createVariantContextFromLikelihoods(ReferenceContext refContext, Allele refAllele, Map GLs) { // no-call everyone for now List noCall = new ArrayList(); @@ -278,7 +307,9 @@ public class UnifiedGenotyperEngine { // estimate our confidence in a reference call and return if ( vc.getNSamples() == 0 ) - return estimateReferenceConfidence(stratifiedContexts, genotypePriors.getHeterozygosity(), false, 1.0); + return (UAC.OutputMode != OUTPUT_MODE.EMIT_ALL_SITES ? + estimateReferenceConfidence(stratifiedContexts, genotypePriors.getHeterozygosity(), false, 1.0) : + new VariantCallContext(generateEmptyContext(tracker, refContext, stratifiedContexts, rawContext), refContext.getBase(), false)); // 'zero' out the AFs (so that we don't have to worry if not all samples have reads at this position) clearAFarray(log10AlleleFrequencyPosteriors.get()); diff --git a/java/test/org/broadinstitute/sting/gatk/walkers/genotyper/UnifiedGenotyperIntegrationTest.java b/java/test/org/broadinstitute/sting/gatk/walkers/genotyper/UnifiedGenotyperIntegrationTest.java index 459851a19..11025c186 100755 --- a/java/test/org/broadinstitute/sting/gatk/walkers/genotyper/UnifiedGenotyperIntegrationTest.java +++ b/java/test/org/broadinstitute/sting/gatk/walkers/genotyper/UnifiedGenotyperIntegrationTest.java @@ -144,7 +144,7 @@ public class UnifiedGenotyperIntegrationTest extends WalkerTest { HashMap e = new HashMap(); e.put( "-sites_only", "71e561ba6fc66bd8b84907252f71ea55" ); e.put( "--output_mode EMIT_ALL_CONFIDENT_SITES", "4ffcb1e1f20ce175783c32c30deef8db" ); - e.put( "--output_mode EMIT_ALL_SITES", "3d98205a31a133c11e518e095dc7ab65" ); + e.put( "--output_mode EMIT_ALL_SITES", "2ca25ad91f7a746f715afdca5d516768" ); for ( Map.Entry entry : e.entrySet() ) { WalkerTest.WalkerTestSpec spec = new WalkerTest.WalkerTestSpec(