Repository: vgteam/sequenceTubeMap (js)

Issue:
Account for deletions in per-node coverage
In https://github.com/vgteam/sequenceTubeMap/pull/324 we don't account for bases of a node being deleted in a read.

We should count the deleted bases and make sure that they aren't counted towards coverage on the node. If one read visits a node and half of the node sequence is deleted in the read, coverage should be 0.5.

Pull request:
Node coverage stats deletions
Fixes #327 

Patch:
diff --git a/src/util/tubemap.js b/src/util/tubemap.js
index 29f1874..7e1322e 100644
--- a/src/util/tubemap.js
+++ b/src/util/tubemap.js
@@ -3241,6 +3241,15 @@ export function coverage(node, allReads) {
     let readPathIndex = readVisit[1];
     let currRead = allReads[readNum];
     let numNodes = currRead.sequenceNew.length;
+    // identify deletion: if there's a deletion, then those bases must be deleted from total base count
+    for (let i = 0; i < currRead.sequenceNew.length; i += 1) {
+      currRead.sequenceNew[i].mismatches.forEach((mm) => {
+        if (mm.type === "deletion") {
+          console.log ("this read has a deletion", currRead)
+          countBases -= mm.length;
+        }
+      })
+    };
     //  if current node is the last node on the read path, add the finalNodeCoverLength number of bases
     if (numNodes === readPathIndex + 1) {
       countBases += currRead.finalNodeCoverLength;
@@ -3255,12 +3264,30 @@ export function coverage(node, allReads) {
     //  indicating read's starting and ending points within the node.
     let readNum = readVisit;
     let currRead = allReads[readNum];
+    // identify deletion: if there's a deletion, then those bases must be deleted from total base count
+    for (let i = 0; i < currRead.sequenceNew.length; i += 1) {
+      currRead.sequenceNew[i].mismatches.forEach((mm) => {
+        if (mm.type === "deletion") {
+          console.log ("this read has a deletion", currRead)
+          countBases -= mm.length;
+        }
+      })
+    };
     countBases += currRead.finalNodeCoverLength - currRead.firstNodeOffset;
   }
   // outgoing reads
   for (let readVisit of node.outgoingReads) {
     let readNum = readVisit[0];
     let currRead = allReads[readNum];
+    // identify deletion: if there's a deletion, then those bases must be deleted from total base count
+    for (let i = 0; i < currRead.sequenceNew.length; i += 1) {
+      currRead.sequenceNew[i].mismatches.forEach((mm) => {
+        if (mm.type === "deletion") {
+          console.log ("this read has a deletion", currRead)
+          countBases -= mm.length;
+        }
+      })
+    };
     // coverage of outgoing read would be the the distance between the end of the node and the
     //  starting point of the read within the node
     countBases += node.sequenceLength - currRead.firstNodeOffset;


Test patch:
diff --git a/src/util/tubemap.test.js b/src/util/tubemap.test.js
index d9d96a0..444bf89 100644
--- a/src/util/tubemap.test.js
+++ b/src/util/tubemap.test.js
@@ -607,4 +607,241 @@ describe('coverage', () => {
     expect(checkNodeExample(node, reads)).toBe(undefined);
     expect(coverage(node, reads)).toBe(6.57);
   })
+  // TEST #5
+  it('it can handle nodes with deletions in 1 read', async () => {
+    const node = {        
+      nodename: '3',    
+      sequenceLength: 6,
+      incomingReads: [],
+      internalReads: [0, 1],
+      outgoingReads: [[2, 0], [3, 0]]
+    }
+    const reads = [
+      {
+        "id": 1,
+        "sequenceNew": [
+          {
+            "nodeName": "1",
+            "mismatches": [
+              {
+                "type": "deletion",
+                "pos": 0,
+                "length": 2
+              }
+            ]
+          }
+        ],
+        "firstNodeOffset": 0,
+        "finalNodeCoverLength": 3,
+      },
+      {
+        "id": 2,
+        "sequenceNew": [
+          {
+            "nodeName": "1",
+            "mismatches": [
+              {
+                "type": "insertion",
+                "pos": 0,
+                "seq": "TCAAGAACAGTCATTCATG"
+              }
+            ]
+          }
+        ],
+        "firstNodeOffset": 1,
+        "finalNodeCoverLength": 5,
+      }, 
+      {
+        "id": 3,
+        "sequenceNew": [
+            {
+              "nodeName": "1",
+              "mismatches": []
+            },
+            {
+              "nodeName": "11",
+              "mismatches": []
+            },
+            {
+              "nodeName": "13",
+              "mismatches": []
+            }
+        ],
+        "firstNodeOffset": 2,
+        "finalNodeCoverLength": 6,
+      },
+      {
+        "id": 4,
+        "sourceTrackID": "1",
+        "sequence": [
+          "1",
+          "12",
+          "13"
+        ],
+        "sequenceNew": [
+            {
+              "nodeName": "1",
+              "mismatches": []
+            },
+            {
+              "nodeName": "12",
+              "mismatches": []
+            },
+            {
+              "nodeName": "13",
+              "mismatches": []
+            }
+        ],
+        "firstNodeOffset": 4,
+        "finalNodeCoverLength": 6,
+      }
+    ];
+    //expect(checkNodeExample(node, reads)).toBe(undefined);
+    //console.log("coverage(node, reads) for test 5:", coverage(node, reads))
+    expect(coverage(node, reads)).toBe(1.83);
+  })
+  // TEST #6
+  it('it can handle larger node with deletions in all reads', async () => {
+    const node = {        
+      nodename: '3',    
+      sequenceLength: 24,
+      incomingReads: [[0, 1]],
+      internalReads: [1, 2],
+      outgoingReads: [[3, 0], [4, 0]]
+    }
+    const reads = [
+      {
+        "id": 1,
+        "sequenceNew": [
+          {
+            "nodeName": "1",
+            "mismatches": [
+              {
+                "type": "deletion",
+                "pos": 0,
+                "length": 1
+              }
+            ]
+          }
+        ],
+        "firstNodeOffset": 1,
+        "finalNodeCoverLength": 10,
+      },
+      {
+        "id": 2,
+        "sequenceNew": [
+          {
+            "nodeName": "1",
+            "mismatches": [
+              {
+                "type": "deletion",
+                "pos": 0,
+                "length": 1
+              }
+            ]
+          }
+        ],
+        "firstNodeOffset": 4,
+        "finalNodeCoverLength": 15,
+      }, 
+      {
+        "id": 3,
+        "sequenceNew": [
+            {
+              "nodeName": "1",
+              "mismatches": [
+                {
+                  "type": "deletion",
+                  "pos": 0,
+                  "length": 1
+                }
+              ]
+            },
+            {
+              "nodeName": "11",
+              "mismatches": [
+                {
+                  "type": "deletion",
+                  "pos": 0,
+                  "length": 2
+                }
+              ]
+            },
+            {
+              "nodeName": "13",
+              "mismatches": [
+                {
+                  "type": "deletion",
+                  "pos": 0,
+                  "length": 1
+                }
+              ]
+            }
+        ],
+        "firstNodeOffset": 2,
+        "finalNodeCoverLength": 8,
+      },
+      {
+        "id": 4,
+        "sourceTrackID": "1",
+        "sequence": [
+          "1",
+          "12",
+          "13"
+        ],
+        "sequenceNew": [
+            {
+              "nodeName": "1",
+              "mismatches": [
+                {
+                  "type": "deletion",
+                  "pos": 0,
+                  "length": 1
+                }
+              ]
+            },
+            {
+              "nodeName": "12",
+              "mismatches": [
+                {
+                  "type": "deletion",
+                  "pos": 0,
+                  "length": 1
+                }
+              ]
+            },
+            {
+              "nodeName": "13",
+              "mismatches": [
+                {
+                  "type": "deletion",
+                  "pos": 0,
+                  "length": 1
+                }
+              ]
+            }
+        ],
+        "firstNodeOffset": 9,
+        "finalNodeCoverLength": 12,
+      },
+      {
+        "id": 5,
+        "sequenceNew": [
+          {
+            "nodeName": "1",
+            "mismatches": [
+              {
+                "type": "deletion",
+                "pos": 0,
+                "length": 2
+              }
+            ]
+          }
+        ],
+        "firstNodeOffset": 2,
+        "finalNodeCoverLength": 14,
+      },
+    ];
+    expect(coverage(node, reads)).toBe(2.79);
+  })
 })