Thursday, 21 April 2016

Multipanel ggplot with subfigure labels (and drawing igraph structures with ggplot2)

I want to make some multipanelled plots in R where, for a given figure:
- each of the subfigures is a ggplot (or other grid::grob object)
- each of the subfigures is given a label in the top-left hand corner

I've been using ggplot to draw graphs (the graph theory kind). It isn't quite set up to do this and indeed, igraph has perfectly good graph-drawing functions. They just don't return grobs (neither does sna::gplot), so they are a bit limited to my mind. Hence, I've set up a couple of simple functions to make ggplot figures out of igraph structures using the igraph::layout functions

Here's an example figure from my current dissertation. We first set up a theme to remove the weird default background for ggplot figures, import some libraries and define a function that returns a label in a rectangle.

library(igraph)
library(grid)
library(gridExtra)
library(ggplot2)

theme_nothing <- function(
  base_size = 12, base_family = "Helvetica"
  ){
  # blank plotting background for ggplot2
  # - for removing all axes etc when plotting igraphs in ggplot2
  # - copied from ggplot2 docs
  theme_bw(base_size = base_size, base_family = base_family
    ) %+replace% theme(
    rect = element_blank(),
    line = element_blank(),
    text = element_blank(),
    legend.position = 'none'
    )
  }


labelGrob <- function(label){
  # function to generate subfigure labels
  require(grid)
  require(gridExtra)
  textGrob(label, y = unit(0.9, 'npc'))
  }
Then we define a function to convert an igraph object into a dataframe of x/y coords so that it can be plotted in ggplot.
igraph_to_ggplot_df <- function( X, coord.func = layout.fruchterman.reingold ){  # Makes a dataframe of vertex / edge positions from an igraph # object stopifnot(is(X, "igraph"))
  # Give the vertices some names, if they don't already have them
  # then use the coord.func to define X/Y positions for each vertex
  v.names  <- names(V(X))
  if(is.null(v.names)){v.names <- paste0('v', 1:length(V(X)))}
  v.coords <- scale(coord.func(X))
  dimnames(v.coords) <- list(v.names, c('x', 'y'))
  # Obtain the edges for the graph, give them names (if they don't
  #   already have them) and then associate the start and end point
  #   of each edge with the X/Y coords of the relevant vertex
  e.list   <- get.edgelist(X)
  stopifnot(nrow(e.list) > 0)
  e.names <- names(E(X))
  if(is.null(e.names)){
    # assume the graph isn't null
    e.names <- paste0('e', 1:nrow(e.list))
    }
  dimnames(e.list) <- list(e.names, c('head', 'tail'))
  stopifnot(length(intersect(v.names, e.names)) == 0)
  # convert the edge / vertex coords into a single dataframe for use in ggplot
  edge_df <- data.frame(
    element.name = rep(e.names, 2),
    element.label = '',
    element.type = 'edge',
    x      = v.coords[as.vector(e.list[, c('head', 'tail')]), 'x'],
    y      = v.coords[as.vector(e.list[, c('head', 'tail')]), 'y']
    )
  vertex_df <- data.frame(
    element.name = rownames(v.coords),
    element.label = rownames(v.coords),
    element.type = 'vertex',
    x = v.coords[, 'x'],
    y = v.coords[, 'y']
    )
  graph_df <- rbind(edge_df, vertex_df)
  graph_df
  }
Then we define a function to turn igraph objects into ggplot graphical objects, so that these can be arranged and plotted out using grid.
igraph_ggplot
<- function( X, coord.func = layout.fruchterman.reingold, offset = 0.5 ){ # Converts a graph into x/y coords, # then plots the graph as straight-line joined coords d.frame <- igraph_to_ggplot_df(X, coord.func) # Note that if I don't change the offset myself, # arrangeGrob cuts off distal bits of my graph x.range <- range(d.frame$x) y.range <- range(d.frame$y) g <- ggplot(data = d.frame, aes(x = x, y = y, group = element.name, label = element.label) ) + geom_line( ) + geom_point(colour = 'darkgoldenrod1', size = 10 ) + geom_point(colour = 'white', size = 7 ) + geom_text(hjust = 0.5, vjust = 0.5, aes(colour = element.type) ) + theme_nothing( ) + xlim(x.range[1] - offset, x.range[2] + offset ) + ylim(y.range[1] - offset, y.range[2] + offset) g }
This is just an example of a figure generated with the above code. We define the incidence
matrix of a graph, and plot out the incidence, adjacency and Laplacian matrices of that graph
(as tables) and a picture of the graph itself. The labels are added using arrangeGrob and
 labelGrob: I define an arrangement of 8 panels (2rows x 4 columns). The labels are 
placed into the first and third columns (which are 1/4 the width of the second and fourth 
columns.)
fig1_2_1 <- function(){
  # Define a graph using it's incidence matrix, then plot
  # Signed incidence matrix
  D <- matrix(c(1, 1, 1, 0,
                -1,0, 0, 1,
                0,-1, 0,-1,
                0, 0,-1, 0,
                0, 0, 0, 0
                ),
              byrow = TRUE, nrow = 5,
              dimnames = list(paste0('v', 1:5),
                              paste0('e', 1:4)))
  B <- abs(D)        # Unsigned incidence matrix
  Q <- D %*% t(D)    # Laplacian
  deg <- diag(Q)     # Degree vector
  A <- diag(deg) - Q # Adjacency matrix
  
  X <- graph_from_adjacency_matrix(
    A, mode = 'undirected'
    )
  g1 <- tableGrob(A)
  g2 <- tableGrob(D)
  g3 <- igraph_ggplot(X, offset = 0.25)
  g4 <- tableGrob(Q)
  ag <- arrangeGrob(labelGrob('a'), g1,
                    labelGrob('b'), g2,
                    labelGrob('c'), g3,
                    labelGrob('d'), g4,
                    nrow = 2,
                    widths = c(1,4,1,4))

  grid.draw(ag)
  NULL
  }

------
Voila!


---------
# TODO - show how to do the same thing with grid.draw(a); grid.draw(b)
 where 'a' contains figures and 'b' contains panel labels
# TODO - eugh! work out why blogger is making my code/paragraphs so ugly?

Wednesday, 4 November 2015

Go from gene symbols to a karyogram showing the positions of those genes

Karyograms can be made quite easily in R using ggbio (see docs; the code below is based on section 8.4 thereof)
I want to make a karyogram that indicates the positions of a bunch of genes on the human genome.

Required packages:
library('ggbio')
library('GenomicRanges')
library('TxDb.Hsapiens.UCSC.hg19.knownGene')
txdb <- TxDb.Hsapiens.UCSC.hg19.knownGene
data(hg19Ideogram, package = "biovizBase")
data(hg19IdeogramCyto, package = "biovizBase")
hg19 <- keepSeqlevels(
    hg19Ideogram,
    paste0("chr", c(1:22, "X", "Y"))
    )
hg19.cyto <- keepSeqlevels(
    hg19IdeogramCyto,
    paste0("chr", c(1:22, "X", "Y"))
    )

A) Obtain genomic coordinates for the genes
Working in R, we want to obtain the genomic coordinates of the exon boundaries for a bunch of human gene symbols. Some of the symbols are deprecated, and some have been updated, some may not have been adopted by EntrezGene etc, but hopefully none have been corrupted by Excel

There are multiple transcript definitions for each gene. So, for each gene, we are going to condense the exons into non-overlapping regions (which won't necessarily represent any specific exon).

To do the above we need
i) to convert gene symbols to entrez gene ids:
# temp - assume this has been done
gene_ids <- c('1000', '1234', '54407')

ii) pull out the exon-boundaries for all transcripts associated with the entrez-gene ids:
exons <- select(
  txdb,
  keys = gene_ids,
  keytype = 'GENEID',
  columns = c('EXONCHROM', 'EXONSTART', 'EXONEND', 'TXID', 'TXSTRAND')
  )

iii) convert the exon boundaries into genomic-ranges objects:
# TODO: use gene symbols instead of GENEID as the names for gr
gr <- GRanges(
  seqnames = Rle(exons$EXONCHROM),
  ranges = with(exons,
    IRanges(start = EXONSTART, end = EXONEND, names = GENEID)),
  strand = Rle(exons$TXSTRAND),
  seqinfo = seqinfo(hg19)
  ) 

iv) collapse the exon boundaries for individual genes and construct a genomicRangesList indexed by the genes (isn't actually needed for the karyogram):
grl <- reduce(split(gr, names(gr)))

v) convert all the above into a function
# TODO...

B) Plot a basic karyogram with the genes overlaid:
# autoplot(unlist(grl), layout = 'karyogram')
p <- ggplot(hg19.cyto) + layout_karyogram(cytoband = TRUE)
p <- p + layout_karyogram(
    gr, geom = 'rect', ylim = c(11, 21), colour = 'blue'
    ) + xlab('')
p




Thursday, 6 August 2015

d3.js book

Noted that Scott Murray's d3.js book 'Interactive Visualisation for the Web' is available for free online. I've yet to find a detailed set of online tutorials for bokeh or d3 (although dashingd3js has a useful intro), so this looks promising.

Saturday, 14 February 2015

biographr intro

This is my first blog post, and indeed my first blog. I'm a cancer genomics guy at Liverpool Uni. I thought I'd collect some notes about genomics, statistics and data visualisation here, so that I can easily find them in the future. So there'll probably be some R/python/jags stuff: I often read papers and try to work out how the authors put their graphs etc together, but when it comes to generating the same types of figures myself I've already forgotten.
.. although I'll probably never post anything beyond this ..
Russ