View Single Post
Old 03-03-2012, 07:42 AM   #3
gringer
David Eccles (gringer)
 
Location: Wellington, New Zealand

Join Date: May 2011
Posts: 843
Default

Here's an R function for producing a heatmap based on the vsd data (together with colour legend and dendrogram):

Code:
make.heatmap <- function(data.vsd, columns = dim(data.vsd)[2], col = colorRampPalette(blues9)(100), log = FALSE, top = dim(data.vsd)[1], method = "canberra"){
  hm.dist <- dist(t(head(data.vsd[,columns],top)), method = method);
  hm.dend <- as.dendrogram(hclust(hm.dist));
  scaled.mat <- as.matrix(hm.dist);
  if(log){
    scaled.mat <- log(scaled.mat);
  }
  diag(scaled.mat) <- NA;
  scaled.mat <- scaled.mat[order.dendrogram(hm.dend),rev(order.dendrogram(hm.dend))];
  scaled.range <- t(as.matrix(seq(range(scaled.mat, na.rm = TRUE)[1],
                                  range(scaled.mat, na.rm = TRUE)[2], length.out = 100)));

  layout(matrix(1:3,1,3),widths = c(1,2,8));
  par(mar = c(10,3,0.5,0.5));
  image(1,scaled.range[1,],scaled.range, col = col,
        xlab = "", ylab = "", xaxt = "n");
  plot(hm.dend, axes = FALSE, horiz = TRUE, yaxs = "i", leaflab = "none");
  par(mar = c(10,0.5,0.5,10));
  image(1:length(rownames(scaled.mat)),1:length(colnames(scaled.mat)),scaled.mat,
        col = col,
        xlab = "", ylab = "", xaxt = "n", yaxt = "n");
  axis(4, 1:length(colnames(scaled.mat)), labels = sub("_seq1","",colnames(scaled.mat)), las = 2,
       line = -0.5, tick = 0, cex.axis = 0.2 + 1/log10(length(colnames(scaled.mat))));
  axis(1, 1:length(rownames(scaled.mat)), labels = rownames(scaled.mat), las = 2,
       line = -0.5, tick = 0, cex.axis = 0.2 + 1/log10(6));
  invisible(scaled.mat);
}
I recommend top <= 1000 if you're looking at all transcripts.
gringer is offline   Reply With Quote