Sequence-to-graph alignment based copy number calling using a network flow formulation

The paper introduces Floco, a novel method that improves copy number calling accuracy on genome graphs by combining node-level read depth probabilities with a network flow formulation solved via integer linear programming, thereby overcoming the limitations of traditional linear reference mapping and simple graph node predictions.

The Gist

Imagine your genome (your body's instruction manual) isn't just a single, straight line of text. Instead, it's a massive, complex choose-your-own-adventure book. Some chapters are identical for everyone, but others have "forks in the road" where different people have different versions of the story. Sometimes, a whole paragraph might be missing in one person's book, or a specific sentence might be repeated three times instead of once.

This is the world of Copy Number Variation (CNV). It's when people have different numbers of copies of specific DNA segments. These differences can cause diseases, like hearing loss, or unique traits.

The Old Problem: The Broken Map

Traditionally, scientists tried to figure out these copy numbers by comparing a person's DNA to a single, standard "reference" map (a linear genome).

Think of this like trying to navigate a city using a map of a different city.

• If your city has a new bridge that isn't on the old map, your GPS gets confused.
• If your city has a street that loops around three times, but the map only shows it once, the GPS thinks you're lost.

Because the "standard map" doesn't show all the possible paths or variations, scientists often got the copy numbers wrong. They would count how many times a piece of DNA appeared in the data (read depth), but if that piece didn't fit the standard map, the count would be messy and inaccurate.

The New Solution: Floco (The Traffic Controller)

The authors of this paper created a new tool called Floco. Instead of using a flat map, Floco uses a Genome Graph.

The Analogy: A Subway System
Imagine the genome graph as a subway system:

• Stations are chunks of DNA.
• Tracks connect the stations.
• Some stations are on a straight line; others are part of a complex web where trains can switch tracks.

When we sequence a person's DNA, we get millions of tiny puzzle pieces (reads). Floco places these pieces onto the subway map.

The Challenge:
If you just count how many trains stop at a specific station, you might get a weird number because:

1. Some trains might have missed the stop (sequencing errors).
2. Some tracks might be confusing, causing trains to pile up or disappear (misalignments).
3. The map itself might have loops that confuse the count.

Floco's Magic: The Network Flow
Floco doesn't just count trains at one station in isolation. It acts like a super-smart traffic controller looking at the entire subway system at once.

1. The "Flow" Concept: Imagine water flowing through pipes. The water represents the DNA. The rules of physics say water can't just appear or disappear in the middle of a pipe; what goes in must come out.
2. Fixing the Errors: Floco uses a mathematical trick (called Integer Linear Programming) to find the most logical "flow" of DNA through the entire graph.
   • If a station looks like it has 0 trains (Copy Number 0) but the stations before and after it are full, Floco realizes, "Wait, that doesn't make sense! The water must be flowing through here." It corrects the error.
   • If a station looks like it has 10 trains, but the pipes leading to it can only handle 2, Floco realizes, "That's a traffic jam caused by a mistake. Let's smooth it out."

Why This Matters

The paper tested Floco on real human and potato data (yes, potatoes have complex genomes too!). Here is what they found:

• Accuracy Boost: By looking at the whole picture (the flow) instead of just individual spots, Floco improved accuracy by up to 43% compared to old methods.
• Consistency: It made sure that if a path exists, the DNA counts make sense along that entire path. It stopped the "glitches" where a path would suddenly disappear and reappear.
• Speed: It's fast enough to run on a regular laptop, which is great for doctors and researchers.

The Bottom Line

Floco is like upgrading from a static, flat map to a dynamic, 3D traffic control system.

Instead of getting confused by the twists, turns, and missing pieces of our complex genetic "choose-your-own-adventure" books, Floco understands the whole story. It ensures that when we count how many copies of a gene a person has, we aren't just guessing based on a broken map, but are calculating the truth based on how the DNA actually flows through the system. This helps us better understand diseases and validate our genetic maps.

Technical Summary

1. Problem Statement

Copy number (CN) variation is a critical genomic feature associated with phenotypic differences and disease. Traditional CN calling relies on mapping sequencing reads to a linear reference genome and estimating CN based on read depth. This approach suffers from two major limitations:

1. Reference Bias: Linear references cannot represent sequences or rearrangements absent from the reference, leading to mapping bias and inaccurate CN estimates for non-reference alleles.
2. Inconsistency in Graphs: While pangenome graphs (representing diverse sequences) solve the bias issue, estimating CN for individual graph nodes based solely on local read depth often yields inconsistent results. Factors like sequencing errors, misalignments, and misassemblies create noise, making it difficult to ensure that CN values across the graph form a valid biological walk (e.g., ensuring flow conservation through the graph).

2. Methodology: Floco

The authors propose Floco, a method that performs CN calling on genome graphs by combining local read-depth statistics with a network flow formulation solved via Integer Linear Programming (ILP).

Input

• A genome graph (in GFA format).
• Sequence-to-graph alignments (in GAF format).

Workflow

1. Preprocessing & Coverage Calculation:

   • Nodes are clipped to remove redundant sequences (similar to graph bluntification).
   • Nodes longer than a threshold (default 100 bp) are split into bins.
   • Coverage ($k_v$) is calculated as the sum of read depths across all base pairs in a node or bin.

2. Probabilistic Modeling:

   • The authors model the distribution of read coverage using a Negative Binomial (NB) distribution for CN > 0.
   • For CN = 0 (deletions), they use a discretized Exponential distribution to model low-level noise caused by alignment errors.
   • Parameters ($r, q$) for the NB distribution are estimated via Maximum Likelihood Estimation (MLE) on filtered nodes, accounting for node length and expected read depth.
   • Node Probabilities: Using Bayes' theorem, the probability $p_{vc}$ of a node $v$ having copy number $c$ is calculated based on the observed coverage.

3. Network Flow Formulation (The Core Innovation):

   • The problem is modeled as a Maximum Probability Flow on a bidirected graph.
   • Graph Construction:
     • Each node $v$ has two ends (Left/Right).
     • A Supersource ($s$) and Supersink ($t$) are added, connected to every node end.
     • Edges between node ends represent the graph topology.
   • Objective Function: Maximize the sum of log-probabilities of the chosen CNs plus the costs of traversing edges.
     $$\text{Maximize } \sum_{v \in V} g(x_v; p_v) + \sum_{e \in E} f_e \cdot \xi_e$$
     Where $x_v$ is the integer CN variable, $g$ is a piecewise linear function mapping CN to log-probability, $f_e$ is flow, and $\xi_e$ is edge cost.
   • Constraints:
     • Flow Conservation: Inflow to a node end must equal outflow from the opposite end.
     • Capacity: Flow must be non-negative integers.
     • Edge Penalties: Edges supported by few reads (below expected depth) incur a high "expensive" cost penalty, discouraging flow through unreliable paths.
   • Solver: The problem is solved using the Gurobi ILP solver to find the globally optimal set of CNs that maximizes probability while respecting graph topology.

3. Key Contributions

• Novel Formulation: First method to apply network flow optimization specifically for CN calling on genome graphs, bridging the gap between local coverage estimation and global topological consistency.
• Robustness to Noise: By enforcing flow conservation, Floco corrects inconsistent CN assignments caused by sequencing errors or misalignments that plague standalone depth-based methods.
• Versatility: The tool works with various graph types (assembly graphs, pangenome graphs) and sequencing technologies (HiFi, ONT, and synthetic chopped assemblies).
• Efficiency: Despite using ILP, the method is computationally efficient, running in minutes on standard hardware for large graphs.

4. Results

The authors evaluated Floco on 15 datasets involving three graphs:

1. HG01114 Assembly Graph (Human diploid).
2. Altus Assembly Graph (Potato tetraploid).
3. HPRC Pangenome Graph (Human, 90 haplotypes).

Key Findings:

• Accuracy Improvement: Adding the network flow formulation increased CN prediction accuracy by up to 43% compared to read-depth estimation alone.
  • Example: On HiFi data for the HPRC graph, node-level accuracy improved from 95.2% to 99.3%.
• Low Coverage Performance: The flow formulation significantly outperformed local methods at low coverage (1×). For ONT data at 1×, accuracy jumped from 53.4% to 96.9% (a ~15-fold reduction in error rate) by aggregating signals from neighboring nodes.
• Concordance: Predictions from three different data sources (HiFi, ONT, and chopped assemblies) showed high concordance:
  • 93.2% concordance for the human assembly graph.
  • 92.9% concordance for the potato assembly graph.
  • 99.9% concordance for the CHM13 in-sample pangenome alignment.
• Quality Control: Floco successfully identified misassembled nodes (assigned CN=0 despite having some coverage) in assembly graphs, flagging bubbles and structural errors that standard tools might miss.
• Performance:
  • Runtime: 8–22 minutes for assembly graphs; 40 mins–3 hours for complex pangenome graphs.
  • Memory: Peak usage remained under 53 GB even for the largest graphs.

5. Significance

• Fills a Critical Gap: Floco provides the first dedicated tool for CN calling on genome graphs, a necessary step as the field shifts from linear references to pangenomes.
• Enables Better Genotyping: It allows for accurate genotyping of complex structural variations (like segmental duplications) that are difficult to resolve with linear references.
• Assembly Validation: It serves as a powerful Quality Control (QC) tool for genome assembly graphs, identifying errors and inconsistencies in the assembly process.
• Future-Proofing: As pangenome references become standard, Floco offers a robust framework for analyzing CN variation across diverse populations without the biases inherent in linear reference mapping.

In summary, Floco transforms CN calling from a local, noise-prone estimation task into a global, topology-aware optimization problem, significantly improving accuracy and consistency for genomic analysis using graph-based references.