RNA-seq visualization

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

An RNA-seq read mapped back to the genome jumps the introns spliced out of the transcript, and the aligner records each jump in the read's CIGAR string. In a stranded library the pair flags also mark which strand the transcript came from. JBrowse draws splice arcs and strand coloring from those two BAM fields, with no extra files or configuration. We read spliced alignments over ACTB, transcript strand at the surfeit locus and a long-read alignment, then load a pipeline's junction table as a track.

Prerequisites

for your own reads, an aligned, sorted and indexed BAM or CRAM from a spliced aligner
a JBrowse instance to load it into: the web quickstart or the desktop quickstart

Where the data comes from

The hg19 alignments, gene models and junction table behind the figures, hosted on jbrowse.org's demo bucket.

Nothing to download: the track configs below read these files by URL.

The files

the paired-end stranded RNA-seq alignments behind every short-read figure on this page, from the sample files listed at RSeQC: https://s3.amazonaws.com/jbrowse.org/genomes/hg19/paired_end_rnaseq/Pairend_StrandSpecific_51mer_Human_hg19.bam

Loading hg19, the stranded reads and RefSeq genes

The alignments are against hg19, so we load that assembly with an alias file that maps the chromosome names in the BAM and the gene models to the assembly's. Then we add the stranded paired-end reads and the RefSeq gene models, each with its index beside it:

Config file

Goes in the assemblies array of config.json. See Assemblies.

{
  "name": "hg19",
  "uri": "https://jbrowse.org/genomes/hg19/fasta/hg19.fa.gz",
  "refNameAliases": {
    "uri": "https://s3.amazonaws.com/jbrowse.org/genomes/hg19/hg19_aliases.txt"
  },
  "cytobands": "https://jbrowse.org/ucsc/hg19/cytoBandIdeo.bed.gz"
}

CLI (add-assembly)

jbrowse add-assembly https://jbrowse.org/genomes/hg19/fasta/hg19.fa.gz \
  --name hg19 \
  --refNameAliases https://s3.amazonaws.com/jbrowse.org/genomes/hg19/hg19_aliases.txt \
  --config '{"cytobands":"https://jbrowse.org/ucsc/hg19/cytoBandIdeo.bed.gz"}'

JBrowse Desktop

hg19 is one of the genomes JBrowse Desktop hosts, with gene tracks already set up: on the start screen click Show all available genomes and pick it. To load the files of this config instead:

In 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:

https://jbrowse.org/genomes/hg19/fasta/hg19.fa.gz
https://jbrowse.org/genomes/hg19/fasta/hg19.fa.gz.fai
https://jbrowse.org/genomes/hg19/fasta/hg19.fa.gz.gzi

JBrowse reads the format off the file name. Then fill in:

Genome name: hg19
refName aliases (under More options): https://s3.amazonaws.com/jbrowse.org/genomes/hg19/hg19_aliases.txt
cytobands (under More options): https://jbrowse.org/ucsc/hg19/cytoBandIdeo.bed.gz

Config file

Goes in the tracks array of config.json. See Tracks.

{
  "trackId": "rnaseq_paired_stranded",
  "name": "Paired-end stranded RNA-seq (RSeQC sample)",
  "uri": "https://s3.amazonaws.com/jbrowse.org/genomes/hg19/paired_end_rnaseq/Pairend_StrandSpecific_51mer_Human_hg19.bam",
  "assemblyNames": ["hg19"]
}

CLI (add-track)

jbrowse add-track https://s3.amazonaws.com/jbrowse.org/genomes/hg19/paired_end_rnaseq/Pairend_StrandSpecific_51mer_Human_hg19.bam \
  --trackId rnaseq_paired_stranded \
  --name "Paired-end stranded RNA-seq (RSeQC sample)" \
  --assemblyNames hg19

In the app

In JBrowse Desktop, or in any running JBrowse Web session, open a view on this track’s assembly, then File → Open track... and, in Add a track from file or URL, enter:

Main file: https://s3.amazonaws.com/jbrowse.org/genomes/hg19/paired_end_rnaseq/Pairend_StrandSpecific_51mer_Human_hg19.bam

Click Next. JBrowse reads the adapter and track type off the file name. Then fill in:

Track name: Paired-end stranded RNA-seq (RSeQC sample)
Assembly: hg19

Click Add.

Open in JBrowse

Open this track in JBrowse Web ↗

Open this track in JBrowse Desktop ↗ (JBrowse Desktop 5.0+)

Config file

Goes in the tracks array of config.json. See Tracks.

{
  "trackId": "ncbi_refseq_hg19",
  "name": "NCBI RefSeq genes (hg19)",
  "uri": "https://s3.amazonaws.com/jbrowse.org/genomes/hg19/ncbi_refseq/GRCh37_latest_genomic.sort.gff.gz",
  "assemblyNames": ["hg19"]
}

CLI (add-track)

jbrowse add-track https://s3.amazonaws.com/jbrowse.org/genomes/hg19/ncbi_refseq/GRCh37_latest_genomic.sort.gff.gz \
  --trackId ncbi_refseq_hg19 \
  --name "NCBI RefSeq genes (hg19)" \
  --assemblyNames hg19

In the app

In JBrowse Desktop, or in any running JBrowse Web session, open a view on this track’s assembly, then File → Open track... and, in Add a track from file or URL, enter:

Main file: https://s3.amazonaws.com/jbrowse.org/genomes/hg19/ncbi_refseq/GRCh37_latest_genomic.sort.gff.gz

Click Next. JBrowse reads the adapter and track type off the file name. Then fill in:

Track name: NCBI RefSeq genes (hg19)
Assembly: hg19

Click Add.

Loading your own RNA-seq data

An aligned, sorted and indexed BAM or CRAM loads as an AlignmentsTrack from one uri, and JBrowse finds the .bai or .crai beside the file:

Config file

Goes in the tracks array of config.json. See Tracks.

{
  "trackId": "my_rnaseq",
  "name": "My RNA-seq",
  "uri": "https://yourhost/rnaseq.bam",
  "assemblyNames": ["hg38"]
}

CLI (add-track)

jbrowse add-track https://yourhost/rnaseq.bam \
  --trackId my_rnaseq \
  --name "My RNA-seq" \
  --assemblyNames hg38

In the app

In JBrowse Desktop, or in any running JBrowse Web session, open a view on this track’s assembly, then File → Open track... and, in Add a track from file or URL, enter:

Main file: https://yourhost/rnaseq.bam

Click Next. JBrowse reads the adapter and track type off the file name. Then fill in:

Track name: My RNA-seq
Assembly: hg38

Click Add.

Your own reads load like the track above once three things are in place:

Assembly. assemblyNames must match an assembly already configured in JBrowse; see the assemblies configuration guide.
Alignment. Align reads with a spliced aligner such as STAR, then samtools sort and samtools index so the .bai sits beside the BAM.
Coverage signal. A precomputed signal, such as a strand-specific BigWig from the aligner, loads separately as a quantitative track.

What RNA-seq looks like in the genome browser

Type chr7:5,567,000-5,570,000 into the location box. The gene there is ACTB, a compact gene with deep, even read coverage.

The thin teal lines jumping across a gap are spliced alignments, where a read maps partly to one exon and partly to the next, skipping the intron between them.

RNA-seq reads over ACTB under the NCBI RefSeq gene model: the coverage histogram with strand-colored splice arcs, then the spliced read pileup.

Read coverage and read height

The histogram counts the reads in the pileup below it at each position.

Read height → Compact in the track menu packs the full read stack into view:

ACTB under compact read height: the whole read stack fits the track, under the per-position coverage histogram and the hg19 NCBI RefSeq gene model.

Spliced reads, CIGAR strings, and splice arcs

A spliced aligner like STAR split-maps a read that crosses an intron and encodes the skip in its CIGAR string, the SAM/BAM field describing how a read aligns to the reference.

A spliced read in the ACTB pileup (reads here are 51 bp) has a CIGAR like this, spaced out for readability:

18M 95N 33M

The CIGAR reads as 18 bp (M, match) aligned to one exon, a 95 bp skip (N) across the intron, and 33 bp (M) aligned to the next.

JBrowse computes the arcs on the fly from the skips in the reads in view, and colors each arc by transcript strand, salmon for forward and purple for reverse. It finds the strand two ways:

Tags. XS and TS record the strand directly. minimap2's ts is relative to the read, so JBrowse combines it with the strand the read aligned to.
Splice motif. A BAM from STAR without --outSAMstrandField intronMotif has no tags, so JBrowse reads the first and last two bases of the intron off the reference: GT-AG on the forward strand reads as CT-AC on the reverse. A junction whose reads disagree, or whose motif is none of GT-AG, GC-AG and AT-AC, draws in the neutral color. Hovering an arc shows the motif beside the read count.

Picking out spliced reads in a deep pileup

In a deep pileup, reads with a skip sit among many more that have none. Three track-menu settings pull the splicing evidence out.

Sort by... → Spliced reads first gives every read whose CIGAR has a skip the lowest rows, so the junction-spanning reads sit together at the top of the pileup.

The ACTB pileup in file order above, and sorted with spliced reads first below. Teal lines mark reads whose CIGAR has a skip. File order scatters them down the stack, and the sort gathers them into the top rows.

Filter by... → Spliced reads → Only spliced reads sets a filter that drops every read without a skip, and the coverage histogram follows, leaving the junction-spanning evidence alone. Hide spliced reads in the same submenu is the complement.

Sashimi arcs → Hide non-canonical junctions drops every arc whose intron does not begin and end with GT-AG, GC-AG or AT-AC on either strand. Raising Sashimi arcs → Min read support also drops thin arcs, along with any real junction few reads support.

Every junction the reads have above, and only the canonical ones below. The motif filter drops the two salmon arcs around the leftmost intron and keeps the purple ones, whose introns the gene model also draws.

Strand-specific RNA-seq

Splice-arc colors give the strand of spliced reads only. A strand-specific library records the transcript strand in which mate of the pair a read is, so every read has it, which separates genes close together or overlapping on opposite strands.

The surfeit locus is a tight cluster of genes on alternating strands (RPL7A, SURF1, SURF2, SURF4). Type chr9:136,214,000-136,229,000 into the location box, then pick this path from the track menu:

Color by... → Paired end → First-of-pair strand

The surfeit locus colored by first-of-pair strand. The pileup splits into two colors, and the switch falls where the genes change strand: RPL7A forward, SURF1 reverse, SURF2 forward.

Grouping the reads by first-of-pair strand shows transcript strand over a whole gene.

Group by... → First-of-pair strand applies a facet, one section per group, each with a coverage histogram from only that group's reads, forward and reverse on one autoscaled axis. Turning off Show... → Show pileup leaves the two histograms alone.

On chr6, NELFE and SKIV2L sit back to back on opposite strands. Type chr6:31,920,500-31,933,000 into the location box:

NELFE and SKIV2L, adjacent and on opposite strands, grouped by first-of-pair strand: each band shows signal over exactly one of the two genes.

Long-read splicing: IsoSeq over ACTB

A long read (PacBio IsoSeq, Nanopore) often spans a whole transcript, aligning across every exon with one N skip per intron, and JBrowse derives the same arcs and connectors from those skips. The IsoSeq alignments are a second BAM on the same assembly:

Config file

Goes in the tracks array of config.json. See Tracks.

{
  "trackId": "alzheimers_isoseq",
  "name": "IsoSeq high-quality isoforms",
  "uri": "https://s3.amazonaws.com/jbrowse.org/genomes/hg19/alzheimers_isoseq/hq_isoforms.fasta.bam",
  "assemblyNames": ["hg19"]
}

CLI (add-track)

jbrowse add-track https://s3.amazonaws.com/jbrowse.org/genomes/hg19/alzheimers_isoseq/hq_isoforms.fasta.bam \
  --trackId alzheimers_isoseq \
  --name "IsoSeq high-quality isoforms" \
  --assemblyNames hg19

In the app

In JBrowse Desktop, or in any running JBrowse Web session, open a view on this track’s assembly, then File → Open track... and, in Add a track from file or URL, enter:

Main file: https://s3.amazonaws.com/jbrowse.org/genomes/hg19/alzheimers_isoseq/hq_isoforms.fasta.bam

Click Next. JBrowse reads the adapter and track type off the file name. Then fill in:

Track name: IsoSeq high-quality isoforms
Assembly: hg19

Click Add.

Open in JBrowse

Open this track in JBrowse Web ↗

Open this track in JBrowse Desktop ↗ (JBrowse Desktop 5.0+)

IsoSeq high-quality isoforms over ACTB. A long read often spans all of a transcript's exons, so one spliced alignment covers the whole transcript.

Loading STAR, regtools and portcullis junction tables

A junction table from an RNA-seq pipeline adds what JBrowse cannot compute from the reads in view: whole-library read counts, an annotated-or-novel flag and portcullis's filtering verdict. One awk line converts each table to BED, with the intron as the interval, the read count as the score, then the strand and the tool's own columns, and the BED draws as arcs on a feature track. Sort and index all three the same way:

sort -k1,1 -k2,2n junctions.bed | bgzip > junctions.bed.gz
tabix -p bed junctions.bed.gz

STAR writes SJ.out.tab with the intron as 1-based inclusive coordinates, then the strand, the motif as a STAR code (0 is non-canonical), the annotated flag and the unique and multi-mapping read counts:

awk -v OFS='\t' '{
  # column 4: 1 = +, 2 = -, 0 = undetermined
  s = $4 == 1 ? "+" : $4 == 2 ? "-" : "."
  # BED is 0-based half-open, so the start moves back one
  print $1, $2 - 1, $3, "junc" NR, $7, s, $5, $6
}' SJ.out.tab > junctions.bed

regtools reads the junctions straight out of the BAM, then annotates each against a GTF with its splice-site motif and whether an annotated transcript joins that donor to that acceptor. The annotated table has a header row and gives the last base of one exon and the first base of the next, so the intron ends one base before end:

# -s FR: a stranded library whose first read runs along the transcript;
#   RF for a dUTP library
regtools junctions extract -s FR rnaseq.bam -o regtools.bed
# annotate reads a plain GTF, not a gzipped one
regtools junctions annotate regtools.bed genome.fa genes.gtf -o annotated.tsv
# column 7 is the motif, column 14 the known-junction flag
awk -F'\t' -v OFS='\t' 'NR > 1 { print $1, $2, $3 - 1, $4, $5, $6, $7, $14 }' \
  annotated.tsv > junctions.bed

portcullis writes a header row naming its columns, so the recipe reads them by name. nb_raw_aln is the raw supporting-read count, canonical_ss is C, S or N for canonical, semi-canonical and non-canonical:

awk -F'\t' -v OFS='\t' '
  NR == 1 { for (i = 1; i <= NF; i++) c[$i] = i; next }
  {
    # start is 0-based and end is inclusive, so end moves forward one
    print $c["refname"], $c["start"], $c["end"] + 1, "junc_" $c["index"],
      $c["nb_raw_aln"], $c["consensus-strand"], $c["canonical_ss"]
  }
' 3-filt/portcullis_filtered.pass.junctions.tab > junctions.bed

This config loads the regtools file built from this page's BAM as a feature track. In grammar-of-graphics terms, the link mark is the shape drawn for each junction, size scales with read count, color is a color scale on the known-junction flag against RefSeq, and filter drops records. columnNames names the extra columns so the color scale can read them:

Config file

Goes in the tracks array of config.json. See Tracks.

{
  "type": "FeatureTrack",
  "trackId": "rnaseq_junctions_hg19",
  "name": "Splice junctions (regtools)",
  "category": ["RNA-seq"],
  "assemblyNames": ["hg19"],
  "adapter": {
    "type": "BedTabixAdapter",
    "uri": "https://jbrowse.org/demos/rnaseq/rnaseq_junctions.bed.gz",
    "columnNames": [
      "chrom",
      "chromStart",
      "chromEnd",
      "name",
      "score",
      "strand",
      "splice_site",
      "known_junction"
    ]
  },
  "displays": [
    {
      "type": "LinearMarkDisplay",
      "displayId": "rnaseq_junctions_hg19-LinearMarkDisplay",
      "transform": [
        {
          "type": "filter",
          "expr": "jexl:feature.score >= 3 && feature.splice_site in ['GT-AG', 'GC-AG', 'AT-AC']"
        }
      ],
      "marks": [
        {
          "mark": "link",
          "encoding": {
            "size": { "field": "score", "scale": "log", "range": [1, 8] },
            "color": {
              "field": "known_junction",
              "scale": "categorical",
              "domain": ["1", "0"],
              "range": ["#377eb8", "#e41a1c"],
              "labels": ["annotated", "novel"],
              "title": "RefSeq"
            }
          }
        },
        { "mark": "text", "encoding": { "text": "score" }, "maxBpPerPx": 50 }
      ]
    }
  ]
}

CLI (add-track-json)

jbrowse add-track-json '{
  "type": "FeatureTrack",
  "trackId": "rnaseq_junctions_hg19",
  "name": "Splice junctions (regtools)",
  "category": ["RNA-seq"],
  "assemblyNames": ["hg19"],
  "adapter": {
    "type": "BedTabixAdapter",
    "uri": "https://jbrowse.org/demos/rnaseq/rnaseq_junctions.bed.gz",
    "columnNames": [
      "chrom",
      "chromStart",
      "chromEnd",
      "name",
      "score",
      "strand",
      "splice_site",
      "known_junction"
    ]
  },
  "displays": [
    {
      "type": "LinearMarkDisplay",
      "displayId": "rnaseq_junctions_hg19-LinearMarkDisplay",
      "transform": [
        {
          "type": "filter",
          "expr": "jexl:feature.score >= 3 && feature.splice_site in ['\''GT-AG'\'', '\''GC-AG'\'', '\''AT-AC'\'']"
        }
      ],
      "marks": [
        {
          "mark": "link",
          "encoding": {
            "size": { "field": "score", "scale": "log", "range": [1, 8] },
            "color": {
              "field": "known_junction",
              "scale": "categorical",
              "domain": ["1", "0"],
              "range": ["#377eb8", "#e41a1c"],
              "labels": ["annotated", "novel"],
              "title": "RefSeq"
            }
          }
        },
        { "mark": "text", "encoding": { "text": "score" }, "maxBpPerPx": 50 }
      ]
    }
  ]
}'

In the app

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": "FeatureTrack",
  "trackId": "rnaseq_junctions_hg19",
  "name": "Splice junctions (regtools)",
  "category": ["RNA-seq"],
  "assemblyNames": ["hg19"],
  "adapter": {
    "type": "BedTabixAdapter",
    "uri": "https://jbrowse.org/demos/rnaseq/rnaseq_junctions.bed.gz",
    "columnNames": [
      "chrom",
      "chromStart",
      "chromEnd",
      "name",
      "score",
      "strand",
      "splice_site",
      "known_junction"
    ]
  },
  "displays": [
    {
      "type": "LinearMarkDisplay",
      "displayId": "rnaseq_junctions_hg19-LinearMarkDisplay",
      "transform": [
        {
          "type": "filter",
          "expr": "jexl:feature.score >= 3 && feature.splice_site in ['GT-AG', 'GC-AG', 'AT-AC']"
        }
      ],
      "marks": [
        {
          "mark": "link",
          "encoding": {
            "size": { "field": "score", "scale": "log", "range": [1, 8] },
            "color": {
              "field": "known_junction",
              "scale": "categorical",
              "domain": ["1", "0"],
              "range": ["#377eb8", "#e41a1c"],
              "labels": ["annotated", "novel"],
              "title": "RefSeq"
            }
          }
        },
        { "mark": "text", "encoding": { "text": "score" }, "maxBpPerPx": 50 }
      ]
    }
  ]
}

Open in JBrowse

Open this track in JBrowse Web ↗

Open this track in JBrowse Desktop ↗ (JBrowse Desktop 5.0+)

The filter applies the sashimi menu's read-support floor to the whole-library counts in the file and keeps the canonical motifs, which drops the duplicate that wrong-strand reads give a junction on the opposite strand. The color's domain lists the known_junction values as strings because the adapter reads extra columns as text.

Advanced → Edit plot... in the track menu edits the same marks, colors and filter on an open track.

Type chr11:49,220,500-49,231,500 into the location box, the 5' end of FOLH1, where the file separates what RefSeq annotates from what this library also splices:

The 5' end of FOLH1 on hg19: RefSeq transcripts above, the library's junctions below, blue where an annotated transcript joins the two ends and red where none does. The large red arc joins a donor and an acceptor RefSeq uses, skipping the exons between them, and it is about as thick as the blue arcs beside it.

To load the other tables:

For the STAR file, the columnNames end in motif and annotated, and the color's field is annotated.
For the portcullis file, the color's field is canonical_ss, its domain ["C", "S", "N"], and its labels canonical, semi-canonical and non-canonical, with a third color in range.

RNA-seq visualization

Prerequisites#

Where the data comes from#

Loading hg19, the stranded reads and RefSeq genes#

Loading your own RNA-seq data#

What RNA-seq looks like in the genome browser#

Read coverage and read height#

Spliced reads, CIGAR strings, and splice arcs#

Picking out spliced reads in a deep pileup#

Strand-specific RNA-seq#

Long-read splicing: IsoSeq over ACTB#

Loading STAR, regtools and portcullis junction tables#

See also#

Prerequisites

Where the data comes from

Loading hg19, the stranded reads and RefSeq genes

Loading your own RNA-seq data

What RNA-seq looks like in the genome browser

Read coverage and read height

Spliced reads, CIGAR strings, and splice arcs

Picking out spliced reads in a deep pileup

Strand-specific RNA-seq

Long-read splicing: IsoSeq over ACTB

Loading STAR, regtools and portcullis junction tables

See also