Strand based annotations will use both reads in an overlapping read pair #5286

takutosato · 2018-10-05T19:29:13Z

This PR includes two changes:

Provide a command line argument to toggle the overlapping base quality correction (i.e. min(bq, 20)) before reassembly, which happens in FragmentUtils. I've found, however, that by the time SomaticGenotypingEngine runs, those the quality of these bases get bumped up to what they used to be, so this may be a no-op. I included it in case I missed something, and to be consistent with the branch @fleharty and @madduran have been using.
Provide a command line argument to count the two reads in an overlapping pair separately in StrandArtfiact and StrandBiasBySample. This feature is only available in Mutect i.e. it won't affect other tools that use StrandBiasBySample

takutosato · 2018-10-05T19:29:59Z

@davidbenjamin will you review?

fleharty · 2018-10-05T19:41:09Z

🍕 🎉 🎈 !

fleharty · 2018-10-09T14:32:15Z

@takutosato Tests appear to be failing, can you take a look at this?

takutosato · 2018-10-09T17:54:59Z

@fleharty fixed the bug, the tests should pass now

fleharty · 2018-10-09T21:16:30Z

@takutosato One of the builds failed, are you sure that fixed the bug?

davidbenjamin

Looks good! You can merge after addressing trivial comments.

davidbenjamin · 2018-10-10T03:32:52Z

src/main/java/org/broadinstitute/hellbender/tools/walkers/annotator/StrandArtifact.java

+        populateMapIfEmpty(refReads);
+
+        // assume that the mate of a forward read is a reverse read and vice versa
+        final int numDiscardedFwdAltMates = (int) altReads.get(Strand.REV).stream().filter(ba -> ba.read.hasAttribute(StrandBiasTest.DISCARDED_MATE_READ_TAG)).count();


Extract a method int countDiscardedMates for this logic even though it's a one-liner.

And when you extract such a method you could put in logic for when there are no reads with the given strand, thereby letting you do away withpopulateMapIfEmpty.

davidbenjamin · 2018-10-10T03:36:05Z

src/main/java/org/broadinstitute/hellbender/tools/walkers/annotator/StrandBiasTest.java

@@ -21,6 +21,7 @@

    protected static final int ARRAY_DIM = 2;
    protected static final int ARRAY_SIZE = ARRAY_DIM * ARRAY_DIM;
+    public static final String DISCARDED_MATE_READ_TAG = "DM";


This constant belongs in the class that applies the tag.

davidbenjamin · 2018-10-10T03:39:23Z

src/main/java/org/broadinstitute/hellbender/tools/walkers/annotator/StrandBiasTest.java

+            // This read's mate got discarded by SomaticGenotypingEngine::clipOverlappingReads()
+            if (read.hasAttribute(DISCARDED_MATE_READ_TAG)){
+                // Note that if this read is forward, then we increment its mate which we assume is reverse
+                table[offset + (isFW ? 1 : 0)]++;


isFW --> isForward

davidbenjamin · 2018-10-10T03:41:56Z

...nstitute/hellbender/tools/walkers/haplotypecaller/AssemblyBasedCallerArgumentCollection.java

@@ -123,4 +124,7 @@
    @Argument(fullName = SMITH_WATERMAN_LONG_NAME, doc = "Which Smith-Waterman implementation to use, generally FASTEST_AVAILABLE is the right choice", optional = true)
    public SmithWatermanAligner.Implementation smithWatermanImplementation = SmithWatermanAligner.Implementation.JAVA;

+    @Argument(fullName = CORRECT_OVERLAPPING_BASE_QUALITIES_LONG_NAME)
+    public boolean correctOverlappingBaseQualities = true;


The convention for boolean arguments is that specifying them yields non-default behavior.

davidbenjamin · 2018-10-10T03:44:41Z

src/main/java/org/broadinstitute/hellbender/tools/walkers/mutect/M2ArgumentCollection.java

@@ -169,4 +171,12 @@ public double getDefaultAlleleFrequency() {
    @Argument(fullName = GET_AF_FROM_AD_LONG_NAME, doc="Use allelic depth to calculate tumor allele fraction; recommended for mitochondrial applications", optional = true)
    public boolean calculateAFfromAD = false;

+    /**
+     * If set to true, count overlapping read pair as two seprate reads instead of one for {@link StrandArtifact} and {@link StrandBiasBySample}
+     * This would break the independence assumption of reads and over-count the alt depth in these annotations. On the other hand it could


This comment can be more direct about the fact that it is, in fact, the objectively correct behavior when it comes to strand artifacts.

davidbenjamin · 2018-10-10T03:45:39Z

src/main/java/org/broadinstitute/hellbender/tools/walkers/mutect/SomaticGenotypingEngine.java

@@ -101,7 +102,7 @@ public CalledHaplotypes callMutations(

        for( final int loc : startPosKeySet ) {
            final List<VariantContext> eventsAtThisLoc = getVariantContextsFromActiveHaplotypes(loc, haplotypes, false);
-            final VariantContext mergedVC = AssemblyBasedCallerUtils.makeMergedVariantContext(eventsAtThisLoc);
+            VariantContext mergedVC = AssemblyBasedCallerUtils.makeMergedVariantContext(eventsAtThisLoc);


Can't this remain final?

davidbenjamin · 2018-10-10T03:49:58Z

src/test/java/org/broadinstitute/hellbender/tools/walkers/mutect/Mutect2IntegrationTest.java

+        final byte altQuality = 50;
+        final byte refQuality = 30;
+
+        final File samFile = File.createTempFile("liquid-biopsy", ".bam");


This setting-up could go in a @BeforeClass method.

I tried to use @BeforeClass and ran into the issue where the @BeforeClass method takes in parameters (numAltPairs and refDepth) that I later use in Assert. Hardcoding them in the @BeforeClass method wouldn't be ideal. I could make them class-wide variables, since there are only two, but widening the scope of variables like that isn't ideal, either. So far I just extracted the sam file creation method without adding the @BeforeClass tag (in my latest commit). We can talk more tomorrow.

@takutosato That works for me, no need to shoehorn it into @BeforeClass.

…tions. 2. Add an option to turn off base quality correction before local reassembly

fleharty · 2018-10-11T18:58:27Z

Thanks Takuto! 🎈

…tions. 2. Add an option to turn off base quality correction before local reassembly (broadinstitute#5286)

takutosato requested a review from davidbenjamin October 5, 2018 19:29

davidbenjamin self-assigned this Oct 5, 2018

davidbenjamin requested changes Oct 10, 2018

View reviewed changes

davidbenjamin approved these changes Oct 11, 2018

View reviewed changes

davidbenjamin assigned takutosato and unassigned davidbenjamin Oct 11, 2018

takutosato force-pushed the ts_read_based_annotations branch from 0dd851a to b800fc4 Compare October 11, 2018 14:49

1. Add an option to count reads, not fragments, in strand bias annota…

d3d5f3e

…tions. 2. Add an option to turn off base quality correction before local reassembly

takutosato force-pushed the ts_read_based_annotations branch from b800fc4 to d3d5f3e Compare October 11, 2018 18:55

takutosato merged commit 4393c86 into master Oct 11, 2018

takutosato deleted the ts_read_based_annotations branch October 11, 2018 18:56

EdwardDixon pushed a commit to EdwardDixon/gatk that referenced this pull request Nov 9, 2018

1. Add an option to count reads, not fragments, in strand bias annota…

eebcb6e

…tions. 2. Add an option to turn off base quality correction before local reassembly (broadinstitute#5286)

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

Strand based annotations will use both reads in an overlapping read pair #5286

Strand based annotations will use both reads in an overlapping read pair #5286

takutosato commented Oct 5, 2018

takutosato commented Oct 5, 2018

fleharty commented Oct 5, 2018

fleharty commented Oct 9, 2018

takutosato commented Oct 9, 2018

fleharty commented Oct 9, 2018

davidbenjamin left a comment

davidbenjamin Oct 10, 2018

davidbenjamin Oct 10, 2018

takutosato Oct 10, 2018

davidbenjamin Oct 10, 2018

takutosato Oct 10, 2018

davidbenjamin Oct 10, 2018

takutosato Oct 10, 2018

davidbenjamin Oct 10, 2018

davidbenjamin Oct 10, 2018

takutosato Oct 10, 2018

takutosato Oct 10, 2018

davidbenjamin Oct 10, 2018

takutosato Oct 10, 2018

takutosato Oct 10, 2018

davidbenjamin Oct 10, 2018

takutosato Oct 10, 2018

davidbenjamin Oct 11, 2018

fleharty commented Oct 11, 2018

Strand based annotations will use both reads in an overlapping read pair #5286

Strand based annotations will use both reads in an overlapping read pair #5286

Conversation

takutosato commented Oct 5, 2018

takutosato commented Oct 5, 2018

fleharty commented Oct 5, 2018

fleharty commented Oct 9, 2018

takutosato commented Oct 9, 2018

fleharty commented Oct 9, 2018

davidbenjamin left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

fleharty commented Oct 11, 2018