refactor: reorganize reconsensus

packages/pangraph/src/pangraph/edits.rs

@@ -1,4 +1,5 @@
 use crate::io::seq::{complement, reverse_complement};
+use crate::make_error;
 use crate::representation::seq::Seq;
 use crate::representation::seq_char::AsciiChar;
 use crate::utils::interval::Interval;
@@ -127,6 +128,93 @@ impl Edit {
     self.subs.is_empty() && self.dels.is_empty() && self.inss.is_empty()
   }
 
+  /// Returns true if this edit contains any insertions or deletions (indels)
+  #[inline]
+  pub fn has_indels(&self) -> bool {
+    self.has_dels() || self.has_inss()
+  }
+
+  /// Returns true if this edit contains any deletions
+  #[inline]
+  pub fn has_dels(&self) -> bool {
+    !self.dels.is_empty()
+  }
+
+  /// Returns true if this edit contains any insertions
+  #[inline]
+  pub fn has_inss(&self) -> bool {
+    !self.inss.is_empty()
+  }
+
+  /// Returns true if this edit contains any substitutions
+  #[inline]
+  pub fn has_subs(&self) -> bool {
+    !self.subs.is_empty()
+  }
+
+  /// Checks if a position is deleted in this edit
+  #[inline]
+  pub fn is_position_deleted(&self, pos: usize) -> bool {
+    self.dels.iter().any(|d| d.contains(pos))
+  }
+
+  /// Adds reversion mutation when no substitutions exist at a position
+  pub fn add_reversion_if_not_deleted(&mut self, pos: usize, original: AsciiChar) {
+    if !self.is_position_deleted(pos) {
+      self.subs.push(Sub::new(pos, original));
+      self.subs.sort_by_key(|s| s.pos);
+    }
+  }
+
+  /// Removes substitution if it matches the new consensus character
+  pub fn remove_matching_substitution(&mut self, substitution: &Sub) {
+    if let Some(existing_sub) = self.subs.iter().find(|s| s.pos == substitution.pos) {
+      if existing_sub.alt == substitution.alt {
+        self
+          .subs
+          .retain(|sub| !(sub.pos == substitution.pos && sub.alt == substitution.alt));
+      }
+    }
+  }
+
+  /// Returns all substitutions at a specific position
+  fn subs_at_position(&self, pos: usize) -> Vec<Sub> {
+    self.subs.iter().filter(|s| s.pos == pos).cloned().collect()
+  }
+
+  /// Reconciles genome alignment when consensus changes due to a substitution.
+  ///
+  /// During reconsensus, when a position in the consensus sequence is updated with a new
+  /// character, this method adjusts the genome's edit to maintain correct alignment.
+  pub fn reconcile_substitution_with_consensus(
+    &mut self,
+    substitution: &Sub,
+    original: AsciiChar,
+  ) -> Result<(), Report> {
+    let subs_count = self.subs.iter().filter(|s| s.pos == substitution.pos).count();
+
+    match subs_count {
+      // If genome has no mutation at this position: adds reversion to original character
+      0 => self.add_reversion_if_not_deleted(substitution.pos, original),
+      // If genome has matching mutation: removes it (now matches new consensus)
+      1 => self.remove_matching_substitution(substitution),
+      // If genome has conflicting mutations: returns error for inconsistent state
+      _ => {
+        let subs_at_pos = self.subs_at_position(substitution.pos);
+        return make_error!(
+          "At position {}: sequence states disagree: {:}",
+          substitution.pos,
+          subs_at_pos
+            .iter()
+            .map(|sub| sub.alt.to_string())
+            .collect_vec()
+            .join(", ")
+        );
+      },
+    }
+    Ok(())
+  }
+
   /// Construct edit which consists of a deletion of length `len`
   pub fn deleted(len: usize) -> Self {
     Self {
@@ -1045,4 +1133,72 @@ mod tests {
     assert_eq!(mean_shift, 3);
     assert_eq!(bandwidth, Some(4));
   }
+
+  #[test]
+  fn test_has_indels() {
+    // Edit with no indels (only substitutions)
+    let edit_no_indels = Edit::new(vec![], vec![], vec![Sub::new(1, 'A')]);
+    assert!(!edit_no_indels.has_indels());
+
+    // Edit with deletions
+    let edit_with_dels = Edit::new(vec![], vec![Del::new(5, 2)], vec![]);
+    assert!(edit_with_dels.has_indels());
+
+    // Edit with insertions
+    let edit_with_inss = Edit::new(vec![Ins::new(10, "ATG")], vec![], vec![]);
+    assert!(edit_with_inss.has_indels());
+
+    // Edit with both insertions and deletions
+    let edit_with_both = Edit::new(vec![Ins::new(10, "ATG")], vec![Del::new(5, 2)], vec![Sub::new(1, 'A')]);
+    assert!(edit_with_both.has_indels());
+
+    // Empty edit
+    let edit_empty = Edit::empty();
+    assert!(!edit_empty.has_indels());
+  }
+
+  #[test]
+  fn test_has_dels() {
+    // Edit with no deletions
+    let edit_no_dels = Edit::new(vec![Ins::new(10, "ATG")], vec![], vec![Sub::new(1, 'A')]);
+    assert!(!edit_no_dels.has_dels());
+
+    // Edit with deletions
+    let edit_with_dels = Edit::new(vec![], vec![Del::new(5, 2)], vec![]);
+    assert!(edit_with_dels.has_dels());
+
+    // Empty edit
+    let edit_empty = Edit::empty();
+    assert!(!edit_empty.has_dels());
+  }
+
+  #[test]
+  fn test_has_inss() {
+    // Edit with no insertions
+    let edit_no_inss = Edit::new(vec![], vec![Del::new(5, 2)], vec![Sub::new(1, 'A')]);
+    assert!(!edit_no_inss.has_inss());
+
+    // Edit with insertions
+    let edit_with_inss = Edit::new(vec![Ins::new(10, "ATG")], vec![], vec![]);
+    assert!(edit_with_inss.has_inss());
+
+    // Empty edit
+    let edit_empty = Edit::empty();
+    assert!(!edit_empty.has_inss());
+  }
+
+  #[test]
+  fn test_has_subs() {
+    // Edit with no substitutions
+    let edit_no_subs = Edit::new(vec![Ins::new(10, "ATG")], vec![Del::new(5, 2)], vec![]);
+    assert!(!edit_no_subs.has_subs());
+
+    // Edit with substitutions
+    let edit_with_subs = Edit::new(vec![], vec![], vec![Sub::new(1, 'A')]);
+    assert!(edit_with_subs.has_subs());
+
+    // Empty edit
+    let edit_empty = Edit::empty();
+    assert!(!edit_empty.has_subs());
+  }
 }

packages/pangraph/src/pangraph/graph_merging.rs

@@ -166,7 +166,7 @@ pub fn self_merge(graph: Pangraph, args: &PangraphBuildArgs) -> Result<(Pangraph
 
   // update consensus and alignment of merged blocks.
   let merge_block_ids = new_blocks_dict.keys().copied().collect_vec();
-  reconsensus_graph(&mut graph, merge_block_ids, args).wrap_err("During reconsensus")?;
+  reconsensus_graph(&mut graph, &merge_block_ids, args).wrap_err("During reconsensus")?;
 
   Ok((graph, true))
 }

packages/pangraph/src/pangraph/pangraph_block.rs

@@ -1,16 +1,18 @@
 use crate::io::fasta::FastaRecord;
 use crate::io::json::{JsonPretty, json_write_str};
 use crate::io::seq::reverse_complement;
-use crate::pangraph::edits::Edit;
+use crate::pangraph::edits::{Del, Edit, Ins, Sub};
 use crate::pangraph::pangraph::Pangraph;
 use crate::pangraph::pangraph_node::NodeId;
 use crate::pangraph::pangraph_path::PathId;
 use crate::representation::seq::Seq;
 use crate::representation::seq_char::AsciiChar;
 use crate::utils::collections::has_duplicates;
+use crate::utils::interval::positions_to_intervals;
 use derive_more::{Display, From};
 use eyre::{Report, WrapErr};
 use getset::{CopyGetters, Getters};
+use itertools::Itertools;
 use maplit::btreemap;
 use schemars::JsonSchema;
 use serde::{Deserialize, Serialize};
@@ -195,6 +197,73 @@ impl PangraphBlock {
       })
     })
   }
+
+  /// Finds all majority edits (substitutions, deletions, insertions) in this block
+  pub fn find_majority_edits(&self) -> Edit {
+    Edit::new(
+      self.find_majority_insertions(),
+      self.find_majority_deletions(),
+      self.find_majority_substitutions(),
+    )
+  }
+
+  /// Helper method to check if a count represents a majority
+  #[inline]
+  pub fn is_majority(&self, count: usize) -> bool {
+    count > self.depth() / 2
+  }
+
+  /// Finds majority substitutions in this block
+  pub fn find_majority_substitutions(&self) -> Vec<Sub> {
+    let mut substitutions: Vec<_> = self
+      .alignments()
+      .values()
+      .flat_map(|edit| &edit.subs)
+      .map(|sub| (sub.pos, sub.alt))
+      .into_group_map()
+      .into_iter()
+      .filter_map(|(pos, alts)| {
+        let (alt, count) = alts.into_iter().counts().into_iter().max_by_key(|(_, count)| *count)?;
+        self.is_majority(count).then_some(Sub::new(pos, alt))
+      })
+      .collect();
+
+    substitutions.sort_by_key(|sub| sub.pos);
+    substitutions
+  }
+
+  /// Finds majority deletions in this block
+  pub fn find_majority_deletions(&self) -> Vec<Del> {
+    let majority_positions: Vec<usize> = self
+      .alignments()
+      .values()
+      .flat_map(|edit| edit.dels.iter().flat_map(|del| del.range()))
+      .counts()
+      .into_iter()
+      .filter_map(|(pos, count)| self.is_majority(count).then_some(pos))
+      .collect();
+
+    positions_to_intervals(&majority_positions)
+      .into_iter()
+      .map(|interval| Del::new(interval.start, interval.len()))
+      .collect()
+  }
+
+  /// Finds majority insertions in this block
+  pub fn find_majority_insertions(&self) -> Vec<Ins> {
+    let mut insertions: Vec<_> = self
+      .alignments()
+      .values()
+      .flat_map(|edit| &edit.inss)
+      .map(|insertion| (insertion.pos, insertion.seq.clone()))
+      .counts()
+      .into_iter()
+      .filter_map(|((pos, seq), count)| self.is_majority(count).then_some(Ins::new(pos, seq)))
+      .collect();
+
+    insertions.sort_by_key(|ins| ins.pos);
+    insertions
+  }
 }
 
 #[derive(Copy, Debug, Clone, PartialEq, Eq, PartialOrd, Ord, Serialize, Deserialize)]