Georgii
Georgii

Reputation: 21

Modifying Continuous Legend in Pheatmap

I am having issues with pheatmap legend annotations. I want to annotate my columns with metadata from my coldata (working with RNAseq, DESeq2).

DESeq2:

dds <- DESeqDataSetFromMatrix(countData = counts, # Raw expression counts.
                                  colData = coldata, # Metadata with age & sex lables.
                                  design = ~ age)


dds <- DESeq(dds, test="LRT", reduced = ~ 1) # LRT test for model fitting.
res <- results(dds)

vsd <- vst(dds, blind = FALSE)

plotPCA(vsd, intgroup = c("age","sex", "batch"), returnData = FALSE)
# Create annotation 
heat_anno <- as.data.frame(colData(vsd)[,c("age", "sex")])

heat_colors <- list(
    sex = c("M" = "#01bfc5", "F" = "salmon"))

Simplified version of the heatmap code:

pheatmap(mat = assay(vsd[genes, ]),
    cluster_rows = TRUE,
    cluster_cols = FALSE,
    annotation_col = heat_anno,
    labels_col = coldata$age)

The result looks like that:

Pheatmap Output

The colors on top of the heatmap are correctly assigned to respective values. However, as you can see in my colnames, the legend for age incorrectly displays the min and max values as 2 and 8 instead of my 0.6 - 9.8. I want to set the legend annotation to 0 and 10 (or at the very least to min and max of my coldata$age). Does anybody know how I could do this? breaks and legend_breaks modifies the color gradient of the heatmap values themselves (here scaled expression), but not additional annotations.

Upvotes: 0

Views: 1441

Answers (1)

Marco Sandri
Marco Sandri

Reputation: 24272

Here is a version of the function pheatmap:::draw_annotation_legend that I have modified.
Copy this function into your script (or into a funct.r file and then call the file with the command source("funct.r")).

draw_annotation_legend <- function (annotation, annotation_colors, border_color, 
                          nbreaks=6, txt="minmax", cex.txt=.9, ...) {
    y = grid:::unit(1, "npc")
    text_height = grid:::unit(1, "grobheight", textGrob("FGH", gp = gpar(...)))
    res = gList()
    for (i in names(annotation)) {
        res[[i]] = textGrob(i, x = 0, y = y, vjust = 1, hjust = 0, 
            gp = gpar(fontface = "bold", ...))
        y = y - 1.5 * text_height
        if (is.character(annotation[[i]]) | is.factor(annotation[[i]])) {
            n = length(annotation_colors[[i]])
            yy = y - (1:n - 1) * 2 * text_height
            res[[paste(i, "r")]] = rectGrob(x = unit(0, "npc"), 
                y = yy, hjust = 0, vjust = 1, height = 2 * text_height, 
                width = 2 * text_height, gp = gpar(col = border_color, 
                  fill = annotation_colors[[i]]))
            res[[paste(i, "t")]] = textGrob(names(annotation_colors[[i]]), 
                x = text_height * 2.4, y = yy - text_height, 
                hjust = 0, vjust = 0.5, gp = gpar(...))
            y = y - n * 2 * text_height
        }
        else {

            # Modified: number of breaks
            yy = y - 8 * text_height + seq(0, 1, length.out=nbreaks+1)[-1] * 
                8 * text_height
            h = 8 * text_height * 1/nbreaks
            res[[paste(i, "r")]] = rectGrob(x = unit(0, "npc"), 
                y = yy, hjust = 0, vjust = 1, height = h, width = 2 * 
                  text_height, gp = gpar(col = NA, 
                  fill = colorRampPalette(annotation_colors[[i]])(nbreaks)))
            res[[paste(i, "r2")]] = rectGrob(x = unit(0, "npc"), 
                y = y, hjust = 0, vjust = 1, height = 8 * text_height, 
                width = 2 * text_height, gp = gpar(col = border_color, 
                  fill = NA))

            # Modified: how to handle colormap labels
            if (is.null(txt)) {
              txt = rev(range(grid.pretty(range(annotation[[i]], 
                  na.rm = TRUE))))
            } else if (txt=="minmax") {
              txt = rev(range(annotation[[i]], na.rm = TRUE))
            } else {
              txt = rev(txt)
            }

            # Modified: label position
            xx <- grid:::unit(1, "npc")*.6
            lab_height <- grid:::unit(1, "grobheight", 
                             textGrob("FGH", gp=gpar(cex=1/cex.txt, ...)))
            yy = grid:::unit(1, "npc") - 1.5 * lab_height - c(0.5, 7.5) * lab_height
            res[[paste(i, "t")]] = textGrob(txt, x = xx, y = yy, 
                                   hjust = 0, vjust = 0.5, gp = gpar(cex=cex.txt,...))

            y = y - 8 * text_height
        }
        y = y - 1.5 * text_height
    }
    res = gTree(children = res)
    return(res)
}

I have added three new parameters to this function:

  • nbreaks, which defines the number of breaks in the color gradient of the heatmap for the annotation legend
  • txt, which determines the labels to be placed at the extremes of the heatmap scale; txt="minmax" sets the labels at the minimum and maximum of the numeric variable (in your example, "age"); txt=NULL uses the labels defined by grid.pretty; you can also define two labels of your choice with txt=c("lab1","lab2")
  • cex.txt, which sets the magnification to be used for the labels.

To minimize modifications to pheatmap, I opted to have these three options defined directly in the defaults of the function draw_annotation_legend. In other words, if you want to modify these options, change their default values in the argument list of the function. It's not an elegant solution, but it's simple from a practical standpoint.
Here is an example demonstrating the use of the modified function with data that I have generated.

mtx <- structure(c(0.73, 0.89, -0.14, -0.57, 0.76, -2.27, 0.79, -0.35, 
-0.36, -1.15, -0.15, 2.27, 0.46, 0.65, -0.94, 1.02, -1.11, -0.41, 
1.4, 0.37, 0.97, 1.82, -0.81, -1.94, -2, 2.26, -0.6, 0.78, 0.77, 
-0.2, 1.01, 2.75, 2.56, 2.04, 0.32, 0.61, -0.41, -2.08, 2.21, 
0.03, 1.41, -2.98, -1.33, 1.17, 1.07, 1.83, -1.77, 0.71, 0.73, 
-1.63, -0.68, 2.43, 0.07, 0.44, -0.84, 0.35, 0.86, 1.03, 2.68, 
-2.93, 0.19, -1.68, 0.69, 1.99, -0.73, -2.29, 1.11, 1.99, 0.65, 
-1.62, 0.07, -0.98, -1.31, 2.91, 1.75, 1.18, 1.03, -1.01, 0.6, 
1.28, 0.81, 1.3, -0.38, 3.1, 1.21, 2.33, 0.84, -0.38, 0.67, 1.03, 
-1.2, -1.07, 2.36, -0.49, -1.23, 0.86, -0.63, 2.7, -0.75, -0.87, 
0.28, -1.45, 0.53, -1.66, 0.18, -0.67, -0.39, 1.95, -0.56, 0.4, 
-1.54, -1.66, 1.58, 1.65, -0.1, -0.63, -0.07, 0.79, 2.73, -0.09, 
1.93, 1.65, 0.4, 1.91, -0.53, -1.45, -2.31, 1.45, -2.65, -1.5, 
2.05, 1.1, 0.66, -1.48, 3.32, -1.31, -1.26, 0.84, 0.16, -0.53, 
-1.43, -1.62, -0.74, -1.88, 0.02, 0.68, -1.93, 1.06, 1.12, 0.17, 
-2.02, 0.69, 0.24, -1.01, -0.14, -0.31, 2.12, -0.43, 0.08, -0.81, 
-0.61, 0.38, -0.3, -0.6, -1.24, -0.35, 0.79, -1.55, 2.21, -0.03, 
0.25, 1.68, 0.05, 1.03, -2.13, 0.6, 3.1, 0.5, 0.27, -2.27, -1.14, 
-0.06, -0.51, 1.41, 1.02, 0.1, 1.82, 0.47, -0.21, -0.63, 0.68, 
-0.63, 0.98, 0.38, 1.64, 1, 1.06, -0.55, -0.56, 0.02, -1.8, -0.79, 
0.3, 1.32, 1.04, 0.13, -2.18, 0.81, 0.12, -0.1, 1.24, -1.07, 
-0.35, 2.58, -0.76, 0.39, 0.83, -2.15, -2.67, 0.09, 1.08, -2.86, 
-0.19, -0.34, 2.24, 0.84, -0.26, 0.02, 1.92, 0.1, 0.1, -0.97, 
0.21, 0.33, 1.11, -0.58, 0.95, -0.8, 0.78, 0.31, -0.88, -0.71, 
-0.33, -1.33, -0.13, 0.96, 3.43, -0.1, 0.68, 0.94, -1.01, 1.25, 
0.57, -1.79, -0.33, 0.8, -0.52, -0.57, -0.99, -1.45, 0.89, 1.58, 
-0.18, -0.64, 1.85, -0.2, 0.05, 0.6, -1.48, -0.83, -0.79, -0.88, 
0.72, -2.64, 0.33, 1.43, 0.02, -0.39, -1.2, 0.59, -1.89, 0.21, 
-1.09, 1.99, 0.81, 0.45, 0.13, -0.84, 1.22, 0.94, -0.54, -0.89, 
-0.24, 0.5, -1.22, 0.23, -2.69, -0.78, -0.96, 0.58, 0.65, 0.01, 
-0.55, 1.5, -0.15, 0.05, 1.49, 0.43, -0.41, 2.09, -1.15, -0.11, 
1.65, 2.16, 2.7, -0.39, -0.72, -1.76, 2.83, -0.96, 0.48, 0.76, 
1.27, 0.59, -2.73, 1.17, 0.59, 0.49, -0.91, 1.23, 1.78, 0.61, 
0.44, 1.08, 0.51, 0.46, -1.9, 1.94, 0.62, 0.58, 2.6, -0.38, 1.19, 
0.67, -0.12, -0.18, -1.48, 0.68, -3.23, 0.97, 0.37, -0.11, -1.83, 
-1.35, 1.32, -0.45, 0.44, 0.04, 0.59, -1.15, -0.47, -2.27, 0.36, 
-0.24), dim = c(7L, 52L), dimnames = list(NULL, c("1", "2", "3", 
"4", "5", "6", "7", "8", "9", "10", "11", "12", "13", "14", "15", 
"16", "17", "18", "19", "20", "21", "22", "23", "24", "25", "26", 
"27", "28", "29", "30", "31", "32", "33", "34", "35", "36", "37", 
"38", "39", "40", "41", "42", "43", "44", "45", "46", "47", "48", 
"49", "50", "51", "52")))

heat_anno <- structure(list(age = c(0.6, 0.6, 0.6, 0.8, 0.8, 0.8, 1, 1, 1, 
1.3, 1.3, 1.3, 2.1, 2.1, 2.1, 2.3, 2.3, 2.5, 2.5, 2.8, 2.8, 2.8, 
3.3, 3.3, 3.6, 3.6, 3.6, 3.8, 3.8, 3.8, 4.2, 4.2, 4.2, 4.8, 4.8, 
5.2, 5.2, 5.2, 5.6, 5.6, 6.6, 6.6, 6.6, 7.5, 7.5, 8, 8, 8, 8.6, 
9.8, 9.8, 9.8), sex = c("M", "M", "M", "F", "F", "F", "F", "M", 
"M", "F", "F", "M", "M", "M", "F", "F", "M", "F", "M", "M", "M", 
"F", "F", "M", "M", "M", "M", "M", "M", "M", "F", "F", "F", "M", 
"F", "F", "F", "F", "F", "F", "F", "F", "F", "M", "F", "M", "M", 
"F", "F", "F", "F", "F")), class = "data.frame", row.names = c(NA, 
52L))

heat_colors <- list(sex = c("M" = "#01bfc5", "F" = "salmon"))

library(pheatmap)
library(gtable)
library(grid)
assignInNamespace(x="draw_annotation_legend", 
                  value=draw_annotation_legend, ns="pheatmap")

pheatmap(mat = mtx,
    cluster_rows = TRUE,
    cluster_cols = FALSE,
    annotation_col = heat_anno[,c("sex", "age")],
    annotation_color=heat_colors)

enter image description here

Upvotes: 1

Related Questions