Skip to main content

Synteny from gene symbols (47 E. coli genomes)

TL;DR: we look at one K-12 operon across forty-six other E. coli and Shigella genomes at once, without aligning any of them. RefSeq's bacterial pipeline gives an orthologous gene the same symbol in every strain it names, so the ortholog table is a join on the gene name over the GFF3 files, and each genome becomes a lane under the K-12 view carrying its own gene models. The join sees the core genome and nothing else: at a cluster that differs between strains the lanes draw their own genes and no ribbons, which is where the page ends.

Prerequisites

Where the data comes from

Forty-seven RefSeq assemblies, each fetched by accession with the datasets CLI: the classic reference strains across phylogroups A, B1, B2, D and E, four Shigella, and complete genomes picked by striding a datasets summary listing. The build script pins every accession; the anchor and the four Shigella are:

A join instead of an alignment

The all-vs-all page describes five E. coli strains to each other by aligning them, which is one minimap2 run per pair and what holds that demo at five. RefSeq's prokaryotic annotation pipeline names a gene by its ortholog (atpA is atpA in every strain that carries it), so the same .blocks table can be filled by matching symbols across the GFF3 files, the route the primate page takes for eight apes. The download is an annotation and a sequence report per genome:

# one RefSeq accession per line; the first is the genome the rows are anchored on
datasets download genome accession --inputfile accessions.txt \
  --include gff3,seq-report --filename genomes.zip
unzip genomes.zip

Each genome's chromosome is the longest sequence in its report, and the plasmids are dropped: a lane follows one contig at a time. The GFF3 filtered to that sequence, sorted, bgzipped and tabix-indexed as in the web quickstart is the genome's gene track.

PGAP writes a gene's locus tag into its Name when it has no symbol for it, so the join has to be told what an unnamed gene looks like, or every hypothetical protein would look named and match nothing:

curl -fO https://raw.githubusercontent.com/GMOD/jbrowse-components/main/scripts/symbols_to_blocks.py
# --anchor names the genome whose genes are the rows; --unnamed is the
# locus-tag shape PGAP falls back to, which joins nothing
python3 symbols_to_blocks.py --anchor MG1655 -o ecoli.blocks --unnamed '_RS[0-9]+$' \
  MG1655=MG1655.gff.gz Sakai=Sakai.gff.gz CFT073=CFT073.gff.gz Sflex301=Sflex301.gff.gz

The helper reports how much of each column it filled, and that number is the screen a strain has to pass. Older PGAP runs named genes by locus tag alone, so a genome can be complete, current and join nothing; of a hundred accessions tried for this page, a quarter of those that still downloaded were like that, and two were phage genomes. Counting the named genes in an annotation says which before anything is built:

gzip -dc strain.gff.gz | awk -F'\t' '$3 == "gene" && $9 ~ /;gene=/' | wc -l

The ortholog track

One SyntenyTrack names all forty-seven assemblies, each of which is a ChromSizesAdapter over its chromosome's length, since the lanes never read sequence. blockAssemblies and bedLocations are positional against the table's columns, in the order the helper printed; the list below is cut to the first four for the page. A lane needs its header's height to draw, so the display's height is set for the whole stack rather than left at the default:

{
  "type": "SyntenyTrack",
  "trackId": "ecoli_orthologs",
  "name": "E. coli orthologs by gene symbol (47 genomes, RefSeq)",
  "assemblyNames": ["MG1655", "Sakai", "CFT073", "Sflex301"],
  "adapter": {
    "type": "MCScanBlocksAdapter",
    "mcscanBlocksLocation": { "uri": "ecoli.blocks" },
    "blockAssemblies": ["MG1655", "Sakai", "CFT073", "Sflex301"],
    "bedLocations": [
      { "uri": "MG1655.bed" },
      { "uri": "Sakai.bed" },
      { "uri": "CFT073.bed" },
      { "uri": "Sflex301.bed" }
    ],
    "assemblyNames": ["MG1655", "Sakai", "CFT073", "Sflex301"]
  },
  "displays": [
    {
      "type": "MultiWaySyntenyDisplay",
      "displayId": "ecoli_orthologs-MultiWaySyntenyDisplay",
      "color": "jexl:feature.name ? randomColor(feature.name) : '#b0b0b0'",
      "height": 1100
    }
  ]
}
jbrowse add-track-json '{
  "type": "SyntenyTrack",
  "trackId": "ecoli_orthologs",
  "name": "E. coli orthologs by gene symbol (47 genomes, RefSeq)",
  "assemblyNames": ["MG1655", "Sakai", "CFT073", "Sflex301"],
  "adapter": {
    "type": "MCScanBlocksAdapter",
    "mcscanBlocksLocation": { "uri": "ecoli.blocks" },
    "blockAssemblies": ["MG1655", "Sakai", "CFT073", "Sflex301"],
    "bedLocations": [
      { "uri": "MG1655.bed" },
      { "uri": "Sakai.bed" },
      { "uri": "CFT073.bed" },
      { "uri": "Sflex301.bed" }
    ],
    "assemblyNames": ["MG1655", "Sakai", "CFT073", "Sflex301"]
  },
  "displays": [
    {
      "type": "MultiWaySyntenyDisplay",
      "displayId": "ecoli_orthologs-MultiWaySyntenyDisplay",
      "color": "jexl:feature.name ? randomColor(feature.name) : '\''#b0b0b0'\''",
      "height": 1100
    }
  ]
}'

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_orthologs",
  "name": "E. coli orthologs by gene symbol (47 genomes, RefSeq)",
  "assemblyNames": ["MG1655", "Sakai", "CFT073", "Sflex301"],
  "adapter": {
    "type": "MCScanBlocksAdapter",
    "mcscanBlocksLocation": { "uri": "ecoli.blocks" },
    "blockAssemblies": ["MG1655", "Sakai", "CFT073", "Sflex301"],
    "bedLocations": [
      { "uri": "MG1655.bed" },
      { "uri": "Sakai.bed" },
      { "uri": "CFT073.bed" },
      { "uri": "Sflex301.bed" }
    ],
    "assemblyNames": ["MG1655", "Sakai", "CFT073", "Sflex301"]
  },
  "displays": [
    {
      "type": "MultiWaySyntenyDisplay",
      "displayId": "ecoli_orthologs-MultiWaySyntenyDisplay",
      "color": "jexl:feature.name ? randomColor(feature.name) : '#b0b0b0'",
      "height": 1100
    }
  ]
}

One operon, forty-seven genomes

Opened on K-12 at the atp operon, the track draws a lane per genome under the K-12 axis. Every gene is colored by its symbol, so a conserved gene is one color running down the whole stack, and a lane's header names its chromosome, where it is looking and [rev] where the strain's chromosome reads the other way. Lanes stack densest first, so the genomes placing the most of the window sit at the top and the reduced Shigella genomes fall toward the bottom without anything naming them.

{
  "defaultSession": {
    "name": "The atp operon across 47 genomes",
    "views": [
      {
        "type": "LinearGenomeView",
        "assembly": "MG1655",
        "loc": "NC_000913.3:3,910,000-3,925,000",
        "tracks": [
          {
            "trackId": "ecoli_orthologs",
            "type": "MultiWaySyntenyDisplay",
            "height": 1100
          }
        ]
      }
    ]
  }
}
jbrowse set-default-session --session - << 'EOF'
{
  "name": "The atp operon across 47 genomes",
  "views": [
    {
      "type": "LinearGenomeView",
      "assembly": "MG1655",
      "loc": "NC_000913.3:3,910,000-3,925,000",
      "tracks": [
        {
          "trackId": "ecoli_orthologs",
          "type": "MultiWaySyntenyDisplay",
          "height": 1100
        }
      ]
    }
  ]
}
EOF

Open this session in JBrowse

The atp operon on K-12 over forty-six E. coli and Shigella lanes from one gene-symbol ortholog track, each lane drawing its own RefSeq gene models. Every gene's color runs the full stack; the lanes reading the operon reversed are the ones whose chromosome was deposited the other way round. Open this view in JBrowse ↗ (large dataset — this session takes a while to load)

The steps behind the figure above, run against your own data.

  1. Open your genome: on the JBrowse Desktop start screen click Open new genome (or Show all available genomes to pick a hosted one).This figure uses MG1655.
  2. Add your own track: File → Open track..., then paste a URL or choose a local file.This figure uses the “ecoli_orthologs” track.
  3. Drag the bar at the bottom edge of the track to resize it.
  4. Type your region of interest into the location box and press Enter.This figure is at NC_000913.3:3,910,000-3,925,000.

Open this view in JBrowse Desktop ↗

Opens JBrowse Desktop (JBrowse Desktop 5.0+) at this view and saves it as a reopenable session — swap in your own files afterwards.

Nothing happens? Paste this link into Desktop's Open JBrowse Web link... (start screen, or File → Session):

https://jbrowse.org/code/jb2/main/?config=https%3A%2F%2Fjbrowse.org%2Fdemos%2Fecoli_orthologs%2Fconfig.json&session=spec-%7B%22views%22%3A%5B%7B%22type%22%3A%22LinearGenomeView%22%2C%22assembly%22%3A%22MG1655%22%2C%22loc%22%3A%22NC_000913.3%3A3%2C910%2C000-3%2C925%2C000%22%2C%22tracks%22%3A%5B%7B%22trackId%22%3A%22ecoli_orthologs%22%2C%22type%22%3A%22MultiWaySyntenyDisplay%22%2C%22height%22%3A1100%7D%5D%7D%5D%7D&sessionName=JBrowse+docs%3A+multiway_synteny%2Fecoli_symbol_atp_operon

This session spec pastes after &session=spec- on any JBrowse Web instance.

{
  "views": [
    {
      "type": "LinearGenomeView",
      "assembly": "MG1655",
      "loc": "NC_000913.3:3,910,000-3,925,000",
      "tracks": [
        {
          "trackId": "ecoli_orthologs",
          "type": "MultiWaySyntenyDisplay",
          "height": 1100
        }
      ]
    }
  ]
}

Loaded against config: https://jbrowse.org/demos/ecoli_orthologs/config.json

The same view with jbrowse-anywidget. Point the adapter at your own file.

from jbrowse_anywidget import LinearGenomeView, make_assembly

view = LinearGenomeView(
    assembly=make_assembly("MG1655", "https://your-server/your-genome.fa.gz"),
    location="NC_000913.3:3,910,000-3,925,000",
)

view.add_track({
    "type": "FeatureTrack",
    "trackId": "my_track",
    "name": "My track",
    "assemblyNames": ["MG1655"],
    "adapter": {"type": "BamAdapter", "uri": "https://your-server/your-file"},
})

view  # display the widget

Hand this to a coding agent. It rebuilds the figure above headlessly, and the session file is then the thing to edit — swap an adapter uri for your own file and rerun.

cat > session.json <<'JSON'
{
  "views": [
    {
      "type": "LinearGenomeView",
      "assembly": "MG1655",
      "loc": "NC_000913.3:3,910,000-3,925,000",
      "tracks": [
        {
          "trackId": "ecoli_orthologs",
          "type": "MultiWaySyntenyDisplay",
          "height": 1100
        }
      ]
    }
  ]
}
JSON

npx @jbrowse/capture --instance https://jbrowse.org/code/jb2/main/ \
  --config https://jbrowse.org/demos/ecoli_orthologs/config.json \
  --session session.json -o figure.png

Using JBrowse with AI agents covers the rest: hosted genomes to point it at, and why waiting for the render is the part that goes wrong.

Where the join stops

The O-antigen cluster between galF and gnd is the locus that differs most between strains: each serotype carries its own set of sugar pathway genes, and those share no symbol across serotypes. The lanes that keep the whole cluster sort to the top, and they are the K-12 derivatives, whose cluster is K-12's gene for gene. Below them the flanking genes chain down every lane, and between the flanks each lane draws whatever its own annotation holds with no ribbon to it. That is the accessory genome as the join sees it.

The O-antigen cluster on K-12 over the same forty-six lanes. The K-12 derivatives at the top of the stack match the cluster gene for gene; in every lane below, the flanking galF and gnd chains run through and the cluster between them is that strain's own, joined to nothing. Open this view in JBrowse ↗ (large dataset — this session takes a while to load)

The steps behind the figure above, run against your own data.

  1. Open your genome: on the JBrowse Desktop start screen click Open new genome (or Show all available genomes to pick a hosted one).This figure uses MG1655.
  2. Add your own track: File → Open track..., then paste a URL or choose a local file.This figure uses the “ecoli_orthologs” track.
  3. Drag the bar at the bottom edge of the track to resize it.
  4. Type your region of interest into the location box and press Enter.This figure is at NC_000913.3:2,095,000-2,115,000.

Open this view in JBrowse Desktop ↗

Opens JBrowse Desktop (JBrowse Desktop 5.0+) at this view and saves it as a reopenable session — swap in your own files afterwards.

Nothing happens? Paste this link into Desktop's Open JBrowse Web link... (start screen, or File → Session):

https://jbrowse.org/code/jb2/main/?config=https%3A%2F%2Fjbrowse.org%2Fdemos%2Fecoli_orthologs%2Fconfig.json&session=spec-%7B%22views%22%3A%5B%7B%22type%22%3A%22LinearGenomeView%22%2C%22assembly%22%3A%22MG1655%22%2C%22loc%22%3A%22NC_000913.3%3A2%2C095%2C000-2%2C115%2C000%22%2C%22tracks%22%3A%5B%7B%22trackId%22%3A%22ecoli_orthologs%22%2C%22type%22%3A%22MultiWaySyntenyDisplay%22%2C%22height%22%3A1100%7D%5D%7D%5D%7D&sessionName=JBrowse+docs%3A+multiway_synteny%2Fecoli_symbol_oantigen

This session spec pastes after &session=spec- on any JBrowse Web instance.

{
  "views": [
    {
      "type": "LinearGenomeView",
      "assembly": "MG1655",
      "loc": "NC_000913.3:2,095,000-2,115,000",
      "tracks": [
        {
          "trackId": "ecoli_orthologs",
          "type": "MultiWaySyntenyDisplay",
          "height": 1100
        }
      ]
    }
  ]
}

Loaded against config: https://jbrowse.org/demos/ecoli_orthologs/config.json

The same view with jbrowse-anywidget. Point the adapter at your own file.

from jbrowse_anywidget import LinearGenomeView, make_assembly

view = LinearGenomeView(
    assembly=make_assembly("MG1655", "https://your-server/your-genome.fa.gz"),
    location="NC_000913.3:2,095,000-2,115,000",
)

view.add_track({
    "type": "FeatureTrack",
    "trackId": "my_track",
    "name": "My track",
    "assemblyNames": ["MG1655"],
    "adapter": {"type": "BamAdapter", "uri": "https://your-server/your-file"},
})

view  # display the widget

Hand this to a coding agent. It rebuilds the figure above headlessly, and the session file is then the thing to edit — swap an adapter uri for your own file and rerun.

cat > session.json <<'JSON'
{
  "views": [
    {
      "type": "LinearGenomeView",
      "assembly": "MG1655",
      "loc": "NC_000913.3:2,095,000-2,115,000",
      "tracks": [
        {
          "trackId": "ecoli_orthologs",
          "type": "MultiWaySyntenyDisplay",
          "height": 1100
        }
      ]
    }
  ]
}
JSON

npx @jbrowse/capture --instance https://jbrowse.org/code/jb2/main/ \
  --config https://jbrowse.org/demos/ecoli_orthologs/config.json \
  --session session.json -o figure.png

Using JBrowse with AI agents covers the rest: hosted genomes to point it at, and why waiting for the render is the part that goes wrong.

For the variable loci the table wants a homology call across the proteomes, an OrthoFinder run, or the all-vs-all alignment that draws the same locus base by base for five strains.

Reproduce it end to end

The script fetches the forty-seven annotations, builds the table and writes the config; see Prerequisites.

curl -fO https://raw.githubusercontent.com/GMOD/jbrowse-components/main/scripts/build_ecoli_orthologs.sh
bash build_ecoli_orthologs.sh

See also

References