Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
14 changes: 7 additions & 7 deletions docs/docs/reference.md
Original file line number Diff line number Diff line change
Expand Up @@ -104,6 +104,12 @@ Align genomes into a multiple sequence alignment graph
Use "-" to write the uncompressed data to standard output (stdout). This is the default, if the argument is not provided.

Default value: `-`
* `-c`, `--circular` — Toggle if input genomes are circular
* `-f`, `--verify` — Sanity check: after construction verifies that the original sequences can be reconstructed exactly from the resulting pangraph. Raises an error otherwise
* `--no-progress-bar` — Toggle to disable progress bar. Notice that the progress bar is only displayed if the output is specified via the `-o` argument
* `--guide-tree <GUIDE_TREE>` — Path to a Newick-format guide tree to use instead of the default neighbor-joining tree.

When provided, the tree's topology drives the bottom-up graph-merging order. Each input FASTA sequence must appear exactly once as a leaf (matched by sequence name), and every internal node must be strictly bifurcating. Branch lengths and internal labels, if present, are ignored. Accepts plain or compressed files (gz, bz2, xz, zst).
* `-l`, `--len <INDEL_LEN_THRESHOLD>` — Minimum block size for alignment graph (in nucleotides)

Default value: `100`
Expand Down Expand Up @@ -134,12 +140,6 @@ Align genomes into a multiple sequence alignment graph
* `--max-alignment-attempts <MAX_ALIGNMENT_ATTEMPTS>` — For within-block alignment: number of times Nextclade will retry alignment with more relaxed results if alignment band boundaries are hit

Default value: `4`
* `-c`, `--circular` — Toggle if input genomes are circular
* `-f`, `--verify` — Sanity check: after construction verifies that the original sequences can be reconstructed exactly from the resulting pangraph. Raises an error otherwise
* `--no-progress-bar` — Toggle to disable progress bar. Notice that the progress bar is only displayed if the output is specified via the `-o` argument
* `--guide-tree <GUIDE_TREE>` — Path to a Newick-format guide tree to use instead of the default neighbor-joining tree.

When provided, the tree's topology drives the bottom-up graph-merging order. Each input FASTA sequence must appear exactly once as a leaf (matched by sequence name), and every internal node must be strictly bifurcating. Branch lengths and internal labels, if present, are ignored. Accepts plain or compressed files (gz, bz2, xz, zst).



Expand Down Expand Up @@ -169,6 +169,7 @@ Merge two pangenome graphs into a single one
Use "-" to write the uncompressed data to standard output (stdout). This is the default, if the argument is not provided.

Default value: `-`
* `-f`, `--verify` — Sanity check: after merging verifies that every genome of the two input graphs can still be reconstructed exactly from the merged graph. Raises an error otherwise
* `-l`, `--len <INDEL_LEN_THRESHOLD>` — Minimum block size for alignment graph (in nucleotides)

Default value: `100`
Expand Down Expand Up @@ -199,7 +200,6 @@ Merge two pangenome graphs into a single one
* `--max-alignment-attempts <MAX_ALIGNMENT_ATTEMPTS>` — For within-block alignment: number of times Nextclade will retry alignment with more relaxed results if alignment band boundaries are hit

Default value: `4`
* `-f`, `--verify` — Sanity check: after merging verifies that every genome of the two input graphs can still be reconstructed exactly from the merged graph. Raises an error otherwise



Expand Down
8 changes: 5 additions & 3 deletions packages/pangraph/src/commands/build/build_args.rs
Original file line number Diff line number Diff line change
Expand Up @@ -29,9 +29,6 @@ pub struct PangraphBuildArgs {
#[clap(value_hint = ValueHint::AnyPath)]
pub output_json: PathBuf,

#[clap(flatten)]
pub merge_params: GraphMergeParams,

/// Toggle if input genomes are circular
#[clap(long, short = 'c')]
pub circular: bool,
Expand All @@ -52,4 +49,9 @@ pub struct PangraphBuildArgs {
/// are ignored. Accepts plain or compressed files (gz, bz2, xz, zst).
#[clap(long, value_hint = ValueHint::FilePath)]
pub guide_tree: Option<PathBuf>,

// Declared last: it opens the "Alignment" help section, and `next_help_heading` applies to every
// argument declared after it.
#[clap(flatten)]
pub merge_params: GraphMergeParams,
}
8 changes: 5 additions & 3 deletions packages/pangraph/src/commands/merge/merge_args.rs
Original file line number Diff line number Diff line change
Expand Up @@ -36,11 +36,13 @@ pub struct PangraphMergeArgs {
#[clap(value_hint = ValueHint::AnyPath)]
pub output_json: PathBuf,

#[clap(flatten)]
pub merge_params: GraphMergeParams,

/// Sanity check: after merging verifies that every genome of the two input graphs can still be
/// reconstructed exactly from the merged graph. Raises an error otherwise.
#[clap(long, short = 'f')]
pub verify: bool,

// Declared last: it opens the "Alignment" help section, and `next_help_heading` applies to every
// argument declared after it.
#[clap(flatten)]
pub merge_params: GraphMergeParams,
}
60 changes: 24 additions & 36 deletions packages/pangraph/src/commands/merge/merge_run.rs
Original file line number Diff line number Diff line change
Expand Up @@ -5,14 +5,13 @@ use crate::make_error;
use crate::pangraph::graph_merging::merge_graphs;
use crate::pangraph::pangraph::Pangraph;
use crate::pangraph::pangraph_path::PangraphPath;
use crate::pangraph::reconstruct::{GenomeCoverage, reconstruct_by_name, verify_graph_sequences};
use crate::representation::seq::Seq;
use crate::pangraph::reconstruct::{path_ids_by_name, verify_graph_against_graphs};
use crate::utils::collections::find_duplicates;
use color_eyre::owo_colors::{AnsiColors, OwoColorize};
use color_eyre::{Help, SectionExt};
use eyre::{Report, WrapErr};
use log::{info, warn};
use std::collections::{BTreeMap, BTreeSet};
use std::collections::BTreeSet;
use std::path::Path;

pub fn merge_run(args: &PangraphMergeArgs) -> Result<(), Report> {
Expand All @@ -27,13 +26,6 @@ pub fn merge_run(args: &PangraphMergeArgs) -> Result<(), Report> {
.make_disjoint_from(&left)
.wrap_err("When making the identifiers of the two input graphs disjoint")?;

// Reconstruct the expected genomes from the (relabeled) inputs, before they are consumed by the
// merger. Keyed by path name: neither path ids nor record order survive a merge.
let expected = args
.verify
.then(|| expected_sequences(args, &left, &right))
.transpose()?;

info!(
"=== Graph merging start: graph sizes {} + {}",
left.paths.len(),
Expand All @@ -50,13 +42,18 @@ pub fn merge_run(args: &PangraphMergeArgs) -> Result<(), Report> {
merged.blocks.len()
);

if let Some(expected) = expected {
if args.verify {
#[cfg(debug_assertions)]
merged.sanity_check().wrap_err("When checking the merged graph")?;

verify_graph_sequences(&merged, &expected, GenomeCoverage::Complete)
// Compared against the inputs one genome at a time: reconstructing both graphs up front would
// hold their entire sequence content in memory for the duration of the check.
verify_graph_against_graphs(&merged, &[&left, &right])
.wrap_err("When verifying the sequences of the merged graph")?;
info!("Merged graph reconstructs all {} input genomes exactly", expected.len());
info!(
"Merged graph reconstructs all {} input genomes exactly",
merged.paths.len()
);
}

json_write_file(&args.output_json, &merged, JsonPretty(true))?;
Expand Down Expand Up @@ -95,6 +92,20 @@ fn merge_cmd_preliminary_checks(args: &PangraphMergeArgs, left: &Pangraph, right
);
}

// Verification matches genomes by name, so it needs every path of both inputs to carry one.
// Checked here rather than at verification time so that it fails before the expensive merge.
if args.verify {
for (graph, filepath) in [(left, &args.left_graph), (right, &args.right_graph)] {
path_ids_by_name(graph)
.wrap_err_with(|| format!("When resolving the genome names of graph '{}'", filepath.display()))
.with_section(|| {
"Verification matches genomes by name. Re-run without `--verify` to skip it."
.color(AnsiColors::Cyan)
.header("Suggestion:")
})?;
}
}

// Circularity is a per-path property, so mixing is structurally fine. It is however most often a
// mistake, since `build --circular` applies to all genomes of a graph at once.
if circularity(left) != circularity(right) {
Expand All @@ -110,26 +121,3 @@ fn merge_cmd_preliminary_checks(args: &PangraphMergeArgs, left: &Pangraph, right
fn circularity(graph: &Pangraph) -> BTreeSet<bool> {
graph.paths().map(PangraphPath::circular).collect()
}

/// Reconstructs the genomes of both input graphs, keyed by genome name.
///
/// Cross-graph name collisions are already rejected by `merge_cmd_preliminary_checks`, so the two
/// sets cannot overwrite each other here.
fn expected_sequences(
args: &PangraphMergeArgs,
left: &Pangraph,
right: &Pangraph,
) -> Result<BTreeMap<String, Seq>, Report> {
let mut expected = BTreeMap::new();
for (graph, filepath) in [(left, &args.left_graph), (right, &args.right_graph)] {
let genomes = reconstruct_by_name(graph)
.wrap_err_with(|| format!("When reconstructing the genomes of graph '{}'", filepath.display()))
.with_section(|| {
"Verification matches genomes by name. Re-run without `--verify` to skip it."
.color(AnsiColors::Cyan)
.header("Suggestion:")
})?;
expected.extend(genomes);
}
Ok(expected)
}
58 changes: 47 additions & 11 deletions packages/pangraph/src/pangraph/pangraph.rs
Original file line number Diff line number Diff line change
Expand Up @@ -12,6 +12,7 @@ use crate::utils::id::id;
use crate::utils::map_merge::{ConflictResolution, map_merge};
use crate::{make_internal_error, make_internal_report};
use eyre::{Report, WrapErr};
use log::warn;
use maplit::btreemap;
use schemars::JsonSchema;
use serde::{Deserialize, Serialize};
Expand Down Expand Up @@ -151,24 +152,35 @@ impl Pangraph {
/// the way to the final graph. Merging therefore requires namespacing one of the two graphs
/// first.
pub fn make_disjoint_from(self, other: &Self) -> Result<Self, Report> {
// Namespace under which the ids of this graph are re-derived.
const SALT: usize = 1;
const MAX_ATTEMPTS: usize = 8;

// Path ids are assigned above every path id of `other`, so they cannot collide by construction.
let path_id_offset = other.paths.keys().map(|pid| pid.0 + 1).max().unwrap_or(0);

let relabeled = self.relabel(SALT, path_id_offset)?;
// The salt must differ from the one used by any earlier merge whose relabeled ids survive into
// `other`, otherwise those ids get re-derived a second time and land on themselves. The path id
// offset provides that: a merger salted with `offset(P)` produces a graph with at least one more
// genome than `P`, so a later merge with that graph on the left uses a strictly larger offset.
let mut relabeled = self;
for attempt in 0..MAX_ATTEMPTS {
relabeled = relabeled.relabel(id((path_id_offset, attempt)), path_id_offset)?;

// Block and node ids are hashes, so a collision with `other` is possible in principle. At
// ~2^-64 per pair it does not happen in practice, but the check is cheap and the alternative is
// one graph silently overwriting a block of the other during the join.
if !relabeled.is_id_disjoint_from(other) {
return make_internal_error!(
"When making graphs id-disjoint: the relabeled graph still shares block or node ids with the other graph. This requires a 64-bit hash collision and should never happen."
);
if relabeled.is_id_disjoint_from(other) {
return Ok(relabeled);
}

// Relabeling is a composition of injective maps, so retrying on top of the previous attempt is
// safe, and path ids are assigned by rank rather than derived from the previous id, so they do
// not drift. Reaching this point needs either a genuine hash collision, or an operation that
// breaks the monotonicity of the offset: `simplify` drops paths without renumbering, so
// `build -> merge -> merge -> simplify -> merge` can bring the offset back to a value already
// used as a salt.
warn!("Identifier collision when relabeling graph ids (attempt {attempt}); retrying");
}

Ok(relabeled)
make_internal_error!(
"When making graphs id-disjoint: no collision-free relabeling of block and node ids found after {MAX_ATTEMPTS} attempts"
)
}

pub fn update(&mut self, u: &GraphUpdate) {
Expand Down Expand Up @@ -643,4 +655,28 @@ mod tests {
assert_eq!(joined.blocks.len(), 4);
assert_eq!(joined.nodes.len(), 4);
}

/// Appending to a graph that already absorbed a relabeled graph. With a constant salt the ids of
/// the third graph were re-derived exactly onto those the second one left behind, so the second
/// append always failed. Every graph here carries the same small ids, which is what `build`
/// assigns to blocks and nodes that never merge.
#[rstest]
fn test_make_disjoint_from_after_a_previous_merge() {
let first = colliding_graph(["a", "b"]);
let second = colliding_graph(["c", "d"]).make_disjoint_from(&first).unwrap();
let joined = crate::pangraph::graph_merging::graph_join(&first, &second);

let third = colliding_graph(["e", "f"]).make_disjoint_from(&joined).unwrap();

assert!(third.is_id_disjoint_from(&joined));
third.sanity_check().unwrap();
assert_eq!(third.path_ids().collect_vec(), vec![PathId(4), PathId(5)]);

// and the three of them can be joined without conflicts
let joined = crate::pangraph::graph_merging::graph_join(&joined, &third);
joined.sanity_check().unwrap();
assert_eq!(joined.paths.len(), 6);
assert_eq!(joined.blocks.len(), 6);
assert_eq!(joined.nodes.len(), 6);
}
}
Loading
Loading