Announcement

Collapse
No announcement yet.

help with a heatmap with deseq - legend?

Collapse
X
 
  • Filter
  • Time
  • Show
Clear All
new posts

  • help with a heatmap with deseq - legend?

    Hi,
    I have made with deseq a heatmap, but it does not show a legend/colorkey and I need to have it?

    I did not see elsewhere how colorkey/legends is made when heatmap is generated from sequencing copy number values, it should not de in log2 values or?
    ------------
    SMART - bioinfo.uni-plovdiv.bg

  • #2
    Dear Vebaev

    DESeq does not do heatmaps. The programming language R, in which DESeq is implemented, has functions to make heatmaps, and I recommend that you have a look at the manual pages of some of these functions, for instance

    ## install (do only once)
    source('http://www.bioconductor.org/biocLite.R')
    biocLite("gplots")

    ## load (do each time you need it)
    library("gplots")
    help("heatmap.2")


    As to transformation - the log2-transformation is problematic due to its singularity at 0 and its large slope for values close to 0. Doing no transformation is also problematic since then the dynamic range of the colour scale tends to be "wasted" on a sparse set of large values. My recommdation is to have a look at Section 6 of the DESeq vignette, Variance stabilizing transformation, where in Fig.13 you also see an example heatmap (but using the more primitive function heatmap that does not offer a color bar).

    Hope this helps,

    Wolfgang.
    Wolfgang Huber
    EMBL

    Comment


    • #3
      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.

      Comment


      • #4
        Thanks both of you, I will try your suggestions!
        ------------
        SMART - bioinfo.uni-plovdiv.bg

        Comment

        Working...
        X