Implements alignedRegionsOnly parameter. This fixes #347

milaboratory · Feb 13, 2018 · 3ca683c · 3ca683c
1 parent 67bbafa
commit 3ca683c
Show file tree

Hide file tree

Showing 4 changed files with 58 additions and 19 deletions.
diff --git a/src/main/java/com/milaboratory/mixcr/assembler/fullseq/FullSeqAssembler.java b/src/main/java/com/milaboratory/mixcr/assembler/fullseq/FullSeqAssembler.java
@@ -6,8 +6,10 @@
 import com.milaboratory.core.alignment.*;
 import com.milaboratory.core.mutations.MutationsBuilder;
 import com.milaboratory.core.sequence.*;
-import com.milaboratory.mixcr.basictypes.*;
-import com.milaboratory.mixcr.cli.ActionExportClonesPretty;
+import com.milaboratory.mixcr.basictypes.Clone;
+import com.milaboratory.mixcr.basictypes.VDJCAlignments;
+import com.milaboratory.mixcr.basictypes.VDJCHit;
+import com.milaboratory.mixcr.basictypes.VDJCPartitionedSequence;
 import com.milaboratory.mixcr.vdjaligners.VDJCAlignerParameters;
 import gnu.trove.impl.Constants;
 import gnu.trove.iterator.TIntIntIterator;
@@ -441,7 +443,9 @@ private Clone buildClone(double count, BranchSequences targets) {
 
             if (range.getTo() < nLeftDummies + lengthV) {
                 boolean floatingLeftBound =
-                        i == 0 && alignerParameters.getVAlignerParameters().getParameters().isFloatingLeftBound();
+                        !parameters.alignedRegionsOnly
+                                && i == 0
+                                && alignerParameters.getVAlignerParameters().getParameters().isFloatingLeftBound();
 
                 // Can be reduced to a single statement
                 if (range.getFrom() < nLeftDummies)
@@ -474,7 +478,9 @@ else if (floatingLeftBound)
                  */
 
                 boolean vFloatingLeftBound =
-                        i == 0 && alignerParameters.getVAlignerParameters().getParameters().isFloatingLeftBound();
+                        !parameters.alignedRegionsOnly
+                                && i == 0
+                                && alignerParameters.getVAlignerParameters().getParameters().isFloatingLeftBound();
 
                 // Can be reduced to a single statement
                 if (range.getFrom() < nLeftDummies)
@@ -507,7 +513,9 @@ else if (vFloatingLeftBound)
                  */
 
                 boolean jFloatingRightBound =
-                        i == targets.ranges.length - 1 && alignerParameters.getJAlignerParameters().getParameters().isFloatingRightBound();
+                        !parameters.alignedRegionsOnly
+                                && i == targets.ranges.length - 1
+                                && alignerParameters.getJAlignerParameters().getParameters().isFloatingRightBound();
 
                 if (range.getTo() >= rightAssemblingFeatureBound + jLength)
                     // This target contain extra non-J nucleotides on the right
@@ -535,7 +543,9 @@ else if (jFloatingRightBound)
                             targets.assemblingFeatureOffset + assemblingFeatureLength, sequence.size() - (targets.assemblingFeatureOffset + assemblingFeatureLength));
             } else if (range.getFrom() > rightAssemblingFeatureBound) {
                 boolean floatingRightBound =
-                        i == targets.ranges.length - 1 && alignerParameters.getJAlignerParameters().getParameters().isFloatingRightBound();
+                        !parameters.alignedRegionsOnly
+                                && i == targets.ranges.length - 1
+                                && alignerParameters.getJAlignerParameters().getParameters().isFloatingRightBound();
 
                 if (range.getTo() >= rightAssemblingFeatureBound + jLength)
                     // This target contain extra non-J nucleotides on the right
@@ -942,39 +952,53 @@ List<PointSequence> toPointSequences(VDJCAlignments alignments, int iTarget) {
 
         List<PointSequence> points = new ArrayList<>();
         if (target.getPartitioning().isAvailable(assemblingFeature.getFirstPoint())) {
+            // This target contains left edge of the assembling feature
             int leftStop = target.getPartitioning().getPosition(assemblingFeature.getFirstPoint());
             if (hasV) {
                 if (vAlignment != null)
                     toPointSequencesByAlignments(points,
                             vAlignment,
                             targetSeq,
-                            new Range(0, leftStop),
+                            new Range(
+                                    parameters.alignedRegionsOnly ? vAlignment.getSequence2Range().getFrom() : 0,
+                                    leftStop),
                             nLeftDummies);
-            } else
+            } else if (!parameters.alignedRegionsOnly)
                 toPointSequencesNoAlignments(points, targetSeq, new Range(0, leftStop), nLeftDummies - leftStop);
         } else if (hasV && vAlignment != null)
+            // This target ends before beginning (left edge) of the assembling feature
             toPointSequencesByAlignments(points,
                     vAlignment,
                     targetSeq,
-                    new Range(0, vAlignment.getSequence2Range().getTo()),
+                    new Range(
+                            parameters.alignedRegionsOnly ? vAlignment.getSequence2Range().getFrom() : 0,
+                            vAlignment.getSequence2Range().getTo()),
                     nLeftDummies);
 
         if (target.getPartitioning().isAvailable(assemblingFeature.getLastPoint())) {
+            // This target contains right edge of the assembling feature
             int rightStart = target.getPartitioning().getPosition(assemblingFeature.getLastPoint());
             if (hasJ) {
                 if (jAlignment != null)
                     toPointSequencesByAlignments(points,
                             jAlignment,
                             targetSeq,
-                            new Range(rightStart, targetSeq.size()),
+                            new Range(rightStart,
+                                    parameters.alignedRegionsOnly
+                                            ? jAlignment.getSequence2Range().getTo()
+                                            : targetSeq.size()),
                             nLeftDummies + lengthV + assemblingFeatureLength - jOffset);
             } else
                 toPointSequencesNoAlignments(points, targetSeq, new Range(rightStart, targetSeq.size()), nLeftDummies + lengthV + assemblingFeatureLength - rightStart);
         } else if (hasJ && jAlignment != null)
+            // This target starts after the end (right edge) of the assembling feature
             toPointSequencesByAlignments(points,
                     jAlignment,
                     targetSeq,
-                    new Range(jAlignment.getSequence2Range().getFrom(), targetSeq.size()),
+                    new Range(jAlignment.getSequence2Range().getFrom(),
+                            parameters.alignedRegionsOnly
+                                    ? jAlignment.getSequence2Range().getTo()
+                                    : targetSeq.size()),
                     nLeftDummies + lengthV + assemblingFeatureLength - jOffset);
 
         return points;
@@ -986,8 +1010,8 @@ void toPointSequencesByAlignments(List<PointSequence> points,
                                       Range seq2Range,
                                       int offset) {
 
-//        if (seq2Range.length() == 0)
-//            return;
+        // if (seq2Range.length() == 0)
+        //     return;
 
         alignment = AlignmentUtils.shiftIndelsAtHomopolymers(alignment);
 

diff --git a/src/main/java/com/milaboratory/mixcr/assembler/fullseq/FullSeqAssemblerParameters.java b/src/main/java/com/milaboratory/mixcr/assembler/fullseq/FullSeqAssemblerParameters.java
@@ -51,6 +51,10 @@ public class FullSeqAssemblerParameters {
      * Region where variants are allowed
      */
     GeneFeature subCloningRegion = GeneFeature.VDJRegion;
+    /**
+     * Assemble only parts of sequences covered by alignments
+     */
+    boolean alignedRegionsOnly = false;
 
     @JsonCreator
     public FullSeqAssemblerParameters(
@@ -60,15 +64,16 @@ public FullSeqAssemblerParameters(
             @JsonProperty("alignedSequenceEdgeDelta") int alignedSequenceEdgeDelta,
             @JsonProperty("alignmentEdgeRegionSize") int alignmentEdgeRegionSize,
             @JsonProperty("minimalNonEdgePointsFraction") double minimalNonEdgePointsFraction,
-            @JsonProperty("subCloningRegion") GeneFeature subCloningRegion) {
+            @JsonProperty("subCloningRegion") GeneFeature subCloningRegion,
+            @JsonProperty("alignedRegionsOnly") boolean alignedRegionsOnly) {
         this.minimalQualityShare = minimalQualityShare;
         this.minimalSumQuality = minimalSumQuality;
         this.decisiveSumQualityThreshold = decisiveSumQualityThreshold;
         this.alignedSequenceEdgeDelta = alignedSequenceEdgeDelta;
         this.alignmentEdgeRegionSize = alignmentEdgeRegionSize;
         this.minimalNonEdgePointsFraction = minimalNonEdgePointsFraction;
         this.subCloningRegion = subCloningRegion;
-
+        this.alignedRegionsOnly = alignedRegionsOnly;
     }
 
     public double getMinimalQualityShare() {
@@ -119,10 +124,19 @@ public void setMinimalNonEdgePointsFraction(double minimalNonEdgePointsFraction)
         this.minimalNonEdgePointsFraction = minimalNonEdgePointsFraction;
     }
 
+    public boolean isAlignedRegionsOnly() {
+        return alignedRegionsOnly;
+    }
+
+    public void setAlignedRegionsOnly(boolean alignedRegionsOnly) {
+        this.alignedRegionsOnly = alignedRegionsOnly;
+    }
+
     @Override
     public FullSeqAssemblerParameters clone() {
         return new FullSeqAssemblerParameters(minimalQualityShare, minimalSumQuality, decisiveSumQualityThreshold,
-                alignedSequenceEdgeDelta, alignmentEdgeRegionSize, minimalNonEdgePointsFraction, subCloningRegion);
+                alignedSequenceEdgeDelta, alignmentEdgeRegionSize, minimalNonEdgePointsFraction, subCloningRegion,
+                alignedRegionsOnly);
     }
 
     private static Map<String, FullSeqAssemblerParameters> knownParameters;

diff --git a/src/main/resources/parameters/full_seq_assembler_parameters.json b/src/main/resources/parameters/full_seq_assembler_parameters.json
@@ -6,6 +6,7 @@
     "alignedSequenceEdgeDelta": 3,
     "alignmentEdgeRegionSize": 7,
     "minimalNonEdgePointsFraction": 0.25,
-    "subCloningRegion": "CDR3"
+    "subCloningRegion": "CDR3",
+    "alignedRegionsOnly": false
   }
 }
diff --git a/src/test/java/com/milaboratory/mixcr/assembler/fullseq/FullSeqAssemblerTest.java b/src/test/java/com/milaboratory/mixcr/assembler/fullseq/FullSeqAssemblerTest.java
@@ -13,7 +13,6 @@
 import com.milaboratory.mixcr.basictypes.VDJCAlignments;
 import com.milaboratory.mixcr.basictypes.VDJCAlignmentsFormatter;
 import com.milaboratory.mixcr.cli.ActionExportClonesPretty;
-import com.milaboratory.mixcr.cli.Main;
 import com.milaboratory.mixcr.util.RunMiXCR;
 import com.milaboratory.mixcr.vdjaligners.VDJCParametersPresets;
 import gnu.trove.set.hash.TIntHashSet;
@@ -35,7 +34,8 @@
 public class FullSeqAssemblerTest {
     static final FullSeqAssemblerParameters DEFAULT_PARAMETERS =
             new FullSeqAssemblerParameters(0.1, 80, 120,
-                    3, 7, 0.25, GeneFeature.VDJRegion);
+                    3, 7, 0.25, GeneFeature.VDJRegion,
+                    false);
 
     static final class MasterSequence {
         final int vPart, cdr3Part, jPart, cPart;