Skip to content

Commit

Permalink
Merge cfee351 into bd08234
Browse files Browse the repository at this point in the history
  • Loading branch information
lorewar2 committed Jan 12, 2023
2 parents bd08234 + cfee351 commit e9c7720
Showing 1 changed file with 59 additions and 7 deletions.
66 changes: 59 additions & 7 deletions src/alignment/poa.rs
Expand Up @@ -195,7 +195,7 @@ impl Traceback {
let mut i = self.last.index() + 1;
let mut j = self.cols;

while i > 0 && j > 0 {
while i > 0 || j > 0 {
// push operation and edge corresponding to (one of the) optimal
// routes
ops.push(self.matrix[i][j].op);
Expand All @@ -212,13 +212,14 @@ impl Traceback {
j -= 1;
}
AlignmentOperation::Match(None) => {
break;
i -= 1;
j -= 1;
}
AlignmentOperation::Del(None) => {
j -= 1;
i -= 1;
}
AlignmentOperation::Ins(None) => {
i -= 1;
j -= 1;
}
}
}
Expand Down Expand Up @@ -432,11 +433,23 @@ impl<F: MatchFunc> Poa<F> {
/// * `aln` - The alignment of the new sequence to the graph
/// * `seq` - The sequence being incorporated
pub fn add_alignment(&mut self, aln: &Alignment, seq: TextSlice) {
let mut prev: NodeIndex<usize> = NodeIndex::new(0);
let head = Topo::new(&self.graph).next(&self.graph).unwrap();
let mut prev: NodeIndex<usize> = NodeIndex::new(head.index());
let mut i: usize = 0;
let mut edge_not_connected: bool = false;
for op in aln.operations.iter() {
match op {
AlignmentOperation::Match(None) => {
let node: NodeIndex<usize> = NodeIndex::new(0);
if (seq[i] != self.graph.raw_nodes()[head.index()].weight) && (seq[i] != b'X') {
let node = self.graph.add_node(seq[i]);
prev = node;
}
if edge_not_connected {
self.graph.add_edge(prev, node, 1);
prev = node;
edge_not_connected = false;
}
i += 1;
}
AlignmentOperation::Match(Some((_, p))) => {
Expand All @@ -452,15 +465,22 @@ impl<F: MatchFunc> Poa<F> {
*self.graph.edge_weight_mut(edge).unwrap() += 1;
}
None => {
// where the previous node was newly added
self.graph.add_edge(prev, node, 1);
if prev.index() != head.index() {
self.graph.add_edge(prev, node, 1);
}
}
}
prev = NodeIndex::new(*p);
}
i += 1;
}
AlignmentOperation::Ins(None) => {
let node = self.graph.add_node(seq[i]);
if edge_not_connected {
self.graph.add_edge(prev, node, 1);
}
prev = node;
edge_not_connected = true;
i += 1;
}
AlignmentOperation::Ins(Some(_)) => {
Expand All @@ -479,6 +499,7 @@ impl<F: MatchFunc> Poa<F> {
mod tests {
use super::*;
use crate::alignment::pairwise::Scoring;
use petgraph::dot::Dot;
use petgraph::graph::NodeIndex;

#[test]
Expand Down Expand Up @@ -587,4 +608,35 @@ mod tests {
.add_to_graph();
assert_eq!(aligner.alignment().score, 10);
}
#[test]
fn test_edge_cases() {
// case 1
let scoring = Scoring::new(-1, 0, |a: u8, b: u8| if a == b { 1i32 } else { -1i32 });
let mut aligner = Aligner::new(scoring, b"BBA");
aligner.global(b"AAA").add_to_graph();
let g = aligner.graph().map(|_, n| (*n) as char, |_, e| *e);
let dot = format!("{:?}", Dot::new(&g));
assert_eq!(dot, "digraph {\n 0 [ label = \"'B'\" ]\n 1 [ label = \"'B'\" ]\n 2 [ label = \"'A'\" ]\n 3 [ label = \"'A'\" ]\n 4 [ label = \"'A'\" ]\n 0 -> 1 [ label = \"1\" ]\n 1 -> 2 [ label = \"1\" ]\n 3 -> 4 [ label = \"1\" ]\n 4 -> 2 [ label = \"1\" ]\n}\n");
// case 2
let scoring = Scoring::new(-1, 0, |a: u8, b: u8| if a == b { 1i32 } else { -1i32 });
let mut aligner = Aligner::new(scoring, b"AAA");
aligner.global(b"ABA").add_to_graph();
let g = aligner.graph().map(|_, n| (*n) as char, |_, e| *e);
let dot = format!("{:?}", Dot::new(&g));
assert_eq!(dot, "digraph {\n 0 [ label = \"'A'\" ]\n 1 [ label = \"'A'\" ]\n 2 [ label = \"'A'\" ]\n 3 [ label = \"'B'\" ]\n 0 -> 1 [ label = \"1\" ]\n 1 -> 2 [ label = \"1\" ]\n 0 -> 3 [ label = \"1\" ]\n 3 -> 2 [ label = \"1\" ]\n}\n");
// case 3
let scoring = Scoring::new(-1, 0, |a: u8, b: u8| if a == b { 1i32 } else { -1i32 });
let mut aligner = Aligner::new(scoring, b"BBBBBAAA");
aligner.global(b"AAA").add_to_graph();
let g = aligner.graph().map(|_, n| (*n) as char, |_, e| *e);
let dot = format!("{:?}", Dot::new(&g));
assert_eq!(dot, "digraph {\n 0 [ label = \"'B'\" ]\n 1 [ label = \"'B'\" ]\n 2 [ label = \"'B'\" ]\n 3 [ label = \"'B'\" ]\n 4 [ label = \"'B'\" ]\n 5 [ label = \"'A'\" ]\n 6 [ label = \"'A'\" ]\n 7 [ label = \"'A'\" ]\n 0 -> 1 [ label = \"1\" ]\n 1 -> 2 [ label = \"1\" ]\n 2 -> 3 [ label = \"1\" ]\n 3 -> 4 [ label = \"1\" ]\n 4 -> 5 [ label = \"1\" ]\n 5 -> 6 [ label = \"2\" ]\n 6 -> 7 [ label = \"2\" ]\n}\n");
// case 4
let scoring = Scoring::new(-1, 0, |a: u8, b: u8| if a == b { 1i32 } else { -1i32 });
let mut aligner = Aligner::new(scoring, b"AAA");
aligner.global(b"BBBBBAAA").add_to_graph();
let g = aligner.graph().map(|_, n| (*n) as char, |_, e| *e);
let dot = format!("{:?}", Dot::new(&g));
assert_eq!(dot, "digraph {\n 0 [ label = \"'A'\" ]\n 1 [ label = \"'A'\" ]\n 2 [ label = \"'A'\" ]\n 3 [ label = \"'B'\" ]\n 4 [ label = \"'B'\" ]\n 5 [ label = \"'B'\" ]\n 6 [ label = \"'B'\" ]\n 7 [ label = \"'B'\" ]\n 0 -> 1 [ label = \"2\" ]\n 1 -> 2 [ label = \"2\" ]\n 3 -> 4 [ label = \"1\" ]\n 4 -> 5 [ label = \"1\" ]\n 5 -> 6 [ label = \"1\" ]\n 6 -> 7 [ label = \"1\" ]\n 7 -> 0 [ label = \"1\" ]\n}\n");
}
}

0 comments on commit e9c7720

Please sign in to comment.