From 3fe011a11cea0c6ddf35f8d9efa46ffadbe74a82 Mon Sep 17 00:00:00 2001
From: Jonn Smith <jonn-smith@users.noreply.github.com>
Date: Fri, 8 Nov 2024 15:54:32 -0500
Subject: [PATCH] Port of `CallableLoci` from GATK3 (#9031)

* This is a port of a tool from GATK which is being used in some active pipelines.
* CallableLoci is a a tool which calculates coverage based stats and identifies regions of the genome which should be "callable" or not based on several filters.
---
 .../tools/walkers/coverage/CallableLoci.java  | 429 ++++++++++++++++++
 .../coverage/CallableLociIntegrationTest.java |  39 ++
 ...lable_loci.testBasicOperation.expected.bed | 422 +++++++++++++++++
 ...ci.testBasicOperation.expected.summary.txt |   7 +
 4 files changed, 897 insertions(+)
 create mode 100644 src/main/java/org/broadinstitute/hellbender/tools/walkers/coverage/CallableLoci.java
 create mode 100644 src/test/java/org/broadinstitute/hellbender/tools/walkers/coverage/CallableLociIntegrationTest.java
 create mode 100644 src/test/resources/callable_loci.testBasicOperation.expected.bed
 create mode 100644 src/test/resources/callable_loci.testBasicOperation.expected.summary.txt

diff --git a/src/main/java/org/broadinstitute/hellbender/tools/walkers/coverage/CallableLoci.java b/src/main/java/org/broadinstitute/hellbender/tools/walkers/coverage/CallableLoci.java
new file mode 100644
index 00000000000..8432db16606
--- /dev/null
+++ b/src/main/java/org/broadinstitute/hellbender/tools/walkers/coverage/CallableLoci.java
@@ -0,0 +1,429 @@
+package org.broadinstitute.hellbender.tools.walkers.coverage;
+
+
+import org.broadinstitute.barclay.argparser.Argument;
+import org.broadinstitute.barclay.argparser.CommandLineProgramProperties;
+import org.broadinstitute.barclay.argparser.Advanced;
+import org.broadinstitute.barclay.help.DocumentedFeature;
+import org.broadinstitute.barclay.argparser.ExperimentalFeature;
+import org.broadinstitute.hellbender.cmdline.StandardArgumentDefinitions;
+import org.broadinstitute.hellbender.cmdline.programgroups.CoverageAnalysisProgramGroup;
+import org.broadinstitute.hellbender.engine.FeatureContext;
+import org.broadinstitute.hellbender.engine.GATKPath;
+import org.broadinstitute.hellbender.engine.LocusWalker;
+import org.broadinstitute.hellbender.engine.ReferenceContext;
+import org.broadinstitute.hellbender.engine.filters.ReadFilter;
+import org.broadinstitute.hellbender.engine.filters.ReadFilterLibrary;
+import org.broadinstitute.hellbender.engine.filters.WellformedReadFilter;
+import org.broadinstitute.hellbender.engine.AlignmentContext;
+import org.broadinstitute.hellbender.exceptions.UserException;
+import org.broadinstitute.hellbender.utils.BaseUtils;
+import org.broadinstitute.hellbender.utils.SimpleInterval;
+import org.broadinstitute.hellbender.utils.pileup.PileupElement;
+import org.broadinstitute.hellbender.exceptions.GATKException;
+import htsjdk.samtools.util.Locatable;
+
+import java.util.ArrayList;
+import java.util.List;
+import java.io.IOException;
+import java.io.OutputStreamWriter;
+import java.util.stream.Collectors;
+import org.broadinstitute.hellbender.utils.Utils;
+import htsjdk.samtools.util.Lazy;
+import java.util.EnumMap;
+
+/**
+ * Collect statistics on callable, uncallable, poorly mapped, and other parts of the genome
+ *
+ * <p>
+ * A very common question about a NGS set of reads is what areas of the genome are considered callable. This tool
+ * considers the coverage at each locus and emits either a per base state or a summary interval BED file that
+ * partitions the genomic intervals into the following callable states:
+ * <dl>
+ * <dt>REF_N</dt>
+ * <dd>The reference base was an N, which is not considered callable the GATK</dd>
+ * <dt>PASS</dt>
+ * <dd>The base satisfied the min. depth for calling but had less than maxDepth to avoid having EXCESSIVE_COVERAGE</dd>
+ * <dt>NO_COVERAGE</dt>
+ * <dd>Absolutely no reads were seen at this locus, regardless of the filtering parameters</dd>
+ * <dt>LOW_COVERAGE</dt>
+ * <dd>There were fewer than min. depth bases at the locus, after applying filters</dd>
+ * <dt>EXCESSIVE_COVERAGE</dt>
+ * <dd>More than --max-depth read at the locus, indicating some sort of mapping problem</dd>
+ * <dt>POOR_MAPPING_QUALITY</dt>
+ * <dd>More than --max-fraction-of-reads-with-low-mapq at the locus, indicating a poor mapping quality of the reads</dd>
+ * </dl>
+ * </p>
+ * <p/>
+ * <h3>Input</h3>
+ * <p>
+ * A BAM file containing <b>exactly one sample</b>.
+ * </p>
+ * <p/>
+ * <h3>Output</h3>
+ * <p>
+ *     A file with the callable status covering each base and a table of callable status x count of all examined bases
+ * </p>
+ * <h3>Usage example</h3>
+ * <pre>
+ *  gatk CallableLoci \
+ *      -I myreads.bam \
+ *      -R myreference.fasta \
+ *      -O callable_status.bed \    
+ *      --summary table.txt 
+ * </pre>
+ * <p/>
+ * would produce a BED file that looks like:
+ * <pre>
+ *     20 10000000 10000864 PASS
+ *     20 10000865 10000985 POOR_MAPPING_QUALITY
+ *     20 10000986 10001138 PASS
+ *     20 10001139 10001254 POOR_MAPPING_QUALITY
+ *     20 10001255 10012255 PASS
+ *     20 10012256 10012259 POOR_MAPPING_QUALITY
+ *     20 10012260 10012263 PASS
+ *     20 10012264 10012328 POOR_MAPPING_QUALITY
+ *     20 10012329 10012550 PASS
+ *     20 10012551 10012551 LOW_COVERAGE
+ *     20 10012552 10012554 PASS
+ *     20 10012555 10012557 LOW_COVERAGE
+ *     20 10012558 10012558 PASS
+ * </pre>
+ * as well as a summary table that looks like:
+ * <p/>
+ * <pre>
+ *                        state nBases
+ *                        REF_N 0
+ *                         PASS 996046
+ *                  NO_COVERAGE 121
+ *                 LOW_COVERAGE 928
+ *           EXCESSIVE_COVERAGE 0
+ *         POOR_MAPPING_QUALITY 2906
+ * </pre>
+ *
+ * @author Mark DePristo / Jonn Smith
+ * @since May 7, 2010 / Nov 1, 2024
+ */
+@ExperimentalFeature
+@DocumentedFeature(groupName = "Coverage Analysis")
+@CommandLineProgramProperties(
+        summary = "Collect statistics on callable, uncallable, poorly mapped, and other parts of the genome",
+        oneLineSummary = "Determine callable status of loci",
+        programGroup = CoverageAnalysisProgramGroup.class
+)
+public class CallableLoci extends LocusWalker {
+
+    @Argument(fullName = StandardArgumentDefinitions.OUTPUT_LONG_NAME,
+            shortName = StandardArgumentDefinitions.OUTPUT_SHORT_NAME,
+            doc = "Output file (BED or per-base format)")
+    private GATKPath outputFile = null;
+
+    /**
+     * Callable loci summary counts will be written to this file.
+     */
+    @Argument(fullName = "summary", 
+            doc = "Name of file for output summary")
+    private GATKPath summaryFile;
+
+    /**
+     * The gap between this value and mmq are reads that are not sufficiently well mapped for calling but
+     * aren't indicative of mapping problems.  For example, if maxLowMAPQ = 1 and mmq = 20, then reads with
+     * MAPQ == 0 are poorly mapped, MAPQ >= 20 are considered as contributing to calling, where
+     * reads with MAPQ >= 1 and < 20 are not bad in and of themselves but aren't sufficiently good to contribute to
+     * calling.  In effect this reads are invisible, driving the base to the NO_ or LOW_COVERAGE states
+     */
+    @Argument(fullName = "max-low-mapq", shortName = "mlmq", minValue = 0, maxValue = 255, optional = true,
+            doc = "Maximum value for MAPQ to be considered a problematic mapped read")
+    private int maxLowMAPQ = 1;
+
+    /**
+     * Reads with MAPQ > minMappingQuality are treated as usable for variation detection, contributing to the PASS
+     * state.
+     */
+    @Argument(fullName = "min-mapping-quality", shortName = "mmq", minValue = 0, maxValue = 255, optional = true,
+            doc = "Minimum mapping quality of reads to count towards depth")
+    private int minMappingQuality = 10;
+
+    /**
+     * Bases with less than minBaseQuality are viewed as not sufficiently high quality to contribute to the PASS state
+     */
+    @Argument(fullName = "min-base-quality", shortName = "mbq", minValue = 0, maxValue = 255, optional = true,
+            doc = "Minimum quality of bases to count towards depth")
+    private int minBaseQuality = 20;
+
+    /**
+     * If the number of QC+ bases (on reads with MAPQ > minMappingQuality and with base quality > minBaseQuality) exceeds this
+     * value and is less than maxDepth the site is considered PASS.
+     */
+    @Advanced
+    @Argument(fullName = "min-depth", shortName = "min-depth", minValue = 0, optional = true,
+            doc = "Minimum QC+ read depth before a locus is considered callable")
+    private int minDepth = 4;
+
+    /**
+     * If the QC+ depth exceeds this value the site is considered to have EXCESSIVE_DEPTH
+     */
+    @Argument(fullName = "max-depth", shortName = "max-depth", optional = true,
+            doc = "Maximum read depth before a locus is considered poorly mapped")
+    private Integer maxDepth = null;
+
+    /**
+     * We don't want to consider a site as POOR_MAPPING_QUALITY just because it has two reads, and one is MAPQ.  We
+     * won't assign a site to the POOR_MAPPING_QUALITY state unless there are at least minDepthForLowMAPQ reads
+     * covering the site.
+     */
+    @Advanced
+    @Argument(fullName = "min-depth-for-low-mapq", shortName = "mdflmq", optional = true,
+            doc = "Minimum read depth before a locus is considered a potential candidate for poorly mapped")
+    private int minDepthLowMAPQ = 10;
+
+    /**
+     * If the number of reads at this site is greater than minDepthForLowMAPQ and the fraction of reads with low mapping quality
+     * exceeds this fraction then the site has POOR_MAPPING_QUALITY.
+     */
+    @Argument(fullName = "max-fraction-of-reads-with-low-mapq", shortName = "frlmq", optional = true,
+            doc = "If the fraction of reads at a base with low mapping quality exceeds this value, the site may be poorly mapped")
+    private double maxLowMAPQFraction = 0.1;
+
+    /**
+     * The output of this tool will be written in this format.  The recommended option is BED.
+     */
+    @Advanced
+    @Argument(fullName = "format", shortName = "format", optional = true,
+            doc = "Output format")
+    private OutputFormat outputFormat = OutputFormat.BED;
+
+    private OutputStreamWriter outputStream = null;
+    private OutputStreamWriter summaryStream = null;
+    
+    public enum OutputFormat {
+        BED,
+        STATE_PER_BASE
+    }
+
+    public enum State {
+        REF_N,
+        CALLABLE,
+        NO_COVERAGE,
+        LOW_COVERAGE,
+        EXCESSIVE_COVERAGE,
+        POOR_MAPPING_QUALITY
+    }
+
+    private BaseState currentState = null;
+    private final StateCounter stateCounter = new StateCounter();
+
+    protected static class BaseState implements Locatable {
+        private Locatable interval;
+        private final State state;
+
+        public BaseState(Locatable interval, State state) {
+            this.interval = interval;
+            this.state = state;
+        }
+
+        public State getState() {
+            return state;
+        }
+
+        @Override
+        public String getContig() {
+            return interval.getContig();
+        }
+
+        @Override
+        public int getStart() {
+            return interval.getStart();
+        }
+
+        @Override
+        public int getEnd() {
+            return interval.getEnd();
+        }
+
+        public void updateInterval(final Locatable newInterval) {
+            this.interval = new SimpleInterval(
+                interval.getContig(),
+                interval.getStart(),
+                newInterval.getEnd()
+            );
+        }
+
+        @Override
+        public String toString() {
+            return toBedString();
+        }
+
+        public String toBedString() {
+            // BED format is 0-based, so subtract 1 from start
+            return String.format("%s\t%d\t%d\t%s", 
+                interval.getContig(), 
+                interval.getStart() - 1, 
+                interval.getEnd(), 
+                state);
+        }
+    }
+
+    @Override
+    public boolean requiresReference() {
+        return true;
+    }
+
+    @Override
+    public boolean emitEmptyLoci() {
+        return true;
+    }
+
+    @Override
+    public boolean includeNs() {
+        return true;
+    }
+
+    @Override
+    // This is the default set of filters for CallableLoci as implemented in the GATK3 version of this tool.
+    public List<ReadFilter> getDefaultReadFilters() {
+        final List<ReadFilter> defaultFilters = new ArrayList<>(6);
+        defaultFilters.add(new ReadFilterLibrary.GoodCigarReadFilter());
+        defaultFilters.add(new ReadFilterLibrary.NotDuplicateReadFilter());
+        defaultFilters.add(new ReadFilterLibrary.PassesVendorQualityCheckReadFilter());
+        defaultFilters.add(new WellformedReadFilter());
+        defaultFilters.add(new ReadFilterLibrary.PrimaryLineReadFilter());
+        defaultFilters.add(new ReadFilterLibrary.MappedReadFilter());
+        return defaultFilters;
+    }
+
+    @Override
+    public void onTraversalStart() {
+        
+        // Validate sample count
+        final List<String> sampleList = getHeaderForReads().getReadGroups().stream()
+            .map(rg -> rg.getSample())
+            .distinct()
+            .collect(Collectors.toList());
+
+        if (sampleList.size() != 1) {
+            throw new UserException.BadInput("CallableLoci only works for a single sample.  Found " + sampleList.size() + " samples (" + String.join(", ", sampleList) + ").");
+        }
+
+        outputStream = new OutputStreamWriter(outputFile.getOutputStream());
+        summaryStream = new OutputStreamWriter(summaryFile.getOutputStream());
+    }
+
+    private State getCurrentState(ReferenceContext referenceContext, AlignmentContext alignmentContext) {
+        State state;
+
+        if (BaseUtils.isNBase(referenceContext.getBase())) {
+            state = State.REF_N;
+        } else {
+            int rawDepth = 0, QCDepth = 0, lowMAPQDepth = 0;
+            
+            for (PileupElement e : alignmentContext.getBasePileup()) {
+                rawDepth++;
+
+                if (e.getMappingQual() <= maxLowMAPQ) {
+                    lowMAPQDepth++;
+                }
+
+                if (e.getMappingQual() >= minMappingQuality && (e.getQual() >= minBaseQuality || e.isDeletion())) {
+                    QCDepth++;
+                }
+            }
+
+            if (rawDepth == 0) {
+                state = State.NO_COVERAGE;
+            } else if ((rawDepth >= minDepthLowMAPQ) && (((double)lowMAPQDepth) / ((double)rawDepth) >= maxLowMAPQFraction)) {
+                state = State.POOR_MAPPING_QUALITY;
+            } else if (QCDepth < minDepth) {
+                state = State.LOW_COVERAGE;
+            } else if (maxDepth != null && rawDepth >= maxDepth) {
+                state = State.EXCESSIVE_COVERAGE;
+            } else {
+                state = State.CALLABLE;
+            }
+        }
+
+        return state;
+    }
+
+    @Override
+    public void apply(AlignmentContext alignmentContext, ReferenceContext referenceContext, FeatureContext featureContext) {
+        
+        final State state = getCurrentState(referenceContext, alignmentContext);
+        BaseState callableState = new BaseState(alignmentContext, state);
+
+        // Update counts using the counter
+        stateCounter.increment(state);
+
+        try {
+            if (outputFormat == OutputFormat.STATE_PER_BASE) {
+                outputStream.write(callableState.toBedString() + "\n");
+            } else {
+                // BED format - integrate adjacent regions with same state
+                if (currentState == null) {
+                    currentState = callableState;
+                } else if (callableState.getStart() != currentState.getEnd() + 1 || currentState.getState() != callableState.getState()) {
+                    outputStream.write(currentState.toBedString() + "\n");
+                    currentState = callableState;
+                } else {
+                    currentState.updateInterval(callableState);
+                }
+            }
+        } catch (IOException e) {
+            throw new GATKException("Error writing to output stream", e);
+        }
+    }
+
+    @Override
+    public Object onTraversalSuccess() {
+
+        try{
+        // Print final state for BED format
+        if (outputFormat == OutputFormat.BED && currentState != null) {
+            outputStream.write(currentState.toBedString() + "\n");
+        }
+
+        // Write summary statistics using the counter
+        summaryStream.write(String.format("%30s %s%n", "state", "nBases"));
+        for (State state : State.values()) {
+            summaryStream.write(String.format("%30s %d%n", state, stateCounter.getCount(state)));
+        }
+
+        } catch (IOException e) {
+            throw new GATKException("Error writing to summary stream", e);
+        }
+        
+        return null;
+    }
+
+    @Override
+    public void closeTool() {
+        try {
+            outputStream.close();
+            summaryStream.close();
+        } catch (IOException e) {
+            throw new GATKException("Error closing output streams", e);
+        }
+    }
+
+    private static class StateCounter {
+        private final EnumMap<State, Long> counts = new EnumMap<>(State.class);
+
+        public StateCounter() {
+            for (State state : State.values()) {
+                counts.put(state, 0L);
+            }
+        }
+
+        public void increment(State state) {
+            counts.put(state, counts.get(state) + 1);
+        }
+
+        public long getCount(State state) {
+            return counts.get(state);
+        }
+
+        public EnumMap<State, Long> getCounts() {
+            return new EnumMap<>(counts);
+        }
+    }
+}
diff --git a/src/test/java/org/broadinstitute/hellbender/tools/walkers/coverage/CallableLociIntegrationTest.java b/src/test/java/org/broadinstitute/hellbender/tools/walkers/coverage/CallableLociIntegrationTest.java
new file mode 100644
index 00000000000..3346387d027
--- /dev/null
+++ b/src/test/java/org/broadinstitute/hellbender/tools/walkers/coverage/CallableLociIntegrationTest.java
@@ -0,0 +1,39 @@
+package org.broadinstitute.hellbender.tools.walkers.coverage;
+
+import org.broadinstitute.hellbender.CommandLineProgramTest;
+import org.broadinstitute.hellbender.testutils.ArgumentsBuilder;
+import org.broadinstitute.hellbender.testutils.IntegrationTestSpec;
+import org.testng.annotations.Test;
+
+import java.io.File;
+import java.util.Arrays;
+
+public class CallableLociIntegrationTest extends CommandLineProgramTest {
+
+    public static final String testDataDir = publicTestDir;
+
+    @Test
+    public void testBasicOperation() throws Exception {
+        final File outputFile = createTempFile("callableLoci", ".bed");
+        final File summaryFile = createTempFile("callableLoci", ".summary.txt");
+
+        final ArgumentsBuilder args = new ArgumentsBuilder();
+        args.addInput(new File(largeFileTestDir + "CEUTrio.HiSeq.WGS.b37.NA12878.20.21.bam"))
+            .addReference(new File(b37_reference_20_21))
+            .addOutput(outputFile)
+            .add("summary", summaryFile)
+            .add("format", "BED");
+
+        runCommandLine(args);
+
+        IntegrationTestSpec.assertEqualTextFiles(
+            outputFile,
+            new File(testDataDir + "callable_loci.testBasicOperation.expected.bed"),
+            "#");
+
+        IntegrationTestSpec.assertEqualTextFiles(
+            summaryFile,
+            new File(testDataDir + "callable_loci.testBasicOperation.expected.summary.txt"),
+            "#");
+    }
+}
\ No newline at end of file
diff --git a/src/test/resources/callable_loci.testBasicOperation.expected.bed b/src/test/resources/callable_loci.testBasicOperation.expected.bed
new file mode 100644
index 00000000000..055d22443c9
--- /dev/null
+++ b/src/test/resources/callable_loci.testBasicOperation.expected.bed
@@ -0,0 +1,422 @@
+20	0	60000	REF_N
+20	60000	9999901	NO_COVERAGE
+20	9999901	9999904	LOW_COVERAGE
+20	9999904	10000862	CALLABLE
+20	10000862	10000983	POOR_MAPPING_QUALITY
+20	10000983	10001000	CALLABLE
+20	10001000	10001003	POOR_MAPPING_QUALITY
+20	10001003	10001006	CALLABLE
+20	10001006	10001010	POOR_MAPPING_QUALITY
+20	10001010	10001074	CALLABLE
+20	10001074	10001254	POOR_MAPPING_QUALITY
+20	10001254	10015586	CALLABLE
+20	10015586	10015599	LOW_COVERAGE
+20	10015599	10098273	CALLABLE
+20	10098273	10098279	LOW_COVERAGE
+20	10098279	10098288	CALLABLE
+20	10098288	10098327	LOW_COVERAGE
+20	10098327	10098345	CALLABLE
+20	10098345	10098495	LOW_COVERAGE
+20	10098495	10098551	NO_COVERAGE
+20	10098551	10098685	LOW_COVERAGE
+20	10098685	10098689	CALLABLE
+20	10098689	10098690	LOW_COVERAGE
+20	10098690	10130790	CALLABLE
+20	10130790	10130886	POOR_MAPPING_QUALITY
+20	10130886	10199097	CALLABLE
+20	10199097	10199102	LOW_COVERAGE
+20	10199102	10199105	CALLABLE
+20	10199105	10199106	LOW_COVERAGE
+20	10199106	10199109	CALLABLE
+20	10199109	10199110	LOW_COVERAGE
+20	10199110	10199111	CALLABLE
+20	10199111	10199112	LOW_COVERAGE
+20	10199112	10199114	CALLABLE
+20	10199114	10199115	LOW_COVERAGE
+20	10199115	10199125	CALLABLE
+20	10199125	10199126	LOW_COVERAGE
+20	10199126	10199133	CALLABLE
+20	10199133	10199138	LOW_COVERAGE
+20	10199138	10199139	CALLABLE
+20	10199139	10199140	LOW_COVERAGE
+20	10199140	10199141	CALLABLE
+20	10199141	10199142	LOW_COVERAGE
+20	10199142	10199146	CALLABLE
+20	10199146	10199147	LOW_COVERAGE
+20	10199147	10199165	CALLABLE
+20	10199165	10199166	LOW_COVERAGE
+20	10199166	10199170	CALLABLE
+20	10199170	10199176	LOW_COVERAGE
+20	10199176	10199190	CALLABLE
+20	10199190	10199191	LOW_COVERAGE
+20	10199191	10199196	CALLABLE
+20	10199196	10199222	LOW_COVERAGE
+20	10199222	10199223	CALLABLE
+20	10199223	10199224	LOW_COVERAGE
+20	10199224	10199228	CALLABLE
+20	10199228	10199229	LOW_COVERAGE
+20	10199229	10199233	CALLABLE
+20	10199233	10199396	LOW_COVERAGE
+20	10199396	10250098	CALLABLE
+20	10250098	10250100	LOW_COVERAGE
+20	10250100	26319569	NO_COVERAGE
+20	26319569	29419569	REF_N
+20	29419569	29653908	NO_COVERAGE
+20	29653908	29803908	REF_N
+20	29803908	34897085	NO_COVERAGE
+20	34897085	34947085	REF_N
+20	34947085	61091437	NO_COVERAGE
+20	61091437	61141437	REF_N
+20	61141437	61213369	NO_COVERAGE
+20	61213369	61263369	REF_N
+20	61263369	62965520	NO_COVERAGE
+20	62965520	63025520	REF_N
+21	0	9411193	REF_N
+21	9411193	9595548	NO_COVERAGE
+21	9595548	9645548	REF_N
+21	9645548	9775437	NO_COVERAGE
+21	9775437	9825437	REF_N
+21	9825437	9999900	NO_COVERAGE
+21	9999900	9999906	LOW_COVERAGE
+21	9999906	10001414	POOR_MAPPING_QUALITY
+21	10001414	10001421	CALLABLE
+21	10001421	10004231	POOR_MAPPING_QUALITY
+21	10004231	10004255	CALLABLE
+21	10004255	10008552	POOR_MAPPING_QUALITY
+21	10008552	10008580	CALLABLE
+21	10008580	10009276	POOR_MAPPING_QUALITY
+21	10009276	10009375	CALLABLE
+21	10009375	10011967	POOR_MAPPING_QUALITY
+21	10011967	10012000	CALLABLE
+21	10012000	10012283	POOR_MAPPING_QUALITY
+21	10012283	10012313	CALLABLE
+21	10012313	10012505	POOR_MAPPING_QUALITY
+21	10012505	10012608	CALLABLE
+21	10012608	10012886	POOR_MAPPING_QUALITY
+21	10012886	10012910	CALLABLE
+21	10012910	10013107	POOR_MAPPING_QUALITY
+21	10013107	10013108	CALLABLE
+21	10013108	10013312	POOR_MAPPING_QUALITY
+21	10013312	10013362	CALLABLE
+21	10013362	10013598	POOR_MAPPING_QUALITY
+21	10013598	10013703	CALLABLE
+21	10013703	10013723	POOR_MAPPING_QUALITY
+21	10013723	10013724	CALLABLE
+21	10013724	10013725	POOR_MAPPING_QUALITY
+21	10013725	10013758	CALLABLE
+21	10013758	10014526	POOR_MAPPING_QUALITY
+21	10014526	10014557	CALLABLE
+21	10014557	10015174	POOR_MAPPING_QUALITY
+21	10015174	10015243	CALLABLE
+21	10015243	10015680	POOR_MAPPING_QUALITY
+21	10015680	10015682	CALLABLE
+21	10015682	10016888	POOR_MAPPING_QUALITY
+21	10016888	10016891	CALLABLE
+21	10016891	10018234	POOR_MAPPING_QUALITY
+21	10018234	10018282	CALLABLE
+21	10018282	10018667	POOR_MAPPING_QUALITY
+21	10018667	10018740	CALLABLE
+21	10018740	10019065	POOR_MAPPING_QUALITY
+21	10019065	10019144	CALLABLE
+21	10019144	10019342	POOR_MAPPING_QUALITY
+21	10019342	10019357	CALLABLE
+21	10019357	10020750	POOR_MAPPING_QUALITY
+21	10020750	10020759	CALLABLE
+21	10020759	10021272	POOR_MAPPING_QUALITY
+21	10021272	10021279	CALLABLE
+21	10021279	10021908	POOR_MAPPING_QUALITY
+21	10021908	10021915	CALLABLE
+21	10021915	10022538	POOR_MAPPING_QUALITY
+21	10022538	10022548	CALLABLE
+21	10022548	10023446	POOR_MAPPING_QUALITY
+21	10023446	10023456	CALLABLE
+21	10023456	10024149	POOR_MAPPING_QUALITY
+21	10024149	10024159	CALLABLE
+21	10024159	10026817	POOR_MAPPING_QUALITY
+21	10026817	10026831	CALLABLE
+21	10026831	10026952	POOR_MAPPING_QUALITY
+21	10026952	10026954	CALLABLE
+21	10026954	10027140	POOR_MAPPING_QUALITY
+21	10027140	10027142	CALLABLE
+21	10027142	10027311	POOR_MAPPING_QUALITY
+21	10027311	10027405	CALLABLE
+21	10027405	10027537	POOR_MAPPING_QUALITY
+21	10027537	10027545	CALLABLE
+21	10027545	10029107	POOR_MAPPING_QUALITY
+21	10029107	10029120	CALLABLE
+21	10029120	10029401	POOR_MAPPING_QUALITY
+21	10029401	10029428	CALLABLE
+21	10029428	10029429	POOR_MAPPING_QUALITY
+21	10029429	10029436	CALLABLE
+21	10029436	10029783	POOR_MAPPING_QUALITY
+21	10029783	10029786	CALLABLE
+21	10029786	10029787	POOR_MAPPING_QUALITY
+21	10029787	10029819	CALLABLE
+21	10029819	10030539	POOR_MAPPING_QUALITY
+21	10030539	10030590	CALLABLE
+21	10030590	10031708	POOR_MAPPING_QUALITY
+21	10031708	10031717	CALLABLE
+21	10031717	10034154	POOR_MAPPING_QUALITY
+21	10034154	10034167	CALLABLE
+21	10034167	10034359	POOR_MAPPING_QUALITY
+21	10034359	10034363	CALLABLE
+21	10034363	10034666	POOR_MAPPING_QUALITY
+21	10034666	10034674	CALLABLE
+21	10034674	10034915	POOR_MAPPING_QUALITY
+21	10034915	10034920	LOW_COVERAGE
+21	10034920	10084920	REF_N
+21	10084920	10085125	CALLABLE
+21	10085125	10085196	POOR_MAPPING_QUALITY
+21	10085196	10085482	CALLABLE
+21	10085482	10085486	POOR_MAPPING_QUALITY
+21	10085486	10085493	CALLABLE
+21	10085493	10085510	POOR_MAPPING_QUALITY
+21	10085510	10085511	CALLABLE
+21	10085511	10085548	POOR_MAPPING_QUALITY
+21	10085548	10086183	CALLABLE
+21	10086183	10086211	POOR_MAPPING_QUALITY
+21	10086211	10086230	CALLABLE
+21	10086230	10086233	POOR_MAPPING_QUALITY
+21	10086233	10086375	CALLABLE
+21	10086375	10086528	POOR_MAPPING_QUALITY
+21	10086528	10086630	CALLABLE
+21	10086630	10087011	POOR_MAPPING_QUALITY
+21	10087011	10087068	CALLABLE
+21	10087068	10088336	POOR_MAPPING_QUALITY
+21	10088336	10088340	LOW_COVERAGE
+21	10088340	10092857	POOR_MAPPING_QUALITY
+21	10092857	10092871	CALLABLE
+21	10092871	10094124	POOR_MAPPING_QUALITY
+21	10094124	10094126	LOW_COVERAGE
+21	10094126	10100485	POOR_MAPPING_QUALITY
+21	10100485	10100491	LOW_COVERAGE
+21	10100491	10100634	POOR_MAPPING_QUALITY
+21	10100634	10100645	LOW_COVERAGE
+21	10100645	10103150	POOR_MAPPING_QUALITY
+21	10103150	10103156	LOW_COVERAGE
+21	10103156	10103157	NO_COVERAGE
+21	10103157	10103160	LOW_COVERAGE
+21	10103160	10107501	POOR_MAPPING_QUALITY
+21	10107501	10107504	LOW_COVERAGE
+21	10107504	10108423	POOR_MAPPING_QUALITY
+21	10108423	10108425	LOW_COVERAGE
+21	10108425	10108428	POOR_MAPPING_QUALITY
+21	10108428	10108429	LOW_COVERAGE
+21	10108429	10110001	POOR_MAPPING_QUALITY
+21	10110001	10110007	LOW_COVERAGE
+21	10110007	10113883	POOR_MAPPING_QUALITY
+21	10113883	10113888	LOW_COVERAGE
+21	10113888	10113891	NO_COVERAGE
+21	10113891	10113909	LOW_COVERAGE
+21	10113909	10114200	POOR_MAPPING_QUALITY
+21	10114200	10114213	LOW_COVERAGE
+21	10114213	10118461	POOR_MAPPING_QUALITY
+21	10118461	10118466	LOW_COVERAGE
+21	10118466	10121719	POOR_MAPPING_QUALITY
+21	10121719	10121744	LOW_COVERAGE
+21	10121744	10122767	POOR_MAPPING_QUALITY
+21	10122767	10122774	LOW_COVERAGE
+21	10122774	10124504	POOR_MAPPING_QUALITY
+21	10124504	10124513	LOW_COVERAGE
+21	10124513	10124772	POOR_MAPPING_QUALITY
+21	10124772	10124774	CALLABLE
+21	10124774	10127843	POOR_MAPPING_QUALITY
+21	10127843	10127895	CALLABLE
+21	10127895	10129890	POOR_MAPPING_QUALITY
+21	10129890	10129926	CALLABLE
+21	10129926	10131172	POOR_MAPPING_QUALITY
+21	10131172	10131185	CALLABLE
+21	10131185	10132125	POOR_MAPPING_QUALITY
+21	10132125	10132147	LOW_COVERAGE
+21	10132147	10134066	POOR_MAPPING_QUALITY
+21	10134066	10134072	CALLABLE
+21	10134072	10136047	POOR_MAPPING_QUALITY
+21	10136047	10136097	CALLABLE
+21	10136097	10136560	POOR_MAPPING_QUALITY
+21	10136560	10136562	CALLABLE
+21	10136562	10136743	POOR_MAPPING_QUALITY
+21	10136743	10136765	CALLABLE
+21	10136765	10137054	POOR_MAPPING_QUALITY
+21	10137054	10137087	CALLABLE
+21	10137087	10137327	POOR_MAPPING_QUALITY
+21	10137327	10137357	CALLABLE
+21	10137357	10137582	POOR_MAPPING_QUALITY
+21	10137582	10137588	CALLABLE
+21	10137588	10137803	POOR_MAPPING_QUALITY
+21	10137803	10137808	CALLABLE
+21	10137808	10138250	POOR_MAPPING_QUALITY
+21	10138250	10138252	CALLABLE
+21	10138252	10138860	POOR_MAPPING_QUALITY
+21	10138860	10138928	CALLABLE
+21	10138928	10140888	POOR_MAPPING_QUALITY
+21	10140888	10140933	CALLABLE
+21	10140933	10141123	POOR_MAPPING_QUALITY
+21	10141123	10141250	CALLABLE
+21	10141250	10141280	POOR_MAPPING_QUALITY
+21	10141280	10141286	CALLABLE
+21	10141286	10141585	POOR_MAPPING_QUALITY
+21	10141585	10141587	CALLABLE
+21	10141587	10144390	POOR_MAPPING_QUALITY
+21	10144390	10144412	CALLABLE
+21	10144412	10145592	POOR_MAPPING_QUALITY
+21	10145592	10145594	LOW_COVERAGE
+21	10145594	10145595	NO_COVERAGE
+21	10145595	10145599	LOW_COVERAGE
+21	10145599	10147275	POOR_MAPPING_QUALITY
+21	10147275	10147277	CALLABLE
+21	10147277	10147285	POOR_MAPPING_QUALITY
+21	10147285	10147304	CALLABLE
+21	10147304	10147916	POOR_MAPPING_QUALITY
+21	10147916	10147917	LOW_COVERAGE
+21	10147917	10148437	POOR_MAPPING_QUALITY
+21	10148437	10148446	CALLABLE
+21	10148446	10150831	POOR_MAPPING_QUALITY
+21	10150831	10150840	CALLABLE
+21	10150840	10151788	POOR_MAPPING_QUALITY
+21	10151788	10151790	CALLABLE
+21	10151790	10153161	POOR_MAPPING_QUALITY
+21	10153161	10153175	CALLABLE
+21	10153175	10153853	POOR_MAPPING_QUALITY
+21	10153853	10153855	LOW_COVERAGE
+21	10153855	10155312	POOR_MAPPING_QUALITY
+21	10155312	10155323	LOW_COVERAGE
+21	10155323	10155683	POOR_MAPPING_QUALITY
+21	10155683	10155688	CALLABLE
+21	10155688	10156051	POOR_MAPPING_QUALITY
+21	10156051	10156055	CALLABLE
+21	10156055	10157124	POOR_MAPPING_QUALITY
+21	10157124	10157125	CALLABLE
+21	10157125	10158275	POOR_MAPPING_QUALITY
+21	10158275	10158281	LOW_COVERAGE
+21	10158281	10158282	NO_COVERAGE
+21	10158282	10158287	LOW_COVERAGE
+21	10158287	10160175	POOR_MAPPING_QUALITY
+21	10160175	10160178	LOW_COVERAGE
+21	10160178	10164945	POOR_MAPPING_QUALITY
+21	10164945	10164950	CALLABLE
+21	10164950	10165347	POOR_MAPPING_QUALITY
+21	10165347	10165361	CALLABLE
+21	10165361	10165647	POOR_MAPPING_QUALITY
+21	10165647	10165650	CALLABLE
+21	10165650	10166922	POOR_MAPPING_QUALITY
+21	10166922	10166927	CALLABLE
+21	10166927	10174185	POOR_MAPPING_QUALITY
+21	10174185	10174189	CALLABLE
+21	10174189	10174826	POOR_MAPPING_QUALITY
+21	10174826	10174830	CALLABLE
+21	10174830	10175880	POOR_MAPPING_QUALITY
+21	10175880	10175885	CALLABLE
+21	10175885	10176078	POOR_MAPPING_QUALITY
+21	10176078	10176079	LOW_COVERAGE
+21	10176079	10177872	POOR_MAPPING_QUALITY
+21	10177872	10177877	CALLABLE
+21	10177877	10178838	POOR_MAPPING_QUALITY
+21	10178838	10178843	LOW_COVERAGE
+21	10178843	10179300	POOR_MAPPING_QUALITY
+21	10179300	10179308	CALLABLE
+21	10179308	10179344	POOR_MAPPING_QUALITY
+21	10179344	10179354	CALLABLE
+21	10179354	10181387	POOR_MAPPING_QUALITY
+21	10181387	10181389	CALLABLE
+21	10181389	10182971	POOR_MAPPING_QUALITY
+21	10182971	10182975	CALLABLE
+21	10182975	10183371	POOR_MAPPING_QUALITY
+21	10183371	10183381	CALLABLE
+21	10183381	10184122	POOR_MAPPING_QUALITY
+21	10184122	10184130	CALLABLE
+21	10184130	10185818	POOR_MAPPING_QUALITY
+21	10185818	10185820	CALLABLE
+21	10185820	10186268	POOR_MAPPING_QUALITY
+21	10186268	10186286	CALLABLE
+21	10186286	10186356	POOR_MAPPING_QUALITY
+21	10186356	10186364	CALLABLE
+21	10186364	10186900	POOR_MAPPING_QUALITY
+21	10186900	10186909	CALLABLE
+21	10186909	10187286	POOR_MAPPING_QUALITY
+21	10187286	10187356	CALLABLE
+21	10187356	10187389	POOR_MAPPING_QUALITY
+21	10187389	10187398	CALLABLE
+21	10187398	10188005	POOR_MAPPING_QUALITY
+21	10188005	10188011	CALLABLE
+21	10188011	10189330	POOR_MAPPING_QUALITY
+21	10189330	10189342	CALLABLE
+21	10189342	10191035	POOR_MAPPING_QUALITY
+21	10191035	10191046	LOW_COVERAGE
+21	10191046	10191652	POOR_MAPPING_QUALITY
+21	10191652	10191654	CALLABLE
+21	10191654	10192714	POOR_MAPPING_QUALITY
+21	10192714	10192718	CALLABLE
+21	10192718	10192814	POOR_MAPPING_QUALITY
+21	10192814	10192815	CALLABLE
+21	10192815	10198849	POOR_MAPPING_QUALITY
+21	10198849	10198850	CALLABLE
+21	10198850	10199745	POOR_MAPPING_QUALITY
+21	10199745	10199748	LOW_COVERAGE
+21	10199748	10199815	POOR_MAPPING_QUALITY
+21	10199815	10199816	CALLABLE
+21	10199816	10199821	POOR_MAPPING_QUALITY
+21	10199821	10199826	CALLABLE
+21	10199826	10200721	POOR_MAPPING_QUALITY
+21	10200721	10200740	CALLABLE
+21	10200740	10200950	POOR_MAPPING_QUALITY
+21	10200950	10201039	CALLABLE
+21	10201039	10201248	POOR_MAPPING_QUALITY
+21	10201248	10201323	CALLABLE
+21	10201323	10201929	POOR_MAPPING_QUALITY
+21	10201929	10201935	CALLABLE
+21	10201935	10202412	POOR_MAPPING_QUALITY
+21	10202412	10202423	CALLABLE
+21	10202423	10203782	POOR_MAPPING_QUALITY
+21	10203782	10203795	CALLABLE
+21	10203795	10205580	POOR_MAPPING_QUALITY
+21	10205580	10205647	CALLABLE
+21	10205647	10208178	POOR_MAPPING_QUALITY
+21	10208178	10208238	CALLABLE
+21	10208238	10208952	POOR_MAPPING_QUALITY
+21	10208952	10208953	CALLABLE
+21	10208953	10208954	POOR_MAPPING_QUALITY
+21	10208954	10208956	CALLABLE
+21	10208956	10209828	POOR_MAPPING_QUALITY
+21	10209828	10209835	CALLABLE
+21	10209835	10209956	POOR_MAPPING_QUALITY
+21	10209956	10210047	CALLABLE
+21	10210047	10210393	POOR_MAPPING_QUALITY
+21	10210393	10210401	CALLABLE
+21	10210401	10212341	POOR_MAPPING_QUALITY
+21	10212341	10212342	CALLABLE
+21	10212342	10212438	POOR_MAPPING_QUALITY
+21	10212438	10212448	CALLABLE
+21	10212448	10212848	POOR_MAPPING_QUALITY
+21	10212848	10212850	CALLABLE
+21	10212850	10214810	POOR_MAPPING_QUALITY
+21	10214810	10214812	CALLABLE
+21	10214812	10215611	POOR_MAPPING_QUALITY
+21	10215611	10215615	CALLABLE
+21	10215615	10215950	POOR_MAPPING_QUALITY
+21	10215950	10215972	LOW_COVERAGE
+21	10215972	10215976	NO_COVERAGE
+21	10215976	10365976	REF_N
+21	10365976	10647896	NO_COVERAGE
+21	10647896	10697896	REF_N
+21	10697896	11188129	NO_COVERAGE
+21	11188129	14338129	REF_N
+21	14338129	33157035	NO_COVERAGE
+21	33157035	33157055	REF_N
+21	33157055	33157379	NO_COVERAGE
+21	33157379	33157389	REF_N
+21	33157389	40285944	NO_COVERAGE
+21	40285944	40285954	REF_N
+21	40285954	42955559	NO_COVERAGE
+21	42955559	43005559	REF_N
+21	43005559	43226828	NO_COVERAGE
+21	43226828	43227328	REF_N
+21	43227328	43249342	NO_COVERAGE
+21	43249342	43250842	REF_N
+21	43250842	44035894	NO_COVERAGE
+21	44035894	44035904	REF_N
+21	44035904	44632664	NO_COVERAGE
+21	44632664	44682664	REF_N
+21	44682664	44888040	NO_COVERAGE
+21	44888040	44888050	REF_N
+21	44888050	48119895	NO_COVERAGE
+21	48119895	48129895	REF_N
diff --git a/src/test/resources/callable_loci.testBasicOperation.expected.summary.txt b/src/test/resources/callable_loci.testBasicOperation.expected.summary.txt
new file mode 100644
index 00000000000..910c5d2771a
--- /dev/null
+++ b/src/test/resources/callable_loci.testBasicOperation.expected.summary.txt
@@ -0,0 +1,7 @@
+                         state nBases
+                         REF_N 16543253
+                      CALLABLE 253137
+                   NO_COVERAGE 94195953
+                  LOW_COVERAGE 800
+            EXCESSIVE_COVERAGE 0
+          POOR_MAPPING_QUALITY 162272