Skip to main content

Reviewing a whole SV callset

JBrowse v5 beta

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.8. Desktop beta builds are coming soon.

jb2export batch renders every junction in a somatic callset as a breakpoint split view, so triage becomes a directory of images. The matched normal, rendered the same way, is the control.

Prerequisites

  • @jbrowse/img, which puts jb2export on your PATH
  • nothing to download: the callset, the tumor reads and the matched normal are all hosted
  • a JBrowse, for the last section only (Web or Desktop); the renders above need none
npm install -g @jbrowse/img

Where the data comes from

COLO829's somatic SV callset is the ONT open-data release's own wf-somatic-variation run (Valle-Inclán et al. 2022), rehosted alongside the cancer SV demo.

COLO829 and its matched normal

COLO829 is a melanoma cell line with a matched normal, COLO829BL, and a community reference for somatic structural-variant calling. The multi-hop tutorial follows one event in this callset all the way down; this page renders every junction at a glance.

The contact sheet

A junction joins two loci, which are the two panels of a breakpoint split view, so a callset renders straight into a review queue:

curl -fO https://jbrowse.org/demos/cancer_sv/COLO829.somatic-sv.vcf.gz
jb2export batch --vcf COLO829.somatic-sv.vcf.gz \
  --config https://jbrowse.org/demos/cancer_sv/config.json --assembly hg38 \
  --track COLO829_tumor_ont height:240 \
  --outDir tumor --flank 600 --width 1100
Warning: skipped 35 record(s), e.g. line 271: names no second locus (INS)
[########################] 100% 100/100
wrote 100/100 images to tumor

Insertions name one locus, so there is no second panel and the warning counts them out. Each breakend pair is written twice and collapses to one. The result is the same 100 junctions, in the same order, that sv_multihop.py chains reports in the multi-hop tutorial, because neither parses the ALT bracket by hand, which goes wrong four ways, none of them raising an error:

  • the replacement string may carry inserted sequence either side of the bracket (GTGATGGATTCA[CHR12:72273112[)
  • callers upper-case the mate contig, and CHR12 is not a region hg38 has
  • END= matches inside CIEND=, and the first hit wins
  • the two records of one breakend pair name the same translocation twice

One image per row, named 1_chr1_33053494-chr6_2919922_gridss12o.png: index first so the directory sorts in callset order, coordinates next, the caller's ID last. A file with no ID column falls back to junction_<n>.

--flank frames the panel, since a breakend is one base. --dryRun prints the file and loci of every row and renders nothing, and --limit 20 renders the first few, to check the framing before the whole callset.

For a long run:

  • --resume skips a row whose image is already in --outDir
  • --manifest writes manifest.tsv beside the images: one row per junction with its file, both loci, its name, and whether it rendered
  • --passOnly drops records the caller filtered out. --limit takes the first N in file order, so on an unfiltered callset the two go together

The reads stream from the hosted CRAM, the module graph loads once for the whole callset, and a --config URL or --hub is fetched once. A row that cannot be rendered is reported and the run continues.

A connector drawn dashed means the read has a segment at a locus the frame does not show. These reads also visit chr10, so this junction wants a third panel, and the control belongs beside it: one render per sample, the same --loc list and --width.

The three loci of COLO829's der(3), chr3 then chr10 then chr12, at the same width in every panel. The tumor nanopore reads carry a solid curve at every breakend and the matched normal carries none. On the right, the same three loci as one reconstructed contig.

featureHeight:super-compact draws reads at 1 px apiece, which keeps six pileups on one screen.

A curve marks two loci as joined; a contig shows the order and orientation. The multi-hop tutorial builds that contig from these reads, and rendering it is another jb2export run with a different --assembly.

The same export over the normal

One directory per track:

jb2export batch --vcf COLO829.somatic-sv.vcf.gz \
  --config https://jbrowse.org/demos/cancer_sv/config.json --assembly hg38 \
  --track COLO829BL_normal_ont height:240 \
  --outDir normal --flank 600 --width 1100

Side by side, the somatic calls are the ones with curves in tumor/ and none in normal/.

Reading the sheet

What each picture shows:

  • a fan of curves at both breakends is the junction as the reads describe it
  • nothing connecting the panels means the reads do not support the caller's coordinates, which is either a false call or a breakpoint placed far enough off that --flank missed it. Re-render that row wider
  • curves in the normal too means germline, not somatic
  • a dense fan in a region of ragged coverage is usually a repeat. The connectors are drawn from what the aligner reported, so a read mismapped into a repeat contributes a confident-looking curve

Opening a call in the browser

Take the coordinates from an image's filename, open the SV inspector on the same VCF, and click through to the breakpoint split view, with the gene track and read details attached.

For a junction that is one hop of something larger, the alignments track menu's Reconstruct derivative allele... groups the reads in view by the route their split alignments describe. COLO829's der(3) is three junctions across three chromosomes, and the multi-hop tutorial follows it the rest of the way.

The junction the figures above are of, taken the rest of the way in the browser: the tumor track menu, the routes the reads describe with the count behind each one, and Breakpoint split view replacing the window with one panel per segment of the route. Open the session this video starts in ↗

Open the session the video starts in ↗, then take the steps it shows:

  1. Let the reads reconstruct the derivative allele
  2. 3 → 10 → 12 (inverted) → 3 (inverted)
  3. Draw as: Breakpoint split view

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 hg38.
  2. Add your own track: File → Open track..., then paste a URL or choose a local file.This figure uses the “ncbi_refseq_hg38” track.
  3. Track menu → Gene glyph → Representative transcript“ncbi_refseq_hg38”.
  4. Drag the bar at the bottom edge of the track to resize it.“ncbi_refseq_hg38”.
  5. Add your own track: File → Open track..., then paste a URL or choose a local file.This figure uses the “COLO829_tumor_ont” track.
  6. Click Force load in the track's "Zoom in to see features or force load" message.Loads the region even past the byte-size limit, which can be slow.“COLO829_tumor_ont”.
  7. Track menu → Show... → Show soft clipping (checked)Reveals clipped bases — the signal that a read spans a breakpoint.“COLO829_tumor_ont”.
  8. Drag the bar at the bottom edge of the track to resize it.“COLO829_tumor_ont”.
  9. Track menu → Read height → Super-compact“COLO829_tumor_ont”.
  10. Type your region of interest into the location box and press Enter.This figure is at chr3:25,357,600-25,361,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%2Fcancer_sv%2Fconfig.json&session=spec-%7B%22views%22%3A%5B%7B%22type%22%3A%22LinearGenomeView%22%2C%22assembly%22%3A%22hg38%22%2C%22loc%22%3A%22chr3%3A25%2C357%2C600-25%2C361%2C000%22%2C%22tracks%22%3A%5B%7B%22trackId%22%3A%22ncbi_refseq_hg38%22%2C%22geneGlyphMode%22%3A%22longestCoding%22%2C%22height%22%3A70%7D%2C%7B%22trackId%22%3A%22COLO829_tumor_ont%22%2C%22type%22%3A%22LinearAlignmentsDisplay%22%2C%22forceLoad%22%3Atrue%2C%22showSoftClipping%22%3Atrue%2C%22height%22%3A130%2C%22featureHeight%22%3A1%7D%5D%7D%5D%7D&sessionName=JBrowse+docs%3A+sv%2Fderivative_allele_route

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

{
  "views": [
    {
      "type": "LinearGenomeView",
      "assembly": "hg38",
      "loc": "chr3:25,357,600-25,361,000",
      "tracks": [
        {
          "trackId": "ncbi_refseq_hg38",
          "geneGlyphMode": "longestCoding",
          "height": 70
        },
        {
          "trackId": "COLO829_tumor_ont",
          "type": "LinearAlignmentsDisplay",
          "forceLoad": true,
          "showSoftClipping": true,
          "height": 130,
          "featureHeight": 1
        }
      ]
    }
  ]
}

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

The same view with jbrowse-anywidget, from the config the figure loads. Swap a track's uri for your own file.

from jbrowse_anywidget import JBrowseApp, fetch_hub

config = fetch_hub("https://jbrowse.org/demos/cancer_sv/config.json")
app = JBrowseApp(
    assemblies=config["assemblies"],
    tracks=config.get("tracks", []),
    views=[
        {
            "type": "LinearGenomeView",
            "assembly": "hg38",
            "loc": "chr3:25,357,600-25,361,000",
            "tracks": [
                {
                    "trackId": "ncbi_refseq_hg38",
                    "geneGlyphMode": "longestCoding",
                    "height": 70
                },
                {
                    "trackId": "COLO829_tumor_ont",
                    "type": "LinearAlignmentsDisplay",
                    "forceLoad": True,
                    "showSoftClipping": True,
                    "height": 130,
                    "featureHeight": 1
                }
            ]
        }
    ],
)
app  # 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": "hg38",
      "loc": "chr3:25,357,600-25,361,000",
      "tracks": [
        {
          "trackId": "ncbi_refseq_hg38",
          "geneGlyphMode": "longestCoding",
          "height": 70
        },
        {
          "trackId": "COLO829_tumor_ont",
          "type": "LinearAlignmentsDisplay",
          "forceLoad": true,
          "showSoftClipping": true,
          "height": 130,
          "featureHeight": 1
        }
      ]
    }
  ]
}
JSON

npx @jbrowse/capture --instance https://jbrowse.org/code/jb2/main/ \
  --config https://jbrowse.org/demos/cancer_sv/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.

Other callers

Anything that writes breakends or symbolic SVs to a VCF goes through the same two commands:

  • cuteSV, Sniffles, pbsv, Delly, Manta, GRIDSS all write a VCF that sv_multihop.py bedpe reads directly
  • LINX publishes clusters and chained links as TSVs. Convert the junction columns to the six BEDPE columns with awk; one --outDir per cluster gives a chromothripsis event as a contact sheet
  • PURPLE copy-number segments are not junctions. Convert the segment TSV to a bedGraph, bedGraphToBigWig it, and add it as a --bigwig so every image carries the copy number under the reads

Ordering breakends into a derivative chromosome needs allele-specific copy number and a centromere constraint. LINX derives both from PURPLE's purity and ploidy.

Reproduce it end to end

Everything on this page is the commands above against hosted files. The figure is three jb2export invocations: two breakpoint renders with one --loc per panel, one per sample, and a plain render of the derivative assembly, which are the sv_review_tumor, sv_review_normal and sv_review_derivative specs in website/scripts/specs/jbrowse-img.ts.

See also

References