Fix Mutect2 IndexOutOfBoundException with germline resource (#7979)

broadinstitute · Aug 9, 2022 · e818963 · e818963
1 parent e1cc573
commit e818963
Show file tree

Hide file tree

Showing 3 changed files with 25 additions and 3 deletions.
diff --git a/src/main/java/org/broadinstitute/hellbender/tools/walkers/mutect/Mutect2Engine.java b/src/main/java/org/broadinstitute/hellbender/tools/walkers/mutect/Mutect2Engine.java
@@ -490,10 +490,9 @@ public ActivityProfileState isActive(final AlignmentContext context, final Refer
 
             for (final VariantContext germlineVC : germline) {
                 final List<Double> germlineAlleleFrequencies = getAttributeAsDoubleList(germlineVC, VCFConstants.ALLELE_FREQUENCY_KEY, 0.0);
-                final List<Allele> germlineAlts = germlineVC.getAlternateAlleles();
                 final Allele germlineRef = germlineVC.getReference();
 
-                for (int germlineAltIdx = 0; germlineAltIdx < germlineAlts.size(); germlineAltIdx++) {
+                for (int germlineAltIdx = 0; germlineAltIdx < germlineAlleleFrequencies.size(); germlineAltIdx++) {
                     if (germlineAlleleFrequencies.get(germlineAltIdx) < MTAC.maxPopulationAlleleFrequency) {
                         continue;
                     }

diff --git a/...main/java/org/broadinstitute/hellbender/tools/walkers/mutect/SomaticGenotypingEngine.java b/...main/java/org/broadinstitute/hellbender/tools/walkers/mutect/SomaticGenotypingEngine.java
@@ -12,6 +12,8 @@
 import org.apache.commons.math3.linear.DefaultRealMatrixChangingVisitor;
 import org.apache.commons.math3.linear.RealMatrix;
 import org.apache.commons.math3.util.FastMath;
+import org.apache.logging.log4j.LogManager;
+import org.apache.logging.log4j.Logger;
 import org.broadinstitute.hellbender.engine.FeatureContext;
 import org.broadinstitute.hellbender.engine.ReferenceContext;
 import org.broadinstitute.hellbender.tools.walkers.annotator.VariantAnnotatorEngine;
@@ -34,6 +36,7 @@
 import java.util.stream.IntStream;
 
 public class SomaticGenotypingEngine implements AutoCloseable {
+    protected static final Logger logger = LogManager.getLogger(SomaticGenotypingEngine.class);
 
     private final M2ArgumentCollection MTAC;
     private final Set<String> normalSamples;
@@ -375,6 +378,10 @@ private static Map<String, Object> getNegativeLogPopulationAFAnnotation(List<Var
     static double[] getGermlineAltAlleleFrequencies(final List<Allele> allAlleles, final Optional<VariantContext> germlineVC, final double afOfAllelesNotInGermlineResource) {
         Utils.validateArg(!allAlleles.isEmpty(), "allAlleles are empty -- there is not even a reference allele.");
         if (germlineVC.isPresent())  {
+            if (! germlineVC.get().hasAttribute(VCFConstants.ALLELE_FREQUENCY_KEY)) {
+                logger.warn("Germline resource variant at " + germlineVC.get().getContig() + ":" + germlineVC.get().getStart() +" missing AF attribute");
+                return Doubles.toArray(Collections.nCopies(allAlleles.size() - 1, afOfAllelesNotInGermlineResource));
+            }
             List<OptionalInt> germlineIndices = GATKVariantContextUtils.alleleIndices(allAlleles, germlineVC.get().getAlleles());
             final List<Double> germlineAltAFs = Mutect2Engine.getAttributeAsDoubleList(germlineVC.get(), VCFConstants.ALLELE_FREQUENCY_KEY, afOfAllelesNotInGermlineResource);
 

diff --git a/...a/org/broadinstitute/hellbender/tools/walkers/mutect/SomaticGenotypingEngineUnitTest.java b/...a/org/broadinstitute/hellbender/tools/walkers/mutect/SomaticGenotypingEngineUnitTest.java
@@ -42,7 +42,7 @@ Object[][] germlineAltAlleleFrequenciesData() {
 
                 //biallelic, same alt allele in different representations
                 { Arrays.asList(Allele.REF_A, Allele.ALT_G), Arrays.asList(ATref, GT), new double[] {0.1}, new double[] {0.1}},
-                { Arrays.asList(ATref, GT), Arrays.asList(Allele.REF_A, Allele.ALT_G), new double[] {0.1}, new double[] {0.1}},
+                { Arrays.asList(ATref, GT), Arrays.asList(Allele.REF_A, Allele.ALT_G), new double[] {0.1}, new double[] {0.1}}
         };
     }
 
@@ -56,4 +56,20 @@ public void testGetGermlineAltAlleleFrequencies(final List<Allele> calledAlleles
         final double[] result = SomaticGenotypingEngine.getGermlineAltAlleleFrequencies(calledAlleles, Optional.of(vc1), DEFAULT_AF);
         Assert.assertEquals(result, expected, 1.0e-10);
     }
+
+    @DataProvider(name = "missingAFData")
+    Object[][] missingAFData() {
+        return new Object[][]{
+                {new VariantContextBuilder("SOURCE", "1", 1, 1, Arrays.asList(Allele.REF_A, Allele.ALT_C))
+                        .attribute(VCFConstants.ALLELE_FREQUENCY_KEY, VCFConstants.MISSING_VALUE_v4).make()},
+                {new VariantContextBuilder("SOURCE", "1", 1, 1, Arrays.asList(Allele.REF_A, Allele.ALT_C))
+                        .make()}
+        };
+    }
+
+    @Test(dataProvider = "missingAFData")
+    public void testGetGermlineAltAlleleFrequenciesWithMissingAF(final VariantContext vc) {
+        final double[] result = SomaticGenotypingEngine.getGermlineAltAlleleFrequencies(vc.getAlleles(), Optional.of(vc), DEFAULT_AF);
+        Assert.assertEquals(result, new double[] {DEFAULT_AF}, 1.0e-10);
+    }
 }