From 6ce4b34bb499c964cee09ba4158985b685897290 Mon Sep 17 00:00:00 2001 From: Daniel Cameron Date: Mon, 13 May 2019 16:40:43 +1000 Subject: [PATCH] #213 Excluding local assembly anchoring from SC calculation; added ASC --- .../idsv/StructuralVariationCallBuilder.java | 9 ++++++++- .../VariantCallingConfiguration.java | 4 ++-- .../edu/wehi/idsv/vcf/VcfInfoAttributes.java | 3 ++- .../StructuralVariationCallBuilderTest.java | 19 +++++++++++++++++++ 4 files changed, 31 insertions(+), 4 deletions(-) diff --git a/src/main/java/au/edu/wehi/idsv/StructuralVariationCallBuilder.java b/src/main/java/au/edu/wehi/idsv/StructuralVariationCallBuilder.java index 05bc8f7ec..68b2f6091 100644 --- a/src/main/java/au/edu/wehi/idsv/StructuralVariationCallBuilder.java +++ b/src/main/java/au/edu/wehi/idsv/StructuralVariationCallBuilder.java @@ -277,10 +277,17 @@ private void updateNominalCallPosition() { supportingSR.stream().flatMap(l -> l.stream()).forEach(e -> processAnchor(allAnchoredBases, e.getSAMRecord())); supportingIndel.stream().flatMap(l -> l.stream()).forEach(e -> processAnchor(allAnchoredBases, e.getSAMRecord())); supportingDP.stream().flatMap(l -> l.stream()).forEach(e -> processAnchor(allAnchoredBases, e.getLocalledMappedRead())); - supportingAS.stream().forEach(e -> processAnchor(allAnchoredBases, e.getSAMRecord())); + //supportingAS.stream().forEach(e -> processAnchor(allAnchoredBases, e.getSAMRecord())); supportingRAS.stream().forEach(e -> processAnchor(allAnchoredBases, e.getSAMRecord())); supportingCAS.stream().forEach(e -> processAnchor(allAnchoredBases, e.getSAMRecord())); attribute(VcfInfoAttributes.SUPPORT_CIGAR, makeCigar(allAnchoredBases, nominalPosition).toString()); + RangeSet assemblyAnchoredBases = TreeRangeSet.create(); + // TODO: implement https://github.com/PapenfussLab/gridss/issues/213 so we can include local assemblies + //supportingAS.stream().forEach(e -> processAnchor(allAnchoredBases, e.getSAMRecord())); + supportingRAS.stream().forEach(e -> processAnchor(assemblyAnchoredBases, e.getSAMRecord())); + supportingCAS.stream().forEach(e -> processAnchor(assemblyAnchoredBases, e.getSAMRecord())); + attribute(VcfInfoAttributes.ASSEMBLY_SUPPORT_CIGAR, makeCigar(assemblyAnchoredBases, nominalPosition).toString()); + } private void updateAssemblySupportTracking() { diff --git a/src/main/java/au/edu/wehi/idsv/configuration/VariantCallingConfiguration.java b/src/main/java/au/edu/wehi/idsv/configuration/VariantCallingConfiguration.java index 0124612a4..b2e141fe6 100644 --- a/src/main/java/au/edu/wehi/idsv/configuration/VariantCallingConfiguration.java +++ b/src/main/java/au/edu/wehi/idsv/configuration/VariantCallingConfiguration.java @@ -170,8 +170,8 @@ public VariantContextDirectedEvidence applyConfidenceFilter(ProcessingContext pr filters.add(VcfFilter.NO_ASSEMBLY.filter()); } else if (v.getBreakpointEvidenceCountLocalAssembly() == 0 || v.getBreakpointEvidenceCountRemoteAssembly() == 0) { filters.add(VcfFilter.SINGLE_ASSEMBLY.filter()); - } else if (v.getBreakpointEvidenceCountLocalAssembly() + v.getBreakpointEvidenceCountRemoteAssembly() > 0 && - v.getBreakpointEvidenceCountReadPair() + v.getBreakpointEvidenceCountSoftClip() == 0) { + } else if (v.getBreakpointEvidenceCountLocalAssembly() + v.getBreakpointEvidenceCountRemoteAssembly() + v.getBreakpointEvidenceCountCompoundAssembly() > 0 && + v.getBreakpointEvidenceCountReadPair() + v.getBreakpointEvidenceCountSoftClip() + v.getBreakpointEvidenceCountIndel() == 0) { filters.add(VcfFilter.ASSEMBLY_ONLY.filter()); } if (variant.getPhredScaledQual() < processContext.getVariantCallingParameters().breakpointLowQuality) { diff --git a/src/main/java/au/edu/wehi/idsv/vcf/VcfInfoAttributes.java b/src/main/java/au/edu/wehi/idsv/vcf/VcfInfoAttributes.java index e524e8b04..2d57c2db0 100644 --- a/src/main/java/au/edu/wehi/idsv/vcf/VcfInfoAttributes.java +++ b/src/main/java/au/edu/wehi/idsv/vcf/VcfInfoAttributes.java @@ -60,7 +60,8 @@ public enum VcfInfoAttributes { SUPPORT_INTERVAL ("SI", VCFHeaderLineCount.UNBOUNDED, VCFHeaderLineType.Integer, "Support interval offsets from breakend position in which at least one supporting read/read pair/assembly is mapped."), REMOTE_SUPPORT_INTERVAL ("RSI", VCFHeaderLineCount.UNBOUNDED, VCFHeaderLineType.Integer, "Support interval offsets of partner breakend."), INEXACT_HOMPOS ("IHOMPOS", 2, VCFHeaderLineType.Integer, "Position of inexact homology"), - SUPPORT_CIGAR ("SC", 1, VCFHeaderLineType.String, "CIGAR for displaying anchoring alignment of any contributing evidence and microhomologies."), + SUPPORT_CIGAR ("SC", 1, VCFHeaderLineType.String, "CIGAR for displaying anchoring alignment of any contributing evidence and microhomologies. Local assemblies are excluded due to https://github.com/PapenfussLab/gridss/issues/213"), + ASSEMBLY_SUPPORT_CIGAR ("ASC", 1, VCFHeaderLineType.String, "CIGAR encoding assembly contig anchoring alignments. Local assemblies are excluded due to https://github.com/PapenfussLab/gridss/issues/213."), BREAKEND_ASSEMBLY_ID ("BEID", VCFHeaderLineCount.UNBOUNDED, VCFHeaderLineType.String, "Identifiers of assemblies supporting the variant."), BREAKEND_ASSEMBLY_ID_LOCAL_CONTIG_OFFSET ("BEIDL", VCFHeaderLineCount.UNBOUNDED, VCFHeaderLineType.Integer, "Local chimeric alignment offset of corresponding BEID assembly."), BREAKEND_ASSEMBLY_ID_REMOTE_CONTIG_OFFSET ("BEIDH", VCFHeaderLineCount.UNBOUNDED, VCFHeaderLineType.Integer, "Remote chimeric alignment offset of corresponding BEID assembly."), diff --git a/src/test/java/au/edu/wehi/idsv/StructuralVariationCallBuilderTest.java b/src/test/java/au/edu/wehi/idsv/StructuralVariationCallBuilderTest.java index e9d73336e..0e0b0f320 100644 --- a/src/test/java/au/edu/wehi/idsv/StructuralVariationCallBuilderTest.java +++ b/src/test/java/au/edu/wehi/idsv/StructuralVariationCallBuilderTest.java @@ -728,6 +728,25 @@ public void anchor_cigar_should_use_local_coordinates() { assertEquals("9M1X", vc.getAttribute(VcfInfoAttributes.SUPPORT_CIGAR.attribute())); } @Test + public void assembly_anchor_cigar_should_use_only_remote_assemblies() { + ProcessingContext pc = getContext(); + AssemblyEvidenceSource aes = AES(pc); + StructuralVariationCallBuilder builder = new StructuralVariationCallBuilder(pc, (VariantContextDirectedEvidence)new IdsvVariantContextBuilder(getContext()) {{ + breakpoint(new BreakpointSummary(0, FWD, 100, 1, BWD, 200), ""); + phredScore(10); + }}.make()); + List support = Lists.newArrayList( + SR(Read(1, 200, "3S7M"), Read(0, 98, "3M7S")) + ); + SAMRecord ass = AssemblyFactory.createAnchoredBreakend(pc, aes, new SequentialIdGenerator("asm"), BWD, support, fullSupport(support), 1, 200, 7, B("NNNNNNNNNN"), B("NNNNNNNNNN")); + SAMRecord beass = withMapq(44, Read(0, 98, "3M"))[0]; + incorporateRealignment(AES(), ass, ImmutableList.of(beass)); + + builder.addEvidence(SingleReadEvidence.createEvidence(AES(), 1, beass).get(0)); + VariantContextDirectedEvidence vc = builder.make(); + assertEquals("2M1X", vc.getAttribute(VcfInfoAttributes.ASSEMBLY_SUPPORT_CIGAR.attribute())); + } + @Test @Ignore("The assembler doesn't know which orientation it was assembling.") public void spanning_assemblies_should_use_original_parent_assembly_direction_to_determine_local_remote_status() { IndelEvidence r = IE(withMapq(40, Read(2, 90, "10M1I10M"))[0]);