Clone read base qualities rather than reference them #4926

DonFreed · 2018-06-20T22:45:18Z

When adjusting the base quality of overlapping read pairs, the modifications are made in place. If the modified reads are later used in another active region, the results from the later active region will be changed by the earlier modification.

This pull request clones the read's base quality, so that in place modifications do not affect later active regions.

codecov-io · 2018-06-20T23:44:21Z

Codecov Report

Merging #4926 into master will decrease coverage by 21.285%.
The diff coverage is 100%.

@@               Coverage Diff                @@
##              master     #4926        +/-   ##
================================================
- Coverage     80.435%   59.149%   -21.285%     
+ Complexity     17792     12568      -5224     
================================================
  Files           1090      1095         +5     
  Lines          64154     64606       +452     
  Branches       10342     10395        +53     
================================================
- Hits           51602     38214     -13388     
- Misses          8503     22145     +13642     
- Partials        4049      4247       +198

Impacted Files	Coverage Δ	Complexity Δ
...kers/haplotypecaller/AssemblyBasedCallerUtils.java	`63.366% <100%> (-6.931%)`	`22 <0> (-3)`
...broadinstitute/hellbender/utils/svd/SimpleSVD.java	`0% <0%> (-100%)`	`0% <0%> (-5%)`
...ogramgroups/ShortVariantDiscoveryProgramGroup.java	`0% <0%> (-100%)`	`0% <0%> (-3%)`
...adinstitute/hellbender/utils/R/RScriptLibrary.java	`0% <0%> (-100%)`	`0% <0%> (-6%)`
...ls/walkers/mutect/M2FiltersArgumentCollection.java	`0% <0%> (-100%)`	`0% <0%> (-1%)`
...odecs/gencode/GencodeGtfSelenocysteineFeature.java	`0% <0%> (-100%)`	`0% <0%> (-4%)`
...number/gcnv/GermlineCNVSegmentVariantComposer.java	`0% <0%> (-100%)`	`0% <0%> (-6%)`
...ctions/OptionalFeatureInputArgumentCollection.java	`0% <0%> (-100%)`	`0% <0%> (-1%)`
...ecaller/graphs/DeadEndKBestSubHaplotypeFinder.java	`0% <0%> (-100%)`	`0% <0%> (-10%)`
...otypecaller/RandomLikelihoodCalculationEngine.java	`0% <0%> (-100%)`	`0% <0%> (-6%)`
... and 614 more

magicDGS · 2018-06-21T05:45:54Z

src/main/java/org/broadinstitute/hellbender/utils/fragments/FragmentUtils.java

@@ -45,9 +45,9 @@ public static void adjustQualsOfOverlappingPairedFragments(final GATKRead clippe
        final int numOverlappingBases = Math.min(clippedFirstRead.getLength() - firstReadStop, clippedSecondRead.getLength());

        final byte[] firstReadBases = clippedFirstRead.getBases();
-        final byte[] firstReadQuals = clippedFirstRead.getBaseQualities();
+        final byte[] firstReadQuals = clippedFirstRead.getBaseQualities().clone();


GATKRead.getBases() and GATKRead.getBaseQualities() are safe methods and return a copy of the array. Thus, it is not really modified in place. I think that this patch is not needed for that reason.

lbergelson · 2018-06-21T19:40:25Z

@DonFreed, I agree with @magicDGS's assessment about it. This feels like a fix that was applied to Gatk3 but doesn't translate to 4? Of course, there could be implementations of GATKRead that don't obey the given contract about copying, but it's worth fixing those since we were more careful to think about copy/no copy when we wrote the new interface.

Of note: if you haven't seen it, GATKRead provides a set of unsafe getBaseQualitiesNoCopy() methods for times when the copy is a performance bottleneck and you can guarantee safe use of the underlying array.

I'm going to close this. Feel free to reopen if you disagree / can provide a unit test that demonstrates the issue still exists.

DonFreed · 2018-06-22T06:39:25Z

@lbergelson and @magicDGS, thank you for correcting me here. You are both right, this PR is based on an outdated understanding of the behavior of GATKRead.getBaseQualities(). Thanks for closing it!

lbergelson · 2018-06-22T14:49:46Z

@DonFreed Thank you for your continued interest and additions to GATK!

DonFreed · 2018-07-03T21:21:20Z

@lbergelson Could we reopen this? We were able to track down the source of the modified BQs across active regions.

At the end of the function adjustQualsOfOverlappingPairedFragments(final GATKRead clippedFirstRead, final GATKRead clippedSecondRead), the base qualities are set in the clipped reads:

clippedFirstRead.setBaseQualities(firstReadQuals);
clippedSecondRead.setBaseQualities(secondReadQuals);

In some cases, the clipped read is actually the original read, modifying the original read's BQ. I've committed a simple fix to my branch that ensures that the clipped read is never the original read.

droazen · 2018-07-06T17:51:54Z

@lbergelson Can you re-review?

droazen

Added my comments

droazen · 2018-07-06T17:55:42Z

...va/org/broadinstitute/hellbender/tools/walkers/haplotypecaller/AssemblyBasedCallerUtils.java

@@ -86,7 +86,7 @@ public static void finalizeRegion(final AssemblyRegion region,
            if ( ! clippedRead.isEmpty() && clippedRead.getCigar().getReadLength() > 0 ) {
                clippedRead = ReadClipper.hardClipToRegion( clippedRead, region.getExtendedSpan().getStart(), region.getExtendedSpan().getEnd() );
                if ( region.readOverlapsRegion(clippedRead) && clippedRead.getLength() > 0 ) {
-                    readsToUse.add(clippedRead);
+                    readsToUse.add((clippedRead == myRead) ? clippedRead.deepCopy() : clippedRead);


I wonder if we can get away with a less-expensive shallow copy (ie., copy()) rather than a costlier deepCopy() here?

Also, can you add a unit test that fails before the fix and passes after it?

copy() BQ instead of deepCopy().

DonFreed · 2018-07-11T17:25:07Z

@droazen, thanks for the comments. I've added a unit test to demonstrate the issue. The test shows that a shallow copy is fine, so that has been changed as well.

droazen · 2018-08-27T19:30:37Z

Re-assigning to @jamesemery

jamesemery

Just two minor comments on the test you wrote and we can push this in.

jamesemery · 2018-08-27T20:45:55Z

...roadinstitute/hellbender/tools/walkers/haplotypecaller/AssemblyBasedCallerUtilsUnitTest.java

+    }
+
+    @Test
+    public void testfinalizeRegion() {


This test looks good, though in this case since there is only one test in the AssemblyBasedCallerUtilsUnitTests I wouldn't use the @BeforeClass to initialize files. Instead just put them into static fields or build the header etc... into this test.

Also this is a good test

jamesemery · 2018-08-27T20:47:33Z

...roadinstitute/hellbender/tools/walkers/haplotypecaller/AssemblyBasedCallerUtilsUnitTest.java

+
+        SampleList sampleList = SampleList.singletonSampleList("tumor");
+        Byte minbq = 9;
+        AssemblyBasedCallerUtils.finalizeRegion(activeRegion, false, false, minbq, header, sampleList);


You should add a comment explaining what this test is doing here, (like explaining why it was failing due to unsafe operations on original reads and read qualities in the engine) so people know what is going on down the line.

DonFreed · 2018-09-04T20:29:41Z

Thanks @jamesemery! I've updated the PR to address your comments.

jamesemery · 2018-09-05T14:23:50Z

Thank you @DonFreed for your contribution!

…valid base qualities (#6943) This fixes the test that was accidentally disabled (when we added a flag for overlapping read adjustment we forgot to set it correctly for the test) and reinstated something like the code from #4926. It should now necessarily be the case that the finalize() method calls read.copy() at least once for every read (though many multiples of once are still quite possible). Fixes #6882

Clone read base qualities rather than reference them

d30ee96

magicDGS reviewed Jun 21, 2018

View reviewed changes

lbergelson closed this Jun 21, 2018

DonFreed deleted the df_cloneReadBQs branch June 22, 2018 06:39

DonFreed restored the df_cloneReadBQs branch July 3, 2018 17:50

Ensure clippedRead is a copy of myRead in finalizeRegion

4440749

lbergelson reopened this Jul 3, 2018

droazen requested a review from lbergelson July 6, 2018 17:52

droazen assigned lbergelson Jul 6, 2018

droazen suggested changes Jul 6, 2018

View reviewed changes

droazen assigned droazen and unassigned lbergelson Jul 6, 2018

Add a unit test for finalizeRegion.

6619e81

copy() BQ instead of deepCopy().

droazen assigned lbergelson and unassigned droazen Jul 12, 2018

droazen assigned jamesemery and unassigned lbergelson Aug 27, 2018

droazen requested a review from jamesemery August 27, 2018 19:13

jamesemery reviewed Aug 27, 2018

View reviewed changes

Remove @BeforeClass and add some explanatory comments

ef0094c

jamesemery merged commit 916d0ab into broadinstitute:master Sep 5, 2018

DonFreed deleted the df_cloneReadBQs branch September 5, 2018 18:39

DonFreed mentioned this pull request Oct 11, 2020

finalizeRegion modifies read base qualities in place #6882

Closed

2 tasks

jamesemery mentioned this pull request Nov 6, 2020

Fixed a bug where overlapping reads in subsequent regions can have invalid base qualities #6943

Merged

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

Clone read base qualities rather than reference them #4926

Clone read base qualities rather than reference them #4926

DonFreed commented Jun 20, 2018

codecov-io commented Jun 20, 2018 •

edited

Loading

magicDGS Jun 21, 2018

lbergelson commented Jun 21, 2018

DonFreed commented Jun 22, 2018 •

edited

Loading

lbergelson commented Jun 22, 2018

DonFreed commented Jul 3, 2018 •

edited

Loading

droazen commented Jul 6, 2018

droazen left a comment

droazen Jul 6, 2018

DonFreed commented Jul 11, 2018

droazen commented Aug 27, 2018

jamesemery left a comment

jamesemery Aug 27, 2018

jamesemery Aug 27, 2018

jamesemery Aug 27, 2018

DonFreed commented Sep 4, 2018

jamesemery commented Sep 5, 2018

Clone read base qualities rather than reference them #4926

Clone read base qualities rather than reference them #4926

Conversation

DonFreed commented Jun 20, 2018

codecov-io commented Jun 20, 2018 • edited Loading

Codecov Report

magicDGS Jun 21, 2018

Choose a reason for hiding this comment

lbergelson commented Jun 21, 2018

DonFreed commented Jun 22, 2018 • edited Loading

lbergelson commented Jun 22, 2018

DonFreed commented Jul 3, 2018 • edited Loading

droazen commented Jul 6, 2018

droazen left a comment

Choose a reason for hiding this comment

droazen Jul 6, 2018

Choose a reason for hiding this comment

DonFreed commented Jul 11, 2018

droazen commented Aug 27, 2018

jamesemery left a comment

Choose a reason for hiding this comment

jamesemery Aug 27, 2018

Choose a reason for hiding this comment

jamesemery Aug 27, 2018

Choose a reason for hiding this comment

jamesemery Aug 27, 2018

Choose a reason for hiding this comment

DonFreed commented Sep 4, 2018

jamesemery commented Sep 5, 2018

codecov-io commented Jun 20, 2018 •

edited

Loading

DonFreed commented Jun 22, 2018 •

edited

Loading

DonFreed commented Jul 3, 2018 •

edited

Loading