Bio Genome Assembly Scaffolding

Version Compatibility

Reference examples tested with: YaHS 1.2+, SALSA2 2.3+, juicer_tools 1.22/2.0+, bwa 0.7.17+, chromap 0.2+, samtools 1.19+, RagTag 2.1+, tidk 0.2.3+, seqkit 2.6+.

Before using code patterns, verify installed versions match. If versions differ:

• CLI: <tool> --version then <tool> --help to confirm flags

hifiasm output filenames and YaHS/SALSA2 enzyme and resolution defaults have changed across releases; confirm output names and -e/-r behaviour against the installed version. The Hi-C-mapping flag set (bwa mem -5SP vs the Arima/VGP per-end pipeline) varies by vendor (Arima, Dovetail/Omni-C, Phase Genomics) - check the current mapping guide before pasting flags. If a command errors, introspect the tool and adapt rather than retrying.

Genome Scaffolding

"Turn my contigs into chromosomes" -> Order and orient contigs by long-range linking signal, insert N-gap spacers, then curate the contact map - scaffolding adds no sequence, so the output is a draft chromosome structure, not finished sequence.

• CLI: yahs contigs.fa hic_to_contigs.bam (Hi-C default), run_pipeline.py -a contigs.fa -b sorted.bed -m yes (SALSA2), ragtag.py scaffold ref.fa contigs.fa (reference-guided), allmaps path (linkage maps)

The Single Most Important Modern Insight -- Automated Scaffolding Produces a Draft; the Genome Is Trustworthy Only After Manual Curation of the Contact Map

YaHS finishes in minutes and emits a file literally named _scaffolds_final.agp. It is not final. Hi-C scaffolding is a contact-frequency inference - it orders/orients contigs by the polymer-physics signal that contact frequency decays with 1D distance - and it is error-prone exactly at the joins. Every reference-grade pipeline (VGP, Darwin Tree of Life, Earth BioGenome) treats the automated AGP as a first pass that a human curator then corrects in PretextView or Juicebox: breaking misjoins, flipping inversions, removing false duplications, assigning chromosomes by eye. Howe 2021 quantified the hidden labor across 111 VGP/DToL assemblies: on average 221 interventions per Gb (67 breaks, 105 joins, 49 false-duplication removals). Three load-bearing consequences:

1. The contact map is the QC, not a decoration. A correct chromosome shows one bright diagonal with smooth off-diagonal decay; off-diagonal blocks, anti-diagonal "bowties," and bleeding between chromosomes are errors to inspect and break (see Reading the Contact Map). An assembly paper whose methods say "scaffolded with [tool]" but never mention curation/Juicebox/PretextView is shipping a draft - downgrade any claim that depends on large-scale order (synteny, fusions, structural variants).

2. Scaffold N50 is not contig N50, and conflating them is the classic reviewer catch. Scaffolding inserts runs of N (estimated, often arbitrary length) and adds zero sequence. A genome can post "scaffold N50 = 60 Mb, chromosome-scale!" while its contig N50 is 200 kb - the contiguity is borne by gaps, not finished sequence. Always report contig N50, scaffold N50, gap count, and total N bases separately.

3. Reordering cannot rescue a bad contig. Scaffolders take contigs as given; a chimeric contig drags the wrong sequence into the wrong chromosome. Chimeras must be broken before scaffolding (see Chimera-Breaking).

Linking-Data Modality Taxonomy

Decision Tree by Available Linking Data

Modern reference-grade recipe (vertebrate/eukaryote): HiFi -> hifiasm contigs -> (optional Bionano hybrid scaffold) -> Hi-C scaffold (YaHS/SALSA2) -> manual curation in PretextView/JBAT -> gap-fill (TGS-GapCloser) -> QC (tidk telomeres, contact-map diagonal, contig-vs-scaffold N50).

Hi-C Read Mapping (the step people botch)

Hi-C reads are not a normal paired-end library: the two ends come from different genomic loci ligated together, so standard PE proper-pair/insert-size logic mis-handles them.

# Map each end SEPARATELY with no mate rescue/pairing; -5 reports the 5' (junction) portion as primary.
bwa index contigs.fa
bwa mem -5SP -T0 -t 16 contigs.fa hic_R1.fq.gz hic_R2.fq.gz | \
    samtools view -@ 8 -b - > aligned.bam
# Mandatory: mark/remove PCR + optical duplicates (Hi-C is duplicate-rich) before scaffolding.
samtools sort -@ 8 -n aligned.bam | samtools fixmate -m - - | \
    samtools sort -@ 8 - | samtools markdup -@ 8 - hic_to_contigs.bam
samtools index hic_to_contigs.bam

• Each end separately, no mate rescue (-5SP): -5 = 5' portion of a chimeric junction read as primary; -S/-P skip mate rescue and pairing. The Arima/VGP per-end pipeline (bwa mem on R1 and R2 independently, filter_five_end.pl, two_read_bam_combiner.pl) is the higher-fidelity equivalent.
• MAPQ filter at scaffolding (YaHS -q, SALSA2 reads MAPQ) to drop repeat-ambiguous reads.
• Enzyme awareness: Arima/Dovetail-DpnII cut at fixed motifs (GATC); the scaffolder uses -e to model legitimate read starts. Omni-C / DNase Hi-C is enzyme-free -> omit -e in YaHS, use -e DNASE in SALSA2.
• chromap (Zhang 2021 Nat Commun 12:6566) is the fast modern alternative: chromap --preset hic -r contigs.fa -x index -1 R1 -2 R2 ... aligns + dedups Hi-C ~10x faster; YaHS accepts its output.

YaHS (the current default)

samtools faidx contigs.fa
# Enzyme: -e GATC (DpnII/Dovetail), -e GATC,GANTC (Arima 2-enzyme); OMIT -e for Omni-C/enzyme-free.
yahs -e GATC -o out contigs.fa hic_to_contigs.bam
# Outputs: out_scaffolds_final.agp + out_scaffolds_final.fa, intermediate out_rNN.agp, out.bin
# Contig error-correction is ON by default (--no-contig-ec to disable; usually a mistake to disable).

# Juicebox prep for curation (.hic + .assembly to load in JBAT):
# NOTE: `juicer` here is the small utility BUNDLED with YaHS (operates on the .bin), NOT Aiden-lab Juicer;
# `juicer_tools.jar` below IS the separate Aiden-lab jar. Do not confuse the two.
juicer pre -a -o out_JBAT out.bin out_scaffolds_final.agp contigs.fa.fai 2>tmp_assembly.log
java -Xmx48G -jar juicer_tools.jar pre out_JBAT.txt out_JBAT.hic <(cat tmp_assembly.log | grep PRE_C_SIZE | awk '{print $2" "$3}')
# After hand-editing in Juicebox, export out_JBAT.review.assembly, then:
juicer post -o out_curated out_JBAT.review.assembly out_JBAT.liftover.agp contigs.fa

YaHS also accepts a BED (bamToBed) or PA5/pairs input. The _scaffolds_final.agp is the editable, reviewable object - curation edits the AGP, then regenerates FASTA. --telo-motif lets YaHS use telomere signal during scaffolding.

SALSA2 (graph-aware, iterative, breaks chimeras)

bamToBed -i aligned.bam > alignment.bed
sort -k4 alignment.bed > sorted.bed                 # MUST sort by read name (column 4)
python run_pipeline.py -a contigs.fa -l contigs.fa.fai -b sorted.bed \
       -e GATC -o salsa_out -m yes -g contigs_graph.gfa -p yes
# -e DNASE for Omni-C/enzyme-free; -m yes enables input-error (chimera) correction; -p yes writes AGP+FASTA per iteration

The signature -m yes step breaks input contigs where the Hi-C signal contradicts the assembler's join before scaffolding. Assembly headers must not contain :. -i sets iterations (default 3), -c min contig length (default 1000).

Reference-Guided Scaffolding (RagTag) -- Power and Peril

ragtag.py correct ref.fa contigs.fa                 # break query at putative misassemblies vs reference FIRST
ragtag.py scaffold ref.fa ragtag_output/ragtag.correct.fasta -t 16   # order/orient by homology (minimap2 default)
# Outputs ragtag.scaffold.agp + .fasta; inserts 100 bp N gaps by default; -C collapses unplaced into chr0.

The peril (load-bearing): RagTag orders contigs to match the reference, so by construction the output looks like the reference. Every real biological rearrangement that distinguishes the target organism - a chromosome fusion, fission, translocation, or inversion polymorphism - is silently erased and replaced by the reference's structure. Never reference-scaffold a genome whose karyotype or large-scale structure is a biological question - the pipeline "discovers" the reference's karyotype because it was imposed (this has happened in the literature). Acceptable uses: a structurally conserved close relative needing only a coordinate system for gene content; a scaffold hypothesis to then test/correct against Hi-C (ragtag.py merge -b hic.bam lets Hi-C arbitrate). On a genome with interesting structure the honest move is Hi-C + curation, full stop.

Genetic/Linkage Maps (ALLMAPS) and Optical Maps (Bionano)

ALLMAPS (Tang 2015 Genome Biol 16:3) integrates one or more genetic maps to order/orient scaffolds: allmaps merge map1.csv map2.csv -o maps.bed then allmaps path maps.bed contigs.fa. It is the gold standard for validating chromosome assignment and for species where Hi-C is hard; weight multiple maps by quality. Bionano Solve hybrid scaffolding (no journal paper; cite the Bionano Solve Theory of Operation: Hybrid Scaffold technical document) aligns DLS optical maps to in-silico contig maps, resolves conflicts, and merges into hybrid scaffolds - run it before Hi-C on large/repetitive genomes where megabase maps bridge repeat arrays Hi-C fumbles. The decision is economic (separate instrument + library), not purely technical.

Reading the Contact Map (the skill nobody writes down)

The map displays its own errors to a trained eye - render with PretextMap -> PretextView (interactive, emits a curated AGP), PretextSnapshot (static PNG), or Juicebox/JBAT.

• Correct chromosome: one bright diagonal, contacts decaying smoothly off-diagonal; inter-chromosome space dim and uniform.
• Misjoin: the diagonal breaks - signal jumps off-diagonal into a separate block at the junction (two diagonals stitched at a corner). Break it.
• Inversion: an anti-diagonal "bowtie"/butterfly - the gradient runs backwards through the segment. The most common curation fix once recognized; flip it.
• Translocation / wrong assignment: off-diagonal bright blocks bleeding between what should be separate chromosomes.
• False duplication (leaked haplotig): a faint off-diagonal stripe parallel to the diagonal contacting the same neighborhood. Remove it (49 of the 221 interventions/Gb).
• Hardest judgment: the map is noisy near the diagonal even when correct. A real weak join is faint but coherent with the right decay shape; mapping noise is structureless speckle. Treat ambiguous near-diagonal signal as a flag to inspect, not a number to trust.