Skip to content

Commit

Permalink
#156 now passing test cases
Browse files Browse the repository at this point in the history
  • Loading branch information
@d-cameron
d-cameron committed Jul 31, 2018
1 parent 3e96595 commit 0c4d318
Show file tree
Hide file tree
Showing 8 changed files with 100 additions and 53 deletions.
12 changes: 7 additions & 5 deletions src/main/java/au/edu/wehi/idsv/AssemblyAttributes.java
View file
Original file line number Diff line number Diff line change
Expand Up @@ -48,12 +48,14 @@ public AssemblyAttributes(SingleReadEvidence record) {
public static void adjustAssemblyAnnotationDueToContigChange(SAMRecord record, int startTruncatedBases) {
int[] intervalStart = record.getSignedIntArrayAttribute(SamTags.ASSEMBLY_EVIDENCE_OFFSET_START);
int[] intervalEnd = record.getSignedIntArrayAttribute(SamTags.ASSEMBLY_EVIDENCE_OFFSET_END);
for (int i = 0; i < intervalStart.length; i++) {
intervalStart[i] -= startTruncatedBases;
intervalEnd[i] -= startTruncatedBases;
if (intervalStart != null && intervalEnd != null) {
for (int i = 0; i < intervalStart.length; i++) {
intervalStart[i] -= startTruncatedBases;
intervalEnd[i] -= startTruncatedBases;
}
record.setAttribute(SamTags.ASSEMBLY_EVIDENCE_OFFSET_START, intervalStart);
record.setAttribute(SamTags.ASSEMBLY_EVIDENCE_OFFSET_END, intervalEnd);
}
record.setAttribute(SamTags.ASSEMBLY_EVIDENCE_OFFSET_START, intervalStart);
record.setAttribute(SamTags.ASSEMBLY_EVIDENCE_OFFSET_END, intervalEnd);
}

/**
Expand Down
38 changes: 21 additions & 17 deletions src/main/java/au/edu/wehi/idsv/SingleReadEvidence.java
View file
Original file line number Diff line number Diff line change
Expand Up @@ -412,23 +412,27 @@ public Range<Integer> getBreakendAssemblyContigBreakpointInterval() {
return Range.closed(rl - anchorBases, rl - anchorBases);
}
}
// Variant has either inserted sequence or homology
Range<Integer> r;
int insertLength = getUntemplatedSequence().length();
if (insertLength > 0) {
if ((location.direction == BreakendDirection.Forward && !record.getReadNegativeStrandFlag()) |
(location.direction == BreakendDirection.Backward && record.getReadNegativeStrandFlag())) {
r = Range.closed(nominalBreakendAfterReadOffset, nominalBreakendAfterReadOffset + insertLength);
} else {
r = Range.closed(nominalBreakendAfterReadOffset - insertLength, nominalBreakendAfterReadOffset);
}
} else {
if (!record.getReadNegativeStrandFlag()) {
r = Range.closed(nominalBreakendAfterReadOffset + (location.start - location.nominal), nominalBreakendAfterReadOffset + (location.end - location.nominal));
} else {
r = Range.closed(nominalBreakendAfterReadOffset - (location.end - location.nominal), nominalBreakendAfterReadOffset - (location.start - location.nominal));
}
}
Range<Integer> r;
if (location instanceof BreakpointSummary) {
// Variant has either inserted sequence or homology
int insertLength = getUntemplatedSequence().length();
if (insertLength > 0) {
if ((location.direction == BreakendDirection.Forward && !record.getReadNegativeStrandFlag()) |
(location.direction == BreakendDirection.Backward && record.getReadNegativeStrandFlag())) {
r = Range.closed(nominalBreakendAfterReadOffset, nominalBreakendAfterReadOffset + insertLength);
} else {
r = Range.closed(nominalBreakendAfterReadOffset - insertLength, nominalBreakendAfterReadOffset);
}
} else {
if (!record.getReadNegativeStrandFlag()) {
r = Range.closed(nominalBreakendAfterReadOffset + (location.start - location.nominal), nominalBreakendAfterReadOffset + (location.end - location.nominal));
} else {
r = Range.closed(nominalBreakendAfterReadOffset - (location.end - location.nominal), nominalBreakendAfterReadOffset - (location.start - location.nominal));
}
}
} else {
r = Range.closed(nominalBreakendAfterReadOffset, nominalBreakendAfterReadOffset);
}
return r;
}
public int getBreakendAssemblyContigOffset() {
Expand Down
6 changes: 4 additions & 2 deletions src/main/java/au/edu/wehi/idsv/StructuralVariationCallBuilder.java
View file
Original file line number Diff line number Diff line change
Expand Up @@ -424,8 +424,10 @@ private <T extends SingleReadEvidence> double[] prorataAssemblyQualBreakdown(
breakdownQual[i] = aa.getSupportingQualScore(offset, ImmutableSet.of(i), null);
}
double breakdownTotal = DoubleStream.of(breakdownQual).sum();
for (int i = 0; i < prorata.length; i++) {
prorata[i] += assQual * (breakdownQual[i] / breakdownTotal);
if (breakdownTotal != 0) { // defensive check to mitigate impact of 0 qual assemblies (#156)
for (int i = 0; i < prorata.length; i++) {
prorata[i] += assQual * (breakdownQual[i] / breakdownTotal);
}
}
totalAssQual += assQual;
}
Expand Down
2 changes: 1 addition & 1 deletion src/test/java/au/edu/wehi/idsv/AssemblyFactoryTest.java
View file
Original file line number Diff line number Diff line change
Expand Up @@ -430,7 +430,7 @@ public void getQual_should_be_sum_of_evidence() {
support.add(NRRP(tes, DP(0, 1, "2M", true, 0, 15, "5M", false)));
support.add(NRRP(tes, DP(0, 2, "2M", true, 0, 16, "5M", false)));
support.add(NRRP(nes, DP(0, 3, "2M", true, 0, 17, "10M", false)));
SAMRecord ass = AssemblyFactory.createAnchoredBreakend(pc, AES(), new SequentialIdGenerator("asm"), BWD, support, null, 0, 10, 5, B("CGTAAAAT"), new byte[] { 0,1,2,3,4,5,6,7});
SAMRecord ass = AssemblyFactory.createAnchoredBreakend(pc, AES(), new SequentialIdGenerator("asm"), BWD, support, fullSupport(support), 0, 10, 5, B("CGTAAAAT"), new byte[] { 0,1,2,3,4,5,6,7});
SoftClipEvidence asse = SoftClipEvidence.create(AES(), BWD, ass);

double[] supportQual = support.stream().mapToDouble(e -> e.getBreakendQual()).toArray();
Expand Down
19 changes: 17 additions & 2 deletions src/test/java/au/edu/wehi/idsv/SingleReadEvidenceTest.java
View file
Original file line number Diff line number Diff line change
Expand Up @@ -220,19 +220,34 @@ public void strandBias_should_match_read_strand() {
}
@Test
public void getBreakendReadOffsetInterval_should_consider_mapping_strand() {
Range<Integer> r = SingleReadEvidence.createEvidence(SES(), 0, Read(0, 1, "2M3S")).get(0).getBreakendAssemblyContigBreakpointInterval();
Range<Integer> r = SingleReadEvidence.createEvidence(SES(), 0, withSequence("NNNNN", Read(0, 1, "2M3S"))[0]).get(0).getBreakendAssemblyContigBreakpointInterval();
// 01234
// MMSSS
assertEquals(2, (int)r.lowerEndpoint());
assertEquals(2, (int)r.upperEndpoint());

r = SingleReadEvidence.createEvidence(SES(), 0, onNegative(Read(0, 1, "1M4S"))[0]).get(0).getBreakendAssemblyContigBreakpointInterval();
r = SingleReadEvidence.createEvidence(SES(), 0, withSequence("NNNNN", onNegative(Read(0, 1, "1M4S")))[0]).get(0).getBreakendAssemblyContigBreakpointInterval();
// 43210
// MSSSS
assertEquals(4, (int)r.lowerEndpoint());
assertEquals(4, (int)r.upperEndpoint());
}
@Test
public void getBreakendReadOffsetInterval_should_untemplated_inserted_Sequence() {
Range<Integer> r = SingleReadEvidence.createEvidence(SES(), 0, withAttr("SA", "polyA,10,+,3S2M,20,0", withSequence("NNNNN", Read(0, 1, "2M3S")))[0]).get(0).getBreakendAssemblyContigBreakpointInterval();
// 01234
// MMSSS
assertEquals(2, (int)r.lowerEndpoint());
assertEquals(3, (int)r.upperEndpoint());

r = SingleReadEvidence.createEvidence(SES(), 0, withAttr("SA", "polyA,10,+,2M3S,20,0", withSequence("NNNNN", onNegative(Read(0, 1, "1M4S"))))[0]).get(0).getBreakendAssemblyContigBreakpointInterval();
// 43210
// MSSSS
// SSSMM
assertEquals(2, (int)r.lowerEndpoint());
assertEquals(4, (int)r.upperEndpoint());
}
@Test
public void getBreakendReadOffsetInterval_should_consider_homology() {
Range<Integer> r = SingleReadEvidence.createEvidence(SES(), 0, withAttr("SA", "polyA,10,+,2S3M,20,0", withSequence("NAAANN", Read(0, 1, "2M4S")))[0]).get(0).getBreakendAssemblyContigBreakpointInterval();
// 01234
Expand Down
33 changes: 28 additions & 5 deletions src/test/java/au/edu/wehi/idsv/SplitReadEvidenceTest.java
View file
Original file line number Diff line number Diff line change
Expand Up @@ -404,8 +404,8 @@ public void should_handle_multiple_overlapping_fragments() {
Assert.assertEquals(e3.get(0).getEvidenceID(), e2.get(1).getRemoteEvidenceID());
}
@Test
public void should_pro_rata_assembly_support() {
SAMRecord primary = withSequence("NNNN", Read(1, 200, "1M3S"))[0];
public void untemplated_inserted_sequence_should_report_max_over_interval() {
SAMRecord primary = withSequence("NNNNNN", Read(1, 200, "1M5S"))[0];
primary.setMappingQuality(100);
SAMRecord r = Read(1, 100, "2M4S");
r.setMappingQuality(100);
Expand All @@ -415,14 +415,37 @@ public void should_pro_rata_assembly_support() {
r.setAttribute(SamTags.ASSEMBLY_DIRECTION, "f");
r.setAttribute(SamTags.ASSEMBLY_EVIDENCE_TYPE, new byte[] { 0, 1, 1});
r.setAttribute(SamTags.ASSEMBLY_EVIDENCE_CATEGORY, new int[] { 0, 0, 1});
r.setAttribute(SamTags.ASSEMBLY_EVIDENCE_OFFSET_START, new int[] { 3, 2, 3});
r.setAttribute(SamTags.ASSEMBLY_EVIDENCE_OFFSET_END, new int[] { 4, 3, 5});
r.setAttribute(SamTags.ASSEMBLY_EVIDENCE_OFFSET_START, new int[] { 0, 2, 3});
r.setAttribute(SamTags.ASSEMBLY_EVIDENCE_OFFSET_END, new int[] { 6, 3, 5});
r.setAttribute(SamTags.ASSEMBLY_EVIDENCE_QUAL, new float[] { 10, 20, 5});
r.setAttribute(SamTags.ASSEMBLY_EVIDENCE_EVIDENCEID, "1 2 3");
r.setAttribute(SamTags.ASSEMBLY_EVIDENCE_FRAGMENTID, "1 2 3");
AssemblyEvidenceSource aes = new MockAssemblyEvidenceSource(getContext(), ImmutableList.of(SES(0), SES(1)), new File("test.bam"));
SplitReadEvidence e = SplitReadEvidence.create(aes, r).get(0);
Assert.assertEquals(15, e.getBreakpointQual(), 0);
// 0 1 2 3 4 5 6
// M S S S primary
// S S S S M M realigned
Assert.assertEquals(10, e.getBreakpointQual(), 0);
}
@Test
public void homology_should_report_min_over_interval() {
SAMRecord primary = withSequence("NAAAAN", Read(1, 200, "1M5S"))[0];
primary.setMappingQuality(100);
SAMRecord r = withSequence("NAAAAN", Read(1, 300, "1S5M"))[0];
r.setMappingQuality(100);
r.setAttribute("SA", new ChimericAlignment(primary).toString());
r.setAttribute(SamTags.IS_ASSEMBLY, 1);
r.setAttribute(SamTags.ASSEMBLY_DIRECTION, "f");
r.setAttribute(SamTags.ASSEMBLY_EVIDENCE_TYPE, new byte[] { 0, 1, 1});
r.setAttribute(SamTags.ASSEMBLY_EVIDENCE_CATEGORY, new int[] { 0, 0, 1});
r.setAttribute(SamTags.ASSEMBLY_EVIDENCE_OFFSET_START, new int[] { 3, 2, 3});
r.setAttribute(SamTags.ASSEMBLY_EVIDENCE_OFFSET_END, new int[] { 4, 3, 5});
r.setAttribute(SamTags.ASSEMBLY_EVIDENCE_QUAL, new float[] { 10, 20, 5});
r.setAttribute(SamTags.ASSEMBLY_EVIDENCE_EVIDENCEID, "1 2 3");
r.setAttribute(SamTags.ASSEMBLY_EVIDENCE_FRAGMENTID, "1 2 3");
AssemblyEvidenceSource aes = new MockAssemblyEvidenceSource(getContext(), ImmutableList.of(SES(0), SES(1)), new File("test.bam"));
SplitReadEvidence e = SplitReadEvidence.create(aes, r).get(0);
Assert.assertEquals(0, e.getBreakpointQual(), 0);
}
@Test
public void FWD_unanchored_assembly_should_pro_rata_support_to_start_of_contig() {
Expand Down
Loading

0 comments on commit 0c4d318

Please sign in to comment.