Skip to content

Commit

Permalink
reverting back
Browse files Browse the repository at this point in the history
  • Loading branch information
@SHuang-Broad
SHuang-Broad committed Jul 23, 2018
1 parent 8c3b4b8 commit 1266952
Show file tree
Hide file tree
Showing 2 changed files with 27 additions and 23 deletions.
35 changes: 20 additions & 15 deletions ...e/hellbender/tools/spark/sv/discovery/alignment/AssemblyContigAlignmentsConfigPicker.java
View file
Original file line number Diff line number Diff line change
Expand Up @@ -10,7 +10,6 @@
import org.apache.spark.api.java.JavaRDD;
import org.broadinstitute.hellbender.exceptions.GATKException;
import org.broadinstitute.hellbender.tools.spark.sv.discovery.SvDiscoverFromLocalAssemblyContigAlignmentsSpark;
import org.broadinstitute.hellbender.tools.spark.sv.utils.SVInterval;
import org.broadinstitute.hellbender.tools.spark.sv.utils.SVUtils;
import org.broadinstitute.hellbender.tools.spark.sv.utils.SvCigarUtils;
import org.broadinstitute.hellbender.utils.Utils;
Expand Down Expand Up @@ -823,13 +822,24 @@ private static void removeDueToShortReadSpan(final List<AlignmentInterval> selec
// subtract the overlap from the distance covered on the contig by the alignment.
// This gives unique read region it explains.
// If this unique read region is "short": shorter than {@code uniqReadLenInclusive}), drop it.
final Map<AlignmentInterval, SVInterval> maxOverlapMap = getMaxOverlapPairs(selectedAlignments);

// each alignment has an entry of a tuple2, one for max overlap maxFront, one for max overlap maxRear,
// max overlap maxFront is a tuple2 registering the index and overlap bases count
final Map<AlignmentInterval, Tuple2<Integer, Integer>> maxOverlapMap = getMaxOverlapPairs(selectedAlignments);
for(Iterator<AlignmentInterval> iterator = selectedAlignments.iterator(); iterator.hasNext();) {
final AlignmentInterval alignment = iterator.next();

final SVInterval uniqReadSpan = maxOverlapMap.get(alignment);

if (uniqReadSpan == null || uniqReadSpan.getLength() < uniqReadLenInclusive) {
final Tuple2<Integer, Integer> maxOverlapFrontAndRear = maxOverlapMap.get(alignment);
final int maxOverlapFront = Math.max(0, maxOverlapFrontAndRear._1);
final int maxOverlapRear = Math.max(0, maxOverlapFrontAndRear._2);

// theoretically this could be negative for an alignment whose maxFront and maxRear sums together bigger than the read span
// but earlier configuration scoring would make this impossible because such alignments should be filtered out already
// considering that it brings more penalty than value, i.e. read bases explained (even if the MQ is 60),
// but even if it is kept, a negative value won't hurt unless a stupid threshold value is passed in
final int uniqReadSpan = alignment.endInAssembledContig - alignment.startInAssembledContig + 1
- maxOverlapFront - maxOverlapRear;
if (uniqReadSpan < uniqReadLenInclusive) {
lowUniquenessMappings.add(alignment);
iterator.remove();
}
Expand Down Expand Up @@ -857,14 +867,12 @@ private static final class TempMaxOverlapInfo {
}
}

@VisibleForTesting
static final int TEMPORARY_SELF_POINTING_SV_INTERVAL_CONTIG = 0;
/**
* Extract the max overlap information, front and back, for each alignment in {@code configuration}.
* For each alignment, the corresponding value is the unique read span by taking away the maximum overlaps.
* For each alignment, the corresponding tuple2 has the max (front, rear) overlap base counts.
*/
@VisibleForTesting
static Map<AlignmentInterval, SVInterval> getMaxOverlapPairs(final List<AlignmentInterval> configuration) {
static Map<AlignmentInterval, Tuple2<Integer, Integer>> getMaxOverlapPairs(final List<AlignmentInterval> configuration) {

final List<TempMaxOverlapInfo> intermediateResult =
new ArrayList<>(Collections.nCopies(configuration.size(), new TempMaxOverlapInfo()));
Expand Down Expand Up @@ -903,13 +911,10 @@ static Map<AlignmentInterval, SVInterval> getMaxOverlapPairs(final List<Alignmen
}
}

final Map<AlignmentInterval, SVInterval> maxOverlapMap = new HashMap<>(configuration.size());
final Map<AlignmentInterval, Tuple2<Integer, Integer>> maxOverlapMap = new HashMap<>(configuration.size());
for (int i = 0; i < configuration.size(); ++i) {
final AlignmentInterval alignment = configuration.get(i);
int uniqueStart = alignment.startInAssembledContig + Math.max(0, intermediateResult.get(i).maxFront._2) - 1;// "-1" because we are going from 1-based to 0-based
int uniqueEnd = alignment.endInAssembledContig - Math.max(0, intermediateResult.get(i).maxRear._2);
final SVInterval uniqReadSpan = uniqueStart >= uniqueEnd ? null : new SVInterval(TEMPORARY_SELF_POINTING_SV_INTERVAL_CONTIG, uniqueStart, uniqueEnd);
maxOverlapMap.put(configuration.get(i), uniqReadSpan);
maxOverlapMap.put(configuration.get(i),
new Tuple2<>(intermediateResult.get(i).maxFront._2, intermediateResult.get(i).maxRear._2));
}

return maxOverlapMap;
Expand Down
15 changes: 7 additions & 8 deletions ...nder/tools/spark/sv/discovery/alignment/AssemblyContigAlignmentsConfigPickerUnitTest.java
View file
Original file line number Diff line number Diff line change
Expand Up @@ -5,7 +5,6 @@
import org.broadinstitute.hellbender.GATKBaseTest;
import org.broadinstitute.hellbender.engine.spark.SparkContextFactory;
import org.broadinstitute.hellbender.tools.spark.sv.discovery.TestUtilsForAssemblyBasedSVDiscovery;
import org.broadinstitute.hellbender.tools.spark.sv.utils.SVInterval;
import org.broadinstitute.hellbender.utils.SimpleInterval;
import org.broadinstitute.hellbender.utils.read.CigarUtils;
import org.testng.Assert;
Expand Down Expand Up @@ -550,18 +549,18 @@ public void testSortingAfterGapSplit(final GoodAndBadMappings input, final boole
new AlignmentInterval(new SimpleInterval("chr18:77373208-77373391"), 1448, 1631, TextCigarCodec.decode("1447S184M980S"), true, 60, AlignmentInterval.NO_NM, AlignmentInterval.NO_AS, ContigAlignmentsModifier.AlnModType.FROM_SPLIT_GAPPED_ALIGNMENT),
new AlignmentInterval(new SimpleInterval("chr18:77373514-77374571"), 1554, 2611, TextCigarCodec.decode("1165H388S1058M"), true, 60, AlignmentInterval.NO_NM, AlignmentInterval.NO_AS, ContigAlignmentsModifier.AlnModType.FROM_SPLIT_GAPPED_ALIGNMENT));

final Map<AlignmentInterval, SVInterval> overlapMap = new HashMap<>(4);
overlapMap.put(alignments.get(0), new SVInterval(AssemblyContigAlignmentsConfigPicker.TEMPORARY_SELF_POINTING_SV_INTERVAL_CONTIG, 0, 1165));
overlapMap.put(alignments.get(1), null);
overlapMap.put(alignments.get(2), new SVInterval(AssemblyContigAlignmentsConfigPicker.TEMPORARY_SELF_POINTING_SV_INTERVAL_CONTIG, 1503, 1553));
overlapMap.put(alignments.get(3), new SVInterval(AssemblyContigAlignmentsConfigPicker.TEMPORARY_SELF_POINTING_SV_INTERVAL_CONTIG, 1631, 2611));
final Map<AlignmentInterval, Tuple2<Integer, Integer>> overlapMap = new HashMap<>(4);
overlapMap.put(alignments.get(0), new Tuple2<>(-1, 282));
overlapMap.put(alignments.get(1), new Tuple2<>(282, 56));
overlapMap.put(alignments.get(2), new Tuple2<>(56, 78));
overlapMap.put(alignments.get(3), new Tuple2<>(78, -1));
data.add(new Object[]{alignments, overlapMap});

return data.toArray(new Object[data.size()][]);
}
@Test(groups = "sv", dataProvider = "forGetMaxOverlapPairs")
public void testGetMaxOverlapPairs(final List<AlignmentInterval> alignments, final Map<AlignmentInterval, SVInterval> expected) {
final Map<AlignmentInterval, SVInterval> actual = AssemblyContigAlignmentsConfigPicker.getMaxOverlapPairs(alignments);
public void testGetMaxOverlapPairs(final List<AlignmentInterval> alignments, final Map<AlignmentInterval, Tuple2<Integer, Integer>> expected) {
final Map<AlignmentInterval, Tuple2<Integer, Integer>> actual = AssemblyContigAlignmentsConfigPicker.getMaxOverlapPairs(alignments);
twoSetsEqualIgnoreOrder(new HashSet<>(actual.keySet()), new HashSet<>(expected.keySet()));
actual.forEach((k,v) -> Assert.assertEquals(v, expected.get(k)));
}
Expand Down

0 comments on commit 1266952

Please sign in to comment.