Synteny visualization (all-vs-all minimap2)
The tutorials target the JBrowse v5 beta, and the v4.3.0 release
on the download page lacks some of what they show. To install the
beta of JBrowse Web, run npm install -g @jbrowse/cli@next, then
jbrowse create jbrowse2 --branch v5.0.0-beta.13. Desktop beta builds are coming soon.
Strains of one bacterial species share most of their chromosome and differ by gene islands one strain picked up and by stretches that flipped. We align five E. coli strains against each other with minimap2, stack them in a linear synteny view, and read where they differ: Sakai's Shiga-toxin prophage, the phenylacetate (paa) operon three strains lack, and IAI39's inversions. The input is one all-vs-all PAF, minimap2's alignment of every genome against every other:
- align the five strains against each other with
minimap2 -Xto build the PAF - load it with
MultiGenomePAFAdapterand stack the five assemblies as rows - add a gene track per strain, then read one strain against the rest in a single pileup
Prerequisites
- a JBrowse to open the files in: Desktop takes a local file by path, Web through Add track
- the NCBI
datasetsCLI minimap2samtools- htslib (
bgzip,tabix) unzipnode, for the JBrowse CLI
Where the data comes from
Five E. coli RefSeq assemblies, each fetched by accession with the datasets
CLI.
The build script takes these files from their URLs, so there is nothing to download by hand.
- K12: ftp.ncbi.nlm.nih.gov/…/GCF_000005845.2_ASM584v2https://ftp.ncbi.nlm.nih.gov/genomes/all/GCF/000/005/845/GCF_000005845.2_ASM584v2/
- Sakai: ftp.ncbi.nlm.nih.gov/…/GCF_000008865.2_ASM886v2https://ftp.ncbi.nlm.nih.gov/genomes/all/GCF/000/008/865/GCF_000008865.2_ASM886v2/
- CFT073: ftp.ncbi.nlm.nih.gov/…/GCF_000007445.1_ASM744v1https://ftp.ncbi.nlm.nih.gov/genomes/all/GCF/000/007/445/GCF_000007445.1_ASM744v1/
- NCTC86: ftp.ncbi.nlm.nih.gov/…/GCF_002007705.1_ASM200770v1https://ftp.ncbi.nlm.nih.gov/genomes/all/GCF/002/007/705/GCF_002007705.1_ASM200770v1/
- IAI39: ftp.ncbi.nlm.nih.gov/…/GCF_000026345.1_ASM2634v1https://ftp.ncbi.nlm.nih.gov/genomes/all/GCF/000/026/345/GCF_000026345.1_ASM2634v1/
Producing an all-vs-all PAF
minimap2 aligning
PanSN-named genomes against
themselves writes an all-vs-all PAF, as does the mapping step of
PGGB, the PanGenome Graph Builder. PanSN
names every sequence sample#haplotype#contig, e.g. K12#1#chr.
The script reduces each of the five assemblies to
one chr record. The five strain FASTAs become the JBrowse assemblies as-is;
the PanSN names exist only in the concatenated copy minimap2 aligns, with the
haplotype 1 throughout:
for strain in K12 Sakai CFT073 NCTC86 IAI39; do
# '>chr' -> '>K12#1#chr'
awk -v s="$strain" '/^>/{print ">" s "#1#chr"; next} {print}' "$strain.fa"
done > all.fa
# -c: emit the base-level CIGAR the linear synteny view needs
# -X: skip self-alignments and the reciprocal half of every pair
minimap2 -c -x asm20 -X all.fa all.fa > all_vs_all.pafSetting up the five assemblies
Each strain FASTA becomes an assembly whose name matches an entry in
assemblyNames on the track:
for strain in K12 Sakai CFT073 NCTC86 IAI39; do
bgzip -f "$strain.fa"
samtools faidx "$strain.fa.gz" # writes the .fai and .gzi JBrowse needs
jbrowse add-assembly "$strain.fa.gz" --name "$strain" --load copy
doneThe same step for one strain, as a config or in JBrowse Desktop:
Goes in the assemblies array of config.json. See Assemblies.
{
"name": "K12",
"uri": "K12.fa.gz"
}jbrowse add-assembly K12.fa.gz \
--name K12 \
--load copyIn JBrowse Desktop, Open new genome on the start screen (or File → Open genome... in a session), then Open from a URL and paste, one per line:
K12.fa.gz
K12.fa.gz.fai
K12.fa.gz.gziJBrowse reads the format off the file name. Then fill in:
- Genome name:
K12
K12.fa.gz is relative to a config.json. Replace it with its URL or its path on this computer.
An assembly whose refNames still have the PanSN prefix draws empty. The assemblies configuration guide has more.
Loading the PAF with MultiGenomePAFAdapter
List every assembly the file covers in assemblyNames; the adapter keeps only
the records whose PanSN prefixes match the pair of rows each band joins:
Goes in the tracks array of config.json. See Tracks.
{
"type": "SyntenyTrack",
"trackId": "ecoli_ava",
"name": "E. coli pangenome (all-vs-all PAF)",
"assemblyNames": ["K12", "Sakai", "CFT073", "NCTC86", "IAI39"],
"adapter": {
"type": "MultiGenomePAFAdapter",
"uri": "all_vs_all.paf",
"assemblyNames": ["K12", "Sakai", "CFT073", "NCTC86", "IAI39"]
}
}jbrowse add-track all_vs_all.paf \
--trackId ecoli_ava \
--name "E. coli pangenome (all-vs-all PAF)" \
--assemblyNames K12,Sakai,CFT073,NCTC86,IAI39 \
--adapterType MultiGenomePAFAdapter \
--load copyIn JBrowse Desktop, or in any running JBrowse Web session, open a view on this track’s assembly, then File → Open track..., choose Add track from pasted JSON, and paste:
{
"type": "SyntenyTrack",
"trackId": "ecoli_ava",
"name": "E. coli pangenome (all-vs-all PAF)",
"assemblyNames": ["K12", "Sakai", "CFT073", "NCTC86", "IAI39"],
"adapter": {
"type": "MultiGenomePAFAdapter",
"uri": "all_vs_all.paf",
"assemblyNames": ["K12", "Sakai", "CFT073", "NCTC86", "IAI39"]
}
}all_vs_all.paf is relative to a config.json. Replace it with its URL or its path on this computer.
That shows it in the linear view. For the synteny view, open Add → Linear synteny view, pick the track under Quick start, and click Launch.
MultiGenomePAFAdapter has to be named; from a .paf extension JBrowse guesses
the pairwise PAFAdapter, which reads only two assembly names. An assembly name
that differs from its PanSN sample prefix needs a map.1
Large files: index with make-pif
MultiGenomePAFAdapter reads the whole PAF into memory. For many samples, index
it with jbrowse make-pif and switch to MultiGenomeIndexedPAFAdapter,
fetching only the region in view:
# produces all_vs_all.pif.gz and all_vs_all.pif.gz.tbi
jbrowse make-pif all_vs_all.pafmake-pif finishes by printing the add-track command for the samples it
found:
jbrowse add-track all_vs_all.pif.gz --adapterType MultiGenomeIndexedPAFAdapter \
-a CFT073,IAI39,K12,NCTC86,Sakai --load copyOnly the adapter block differs from the un-indexed version:
Goes in the tracks array of config.json. See Tracks.
{
"type": "SyntenyTrack",
"trackId": "ecoli_ava_indexed",
"name": "E. coli all-vs-all (indexed)",
"assemblyNames": ["K12", "Sakai", "CFT073", "NCTC86", "IAI39"],
"adapter": {
"type": "MultiGenomeIndexedPAFAdapter",
"uri": "all_vs_all.pif.gz",
"assemblyNames": ["K12", "Sakai", "CFT073", "NCTC86", "IAI39"]
}
}jbrowse add-track all_vs_all.pif.gz \
--trackId ecoli_ava_indexed \
--name "E. coli all-vs-all (indexed)" \
--assemblyNames K12,Sakai,CFT073,NCTC86,IAI39 \
--adapterType MultiGenomeIndexedPAFAdapter \
--load copyIn JBrowse Desktop, or in any running JBrowse Web session, open a view on this track’s assembly, then File → Open track..., choose Add track from pasted JSON, and paste:
{
"type": "SyntenyTrack",
"trackId": "ecoli_ava_indexed",
"name": "E. coli all-vs-all (indexed)",
"assemblyNames": ["K12", "Sakai", "CFT073", "NCTC86", "IAI39"],
"adapter": {
"type": "MultiGenomeIndexedPAFAdapter",
"uri": "all_vs_all.pif.gz",
"assemblyNames": ["K12", "Sakai", "CFT073", "NCTC86", "IAI39"]
}
}all_vs_all.pif.gz is relative to a config.json. Replace it with its URL or its path on this computer.
That shows it in the linear view. For the synteny view, open Add → Linear synteny view, pick the track under Quick start, and click Launch.
make-pif also writes a coarse copy of the alignments for zoomed-out
views.2
Stacking the five strains
Stacking the strains from the import form
- Add → Linear synteny view opens the form in Quick start.
- Choose
ecoli_ava. Its five assemblies each become a row. - Click Launch.
Manual mode builds the stack by hand: Add row per strain, and the connector button between each pair to pick its track.
Stacking the strains in a defaultSession
A defaultSession holding a LinearSyntenyView opens the stack on load:
Goes at the top level of config.json, replacing any defaultSession there. See Default session.
{
"defaultSession": {
"name": "E. coli 5-strain pangenome",
"views": [
{
"type": "LinearSyntenyView",
"views": [
{ "assembly": "K12" },
{ "assembly": "Sakai" },
{ "assembly": "CFT073" },
{ "assembly": "NCTC86" },
{ "assembly": "IAI39" }
],
"tracks": [["ecoli_ava"], ["ecoli_ava"], ["ecoli_ava"], ["ecoli_ava"]],
"minAlignmentLength": 10000,
"collapseEmptyRows": true
}
]
}
}jbrowse set-default-session --session - << 'EOF'
{
"name": "E. coli 5-strain pangenome",
"views": [
{
"type": "LinearSyntenyView",
"views": [
{ "assembly": "K12" },
{ "assembly": "Sakai" },
{ "assembly": "CFT073" },
{ "assembly": "NCTC86" },
{ "assembly": "IAI39" }
],
"tracks": [["ecoli_ava"], ["ecoli_ava"], ["ecoli_ava"], ["ecoli_ava"]],
"minAlignmentLength": 10000,
"collapseEmptyRows": true
}
]
}
EOFtracksis one entry per band, so five rows take four:tracks[0]connects rows 0-1,tracks[1]rows 1-2, and so oncollapseEmptyRowsgives each trackless row a bare scale bar
Row order is a free choice with an all-vs-all PAF, since the file aligns every pair.
The gaps between ribbons mark where the strains differ.
Adding gene tracks to see what a gap holds
Each strain's GFF becomes a gene track attached to that strain's assembly:
for strain in K12 Sakai CFT073 NCTC86 IAI39; do
jbrowse sort-gff "$strain.gff" | bgzip > "$strain.gff.gz"
tabix "$strain.gff.gz"
# -a "$strain": attach the track to that one assembly, so it stays with that row
jbrowse add-track "$strain.gff.gz" -a "$strain" --name "$strain genes" --load copy
doneNavigate the K-12 row to chr:1,026,000-1,126,000 and the Sakai row to
chr:1,205,000-1,305,000. The gap right of the ribbon holds stx2A and
stx2B, the Shiga-toxin subunits, with no alignment to K-12.
One strain against all the others
In a plain linear genome view, with no second row to name a target assembly, the
ecoli_ava track draws the strain you're viewing against every other strain in
the file. Clicking an alignment can launch a synteny view against its mate, the
strain it aligns to.
Three track-menu items set the pileup up for reading strain by strain:
- Group by... → Mate assembly gives each strain its own lane. Untick Show... → Collapse groups to one row to stack every lane, or expand one from its label.
- Group by... → Hide self-alignment lane drops the lane for the strain you're viewing. The figures below have it ticked.
- Show... → Show coverage adds a histogram of how many other strains cover each base.
The figure below shows K-12's phenylacetate (paa) operon, where three strains stop at the shaded edge and NCTC86 runs through, with a track of the pangenome graph's segments, built with minigraph in the E. coli pangenome tutorial, above the lanes.
At whole-chromosome zoom, the same lanes also fit on the K-12 row of the five-strain stack:
- List IAI39 second in the session's
views. - Add
ecoli_avato the K-12 row from that row's track selector. - Pick Strand from the palette button, so an inversion is blue in both halves.
Percent identity per strain
The ecoli_ava track can also draw each alignment as a line at its identity,
one row per strain, the percent identity plot
PipMaker drew for a pair of genomes. In
the track menu, Display types → Marks draws it with nothing to configure.
The config below is that plot written out, with four details to read:
identitycomes from minimap2'sdedivergence tagrowsgives eachmate.assemblyNameone rowfilterdrops K-12's alignments to itself, which otherwise take a row of their ownscales.y.zeroisfalse, so the autoscaled axis spans the identities alone and does not reach 0
Goes in the tracks array of config.json. See Tracks.
{
"type": "SyntenyTrack",
"trackId": "ecoli_ava_identity",
"name": "K-12 against each strain, identity",
"assemblyNames": ["K12", "Sakai", "CFT073", "NCTC86", "IAI39"],
"adapter": {
"type": "MultiGenomePAFAdapter",
"uri": "https://jbrowse.org/demos/ecoli_pangenome/all_vs_all.paf.gz",
"assemblyNames": ["K12", "Sakai", "CFT073", "NCTC86", "IAI39"]
},
"displays": [
{
"type": "LinearMarkDisplay",
"displayId": "ecoli_ava_identity-LinearMarkDisplay",
"rows": "mate.assemblyName",
"filter": ["jexl:feature.mate.assemblyName != 'K12'"],
"scales": { "y": { "zero": false } },
"marks": [{ "mark": "rule", "encoding": { "y": "identity", "size": 2 } }]
}
]
}jbrowse add-track-json '{
"type": "SyntenyTrack",
"trackId": "ecoli_ava_identity",
"name": "K-12 against each strain, identity",
"assemblyNames": ["K12", "Sakai", "CFT073", "NCTC86", "IAI39"],
"adapter": {
"type": "MultiGenomePAFAdapter",
"uri": "https://jbrowse.org/demos/ecoli_pangenome/all_vs_all.paf.gz",
"assemblyNames": ["K12", "Sakai", "CFT073", "NCTC86", "IAI39"]
},
"displays": [
{
"type": "LinearMarkDisplay",
"displayId": "ecoli_ava_identity-LinearMarkDisplay",
"rows": "mate.assemblyName",
"filter": ["jexl:feature.mate.assemblyName != '\''K12'\''"],
"scales": { "y": { "zero": false } },
"marks": [{ "mark": "rule", "encoding": { "y": "identity", "size": 2 } }]
}
]
}'In JBrowse Desktop, or in any running JBrowse Web session, open a view on this track’s assembly, then File → Open track..., choose Add track from pasted JSON, and paste:
{
"type": "SyntenyTrack",
"trackId": "ecoli_ava_identity",
"name": "K-12 against each strain, identity",
"assemblyNames": ["K12", "Sakai", "CFT073", "NCTC86", "IAI39"],
"adapter": {
"type": "MultiGenomePAFAdapter",
"uri": "https://jbrowse.org/demos/ecoli_pangenome/all_vs_all.paf.gz",
"assemblyNames": ["K12", "Sakai", "CFT073", "NCTC86", "IAI39"]
},
"displays": [
{
"type": "LinearMarkDisplay",
"displayId": "ecoli_ava_identity-LinearMarkDisplay",
"rows": "mate.assemblyName",
"filter": ["jexl:feature.mate.assemblyName != 'K12'"],
"scales": { "y": { "zero": false } },
"marks": [{ "mark": "rule", "encoding": { "y": "identity", "size": 2 } }]
}
]
}That shows it in the linear view. For the synteny view, open Add → Linear synteny view, pick the track under Quick start, and click Launch.
Open this track in JBrowse Web ↗
Open this track in JBrowse Desktop ↗ (JBrowse Desktop 5.0+)
Redrawing each lane in its own strain's coordinates
On the K-12 axis, a strain with no alignment to the backbone is a white gap. Display types → Multi-way synteny display redraws each lane in the coordinates of the strain it shows, with each PAF record as one ribbon; the ortholog-table tutorial has a figure of the display.
Launching a stacked view at one locus
Drag-select a region on the scale bar and pick Launch → Linear synteny view. The dialog lists every assembly aligning to that region as a panel, top to bottom, with arrows to reorder them and a checkbox to drop one; Replace current view or Open in new view launches the stack. Ribbons draw between neighbouring rows only.
Right-clicking a single alignment offers three routes under Launch:
- Linear synteny view with Sakai (or whichever strain the alignment names) opens that one pair.
- Linear synteny view, all assemblies here opens the same multi-strain dialog.
- Open Sakai at the matching region opens a linear genome view of Sakai at that region.
A launched view is a few kilobases wide, and the CIGAR minimap2 -c wrote draws
each insertion and deletion where it falls. CIGAR indels in the settings
menu switches between colored, transparent and none.
Checking the Sakai gap against the PAF
Print the alignments between one strain and another near a gap, here Sakai
against K-12 near the stx2 island. -X emits each pair once in either
direction, so the coordinates come from whichever column the strain landed in.
Set s and m to the two strain prefixes and lo and hi to the window:
awk -F'\t' -v OFS='\t' -v s=Sakai -v m=K12 '
$1 ~ "^" s "#" && $6 ~ "^" m "#" { print $3, $4; next }
$1 ~ "^" m "#" && $6 ~ "^" s "#" { print $8, $9 }
' all_vs_all.paf | sort -n | awk -v lo=1200000 -v hi=1300000 '$1 < hi && $2 > lo'1207288 1207877
1210882 1246166
1251954 1252260
1274685 1275548The second line is the shared backbone the ribbon draws. The third is short, and the fourth starts at 1,274,685, so stx2A and stx2B fall in a stretch with no K-12 counterpart.
Reproduce it end to end
build_ecoli_pangenome_synteny.sh
runs everything on this page, download and preparation included, and needs the
tools under Prerequisites:
- Download the five RefSeq assemblies with their annotation, and keep each
one's chromosome under the name
chr, so the plasmids drop out and every strain's row reads the same name. - Concatenate the strains under PanSN names and align them with
minimap2 -X, so each pair is aligned once and no strain aligns to itself. - Keep each GFF's chromosome features under the same
chrname, so the genes load on that strain's row. - Write the config with the five assemblies, the gene tracks, the all-vs-all track and the five-row session.
curl -fO https://raw.githubusercontent.com/GMOD/jbrowse-components/main/scripts/build_ecoli_pangenome_synteny.sh
bash build_ecoli_pangenome_synteny.sh # builds ./ecoli_pangenome_build/jbrowse2
npx --yes serve ecoli_pangenome_build/jbrowse2 # then open the printed URLSee also
- Pangenome (pggb)
- Synteny from gene symbols (44 E. coli genomes)
- Pangenome (Minigraph-Cactus)
- Synteny visualization (pairwise minimap2)
- Synteny from an ortholog table (grape, peach, cacao)
- Dotplot view
- Synteny track
- MultiGenomePAFAdapter
- MultiGenomeIndexedPAFAdapter
- PIF (Pairwise Indexed Format)
- jbrowse-anywidget
- JBrowseR
Notes
-
assemblyNameToPanSNmaps such a name, e.g.{ "Ecoli_K12": "K12" }. A haplotype-resolved pangenome can map each haplotype to a separate assembly with asample#haplotypeprefix; see PanSN depth. ↩ -
--coarsesets how many bp a coarse row may stray from the real alignment; with a larger value, raisecoarseBpPerPxThresholdon the adapter too.--csiswaps the TBI index for sequences over ~512 Mb. ↩
Feedback on this tutorial is welcome: contact us.