R package

Multiple sequence alignments in R, as the interactive JBrowseMSA viewer and as ggplot2 figures.

• In: an alignment (FASTA, Stockholm, Clustal, A3M, or a Biostrings object), a tree (Newick, ape, treeio or ggtree), and the annotations to draw on them, such as domains, genes, highlights, clades or your own column tracks
• Out: an htmlwidget for RStudio, R Markdown, Quarto and Shiny; a ggplot2 figure drawn from the viewer’s own layout, with msaview_plot(); or the viewer’s own SVG export, with msaview_save()
• Layers: compose with +, the way ggplot2 and ggtree add a geom

remotes::install_github("GMOD/ggmsaview")

Gallery

• Each figure is drawn twice from data installed with the package: JBrowseMSA’s own SVG export (msaview_save()) on the left, msaview_plot() on the right.

• Each recipe runs as written after

  library(ggmsaview)
  e <- function(f) system.file("extdata", f, package = "ggmsaview")

Preproinsulin across nine vertebrates, beside its tree

w <- msaview(msa = e("insulin.aln"), tree = e("insulin.nh"))
autoplot(w)
msaview_save(w, "figure.svg")

The same alignment drawn as its differences from human

w <- msaview(msa = e("insulin.aln"), tree = e("insulin.nh"),
  relative_to = "Human")
autoplot(w)
msaview_save(w, "figure.svg")

Clade brackets over the insulin tree

w <- msaview(msa = e("insulin.aln"), tree = e("insulin.nh")) +
  geom_msa_clade(range = c("Pig", "Mouse"), tips = 5, mark = "bracket",
    label = "Mammals", color = "#1565c0") +
  geom_msa_clade(range = c("Xenopus", "Chicken"), tips = 3, mark = "bracket",
    label = "Non-mammals")
autoplot(w)
msaview_save(w, "figure.svg")

Ten Src-family kinases with their InterPro domains

w <- msaview(msa = e("kinase.aln"), tree = e("kinase.nh"),
  gff = e("kinase-domains.gff"), col_width = 1.5)
autoplot(w)
msaview_save(w, "figure.svg")

tRNA with its consensus secondary structure

w <- msaview(msa = e("trna.stock"))
autoplot(w)
msaview_save(w, "figure.svg")

A gene cluster across six genomes, with each genome’s clade

clades <- data.frame(name = paste0("Genome_", 1:6),
  clade = c("A", "A", "B", "B", "C", "C"))
w <- msaview(msa = e("gene-cluster.stock"), gff = e("gene-cluster.gff"),
  col_width = 1) +
  geom_msa_rowdata(clades) +
  geom_msa_strip("clade", header = "Clade") +
  geom_msa_features(x = "position", color = "Name", label = "Name",
    header = "Genes")
autoplot(w)
msaview_save(w, "figure.svg")

Installation

remotes::install_github("GMOD/ggmsaview")

The package ships its JavaScript, so it installs with no Node. Only msaview_save() runs Node, through react-msaview-cli.

For Bioconductor interop:

BiocManager::install(c("Biostrings", "ggtree", "treeio"))

Quick start

library(ggmsaview)

widget <- msaview(
  msa = ">human\nMVLSPADKTNVKAAWGKVGAHAGEYGAEALERMFLSFPTTKTYFPHFDLSH\n>mouse\nMVLSGEDKSNIKAAWGKIGGHGAEYGAEALERMFASFPTTKTYFPHFDVSH\n>goat\nMSLTRTERTIILSLWSKISTQADVIGTETLERLFSCYPQAKTYFPHFDLHS",
  tree = "((human:0.1,mouse:0.2):0.05,goat:0.3);"
)
widget
# htmlwidgets::saveWidget(widget, "alignment.html")  # self-contained HTML
# msaview_plot(widget)                                # ggplot2 figure, see below

ggplot2 figures

msaview_plot() draws the same widget with ggplot2, so ggsave(), themes and patchwork apply:

• The viewer’s own model runs in an embedded QuickJS engine (the QuickJSR package), so collapsed clades, hidden gappy columns, color schemes, domains, highlights and relative_to all come out as they do on screen. No Node or browser runs
• The figure draws the viewer’s row panels (strips and gene arrows), its column tracks (conservation, the sequence logo, arcs, Stockholm text rows and your own geom_msa_track()s)
• Domains, gene arrows and strips map fill to their category through a scale_fill_manual() seeded with the viewer’s colors, so a scale_fill_*() you add replaces it, and patchwork collects the legends
• msaview_layout() returns the tables behind the figure (rows, letter colors, tree segments, feature shapes, row panels, highlights, clades, column tracks) for your own layers
• A widget knitted into a PDF or Word document from R Markdown or Quarto prints as this figure

library(ggplot2)
p <- msaview_plot(widget, tracks = c("conservation", "sequence-logo"))
ggsave("alignment.pdf", p, width = 8, height = 2)

SVG from the viewer’s own renderer

msaview_save() writes the SVG the viewer exports, or a PNG or PDF of it through the rsvg package:

• It needs Node.js, and runs npx react-msaview-cli unless options(msaview.cli = ...) names another copy
• No browser runs, so it works in a batch job or on a server
• width sets the view’s pixel width; the whole alignment is drawn whatever its length
• The CLI refuses an alignment of more than 500,000 residue letters, since each letter is an SVG element; a col_width below 5 or a row_height below 8 draws it with no letters

msaview_save(widget, "alignment.svg")
msaview_save(widget, "alignment.png", scale = 3)

Examples

A protease family in R is the long version: fourteen UniProt accessions cut to their peptidase S1 domain, aligned with DECIPHER, and drawn with a BLOSUM62 track, disulfide arcs and a band on each catalytic residue.

Layers, composed with +

Every argument of msaview() is also a layer, added the way ggplot2 and ggtree add a geom. Each layer reads a data frame, so a table of variants or per-column numbers goes in as it already sits in your session.

library(ggmsaview)

variants <- data.frame(
  row   = c("Human", "Human"),
  start = c(248, 273),
  end   = c(248, 273),
  label = c("R248Q", "R273H")
)

msaview(msa = "p53.aln", tree = "p53.nh") +
  geom_msa_domains("p53-domains.gff") +
  geom_msa_highlight(variants) +
  geom_msa_track(data.frame(value = conservation), name = "Conservation") +
  scale_residue_color("clustal", encoding = "color") +
  stat_msa_diff("Human") +
  coord_msa(100, 288, row = "Human") +
  theme_msa("dark", col_width = 14)

geom_msa_highlight(), geom_msa_clade(), geom_msa_track(), geom_msa_structure(), geom_msa_strip(), geom_msa_features() and scale_row_color() accumulate, so calling one twice adds two. The others replace what an earlier layer set.

A layer builds the same props the matching msaview() argument does, so the two styles mix freely and produce the same viewer. msa and tree are the viewer itself and stay arguments; every other argument has a layer, which test-layer-coverage.R checks.

Named character vector

seqs <- c(
  human = "MVLSPADKTNVKAAWGKVGAHAGEYGAEALERMFLSFPTTKTYFPHFDLSH",
  mouse = "MVLSGEDKSNIKAAWGKIGGHGAEYGAEALERMFASFPTTKTYFPHFDVSH",
  goat  = "MSLTRTERTIILSLWSKISTQADVIGTETLERLFSCYPQAKTYFPHFDLHS"
)
msaview(msa = seqs, tree = "((human:0.1,mouse:0.2):0.05,goat:0.3);")

From files

msaview(msa = "alignment.stock")
msaview(msa = "alignment.fa", tree = "tree.nwk")

With ape

library(ape)

tree <- rtree(15)
seqs <- setNames(
  replicate(15, paste0(sample(c("A","C","G","T"), 300, TRUE), collapse = "")),
  tree$tip.label
)
msaview(msa = seqs, tree = tree, color_scheme = "nucleotide")

Cladogram (no branch lengths)

msaview(msa = seqs, tree = tree, show_branch_len = FALSE)

With Biostrings

library(Biostrings)

# DNAStringSet
dna <- DNAStringSet(c(
  seq1 = "ATGCGATCGATCGATCG--ATCG",
  seq2 = "ATGCGATCGATCGATCGATCGATCG",
  seq3 = "ATGCG--CGATCGATCGATCGATCG",
  seq4 = "ATGCGATCGATCGATCG--ATCG"
))
msaview(msa = dna, color_scheme = "nucleotide")

# AAStringSet
aa <- AAStringSet(c(
  human = "MVLSPADKTNVKAAWGKVGAHAGEYGAEALERMFLSFPTTKTYFPHFDLSH",
  mouse = "MVLSGEDKSNIKAAWGKIGGHGAEYGAEALERMFASFPTTKTYFPHFDVSH",
  goat  = "MSLTRTERTIILSLWSKISTQADVIGTETLERLFSCYPQAKTYFPHFDLHS"
))
msaview(msa = aa, color_scheme = "clustal")

# DNAMultipleAlignment
# Biostrings rejects a MultipleAlignment whose rows differ in length
aln <- DNAMultipleAlignment(c(
  seq1 = "ATGCGATCGATCGATCG--ATCG",
  seq2 = "ATGCGATCGATCGATCGATCGAC",
  seq3 = "ATGCG--CGATCGATCGATCGAC"
))
msaview(msa = aln)

After running the msa package

library(msa)
library(Biostrings)

unaligned <- readAAStringSet("proteins.fasta")
aligned <- msa(unaligned, method = "ClustalOmega")
msaview(msa = as(aligned, "AAStringSet"))

Domain annotations

Overlay InterProScan-style protein domains as boxes on the alignment. The gff argument accepts a file path, a GFF3 string, or a data frame. This pairs with the CLI, which writes domains as GFF3:

react-msaview-cli interpro accessions.tsv -o domains.gff

seqs <- c(
  GPCR_human  = "MNGTEGPNFYVPFSNATGVVRSPFEYPQYYLAEPWQFSMLAAYMFLLIVLGFPINFLTLYVTVQHKKLR",
  GPCR_mouse  = "MNGTEGPNFYVPFSNKTGVVRSPFEYPQYYLAEPWQFSMLAAYMFLLIMLGFPINFLTLYVTVQHKKLR",
  GPCR_bovine = "MNGTEGPNFYVPFSNATGVVRSPFEYPQYYLAEPWQFSMLAAYMFLLIVLGFPINFLTLYVTVQHKKLR"
)

# from a CLI-generated GFF file
msaview(msa = seqs, gff = "domains.gff")

# or from a data frame: row (or seqname), start, end, and any fields,
# such as a name, a description, a strand or a column an encoding reads
domains <- data.frame(
  seqname     = names(seqs),
  start       = 6,
  end         = 62,
  name        = "PF00001",
  description = "7tm receptor (rhodopsin family)"
)
msaview(msa = seqs, gff = domains, color_scheme = "clustalx_protein_dynamic")

Labeled highlights

Mark a residue, a column range, or a set of rows, with a label saved in the widget. Coordinates are 1-based and inclusive; row makes start and end residues of that sequence, projected through the alignment’s gaps.

msaview(
  msa = seqs,
  highlights = list(
    list(row = "GPCR_human", start = 20, end = 26, label = "ligand pocket"),
    list(rows = c("GPCR_mouse"), label = "knockout", color = "rgba(0,120,255,0.2)")
  )
)

A track from your own numbers

column_tracks draws values you compute per column, or per residue of one sequence, as a track above the alignment. The viewer scales bars by max (default 1). When row names a sequence, the viewer places each value on that sequence’s residues and skips the columns where the sequence has a gap.

hydropathy <- c(I = 4.5, V = 4.2, L = 3.8, F = 2.8, C = 2.5, M = 1.9, A = 1.8,
                G = -0.4, T = -0.7, S = -0.8, W = -0.9, Y = -1.3, P = -1.6,
                H = -3.2, E = -3.5, Q = -3.5, D = -3.5, N = -3.5, K = -3.9,
                R = -4.5)
residues <- strsplit(gsub("-", "", seqs[["GPCR_human"]]), "")[[1]]

msaview(
  msa = seqs,
  column_tracks = list(
    list(id = "kd", name = "Hydropathy (human)", kind = "bar",
         values = hydropathy[residues] + 4.5,
         scales = list(y = list(domainMax = 9)), row = "GPCR_human",
         color = "#6a51a3")
  )
)

A kind = "text" track takes data, one character per column, and an optional encoding = list(color = list(domain =, range =)), a color per character. Every field is listed in the layers reference.

With ggtree

Pass a ggtree plot object directly as the tree argument, and msaview reads the tree out of it.

library(ggtree)
library(ape)

tree <- rtree(20)
p <- ggtree(tree) + geom_tiplab()

seqs <- setNames(
  replicate(20, paste0(sample(c("A","C","G","T"), 200, TRUE), collapse = "")),
  tree$tip.label
)

# pass the ggtree plot object directly
msaview(msa = seqs, tree = p, color_scheme = "nucleotide")

ggtree with annotations

Annotated ggtree plots work too. msaview takes only the tree; the colors, labels and metadata stay in the ggtree plot.

metadata <- data.frame(
  label = tree$tip.label,
  group = sample(c("A", "B", "C"), 20, replace = TRUE)
)
p2 <- ggtree(tree) %<+% metadata +
  geom_tiplab(aes(color = group)) +
  theme(legend.position = "right")

msaview(msa = seqs, tree = p2)

With treeio

library(treeio)

# treedata objects from read.beast, read.raxml, etc.
beast_tree <- read.beast("beast_output.tree")
msaview(tree = beast_tree)

# or convert from phylo
td <- as.treedata(tree)
msaview(msa = seqs, tree = td)

In Shiny

The widget sets three inputs named after its output id. input$<id>_click holds the cell a click pinned: column is the 1-based column of the file, row the row name, residue the 1-based position in that row’s sequence (absent on a gap), and letter the character. It is NULL after a click clears it. input$<id>_viewport holds the columns on screen as startColumn and endColumn. input$<id>_selection holds the block a shift-drag selects: start and end columns of the file, and rows, the selected row names, absent when the block spans every row. It is NULL after Escape clears it. A re-render that keeps the same alignment keeps the reader’s scroll and zoom, so changing color_scheme below restyles the view in place.

library(shiny)
library(ggmsaview)

ui <- fluidPage(
  titlePanel("MSA Viewer"),
  sidebarLayout(
    sidebarPanel(
      fileInput("msa_file", "Upload alignment",
                accept = c(".fa", ".fasta", ".stock", ".sto", ".aln")),
      fileInput("tree_file", "Upload tree (optional)",
                accept = c(".nwk", ".nh", ".newick")),
      selectInput("color_scheme", "Color scheme",
                  choices = c("maeditor", "clustal", "lesk", "cinema",
                              "flower", "nucleotide", "jbrowse_dna",
                              "clustalx_protein_dynamic",
                              "percent_identity_dynamic"),
                  selected = "maeditor"),
      actionButton("example_btn", "Load example data")
    ),
    mainPanel(
      msaviewOutput("msa_viewer", height = "600px"),
      verbatimTextOutput("clicked")
    )
  )
)

server <- function(input, output, session) {
  msa_data <- reactiveVal(NULL)
  tree_data <- reactiveVal(NULL)

  observeEvent(input$msa_file, {
    msa_data(input$msa_file$datapath)
  })

  observeEvent(input$tree_file, {
    tree_data(input$tree_file$datapath)
  })

  observeEvent(input$example_btn, {
    msa_data(paste0(
      ">human\nMVLSPADKTNVKAAWGKVGAHAGEYGAEALERMFLSFPTTKTYFPHFDLSH\n",
      ">mouse\nMVLSGEDKSNIKAAWGKIGGHGAEYGAEALERMFASFPTTKTYFPHFDVSH\n",
      ">goat\nMSLTRTERTIILSLWSKISTQADVIGTETLERLFSCYPQAKTYFPHFDLHS"
    ))
    tree_data("((human:0.1,mouse:0.2):0.05,goat:0.3);")
  })

  output$msa_viewer <- renderMsaview({
    req(msa_data())
    msaview(
      msa = msa_data(),
      tree = tree_data(),
      color_scheme = input$color_scheme
    )
  })

  output$clicked <- renderPrint(input$msa_viewer_click)
}

shinyApp(ui, server)

Supported inputs

Color schemes

Protein: maeditor, clustal, lesk, cinema, flower, clustalx_protein, jalview_taylor, jalview_zappo, jalview_hydrophobicity, jalview_buried, jalview_prophelix, jalview_propstrand, jalview_propturn

Nucleotide: nucleotide, jbrowse_dna, rainbow_dna, clustalx_dna

Dynamic (computed per-column): clustalx_protein_dynamic, percent_identity_dynamic

A named vector of colors keyed by residue letter colors those letters, in either case, and leaves every other letter uncolored:

msaview(msa = aa, color_scheme = c(K = "#1f77b4", R = "#1f77b4", D = "#d62728"))