「熱圖」ComplexHeatmap展示單細胞聚類

實用Seurat自帶的熱圖函數(shù)DoHeatmap繪制的熱圖，感覺有點不上檔次，于是我嘗試使用ComplexHeatmap這個R包來對結(jié)果進行展示。

個人覺得好的熱圖有三個要素

  • 聚類: 能夠讓別人一眼就看到模式
  • 注釋: 附加注釋能提供更多信息
  • 配色: 要符合直覺，比如說大部分都會認為紅色是高表達，藍色是低表達

在正式開始之前，我們需要先獲取一下pbmc的數(shù)據(jù)，Seurat提供了R包SeuratData專門用于獲取數(shù)據(jù)

devtools::install_github('satijalab/seurat-data')
library(SeuratData)
InstallData("pbmc3k")

加載數(shù)據(jù)并進行數(shù)據(jù)預處理，獲取繪制熱圖所需的數(shù)據(jù)

library(SeuratData)
library(Seurat)
data("pbmc3k")
pbmc <- pbmc3k
pbmc[["percent.mt"]] <- PercentageFeatureSet(pbmc, pattern = "^MT-")
pbmc <- subset(pbmc, subset = nFeature_RNA > 200 & nFeature_RNA < 2500 & percent.mt < 5)
pbmc <- NormalizeData(pbmc, normalization.method = "LogNormalize", scale.factor = 10000)
pbmc <- FindVariableFeatures(pbmc, selection.method = "vst", nfeatures = 2000)
all.genes <- rownames(pbmc)
pbmc <- ScaleData(pbmc, features = all.genes)
pbmc <- RunPCA(pbmc, features = VariableFeatures(object = pbmc))
pbmc <- FindNeighbors(pbmc, dims = 1:10)
pbmc <- FindClusters(pbmc, resolution = 0.5)

pbmc.markers <- FindAllMarkers(pbmc, 
                               only.pos = TRUE, 
                               min.pct = 0.25, 
                               logfc.threshold = 0.25)

先感受下Seurat自帶熱圖

top10 <- pbmc.markers %>% group_by(cluster) %>% top_n(n = 10, wt = avg_logFC)
DoHeatmap(pbmc, features = top10$gene) + NoLegend()
Seurat-heatmap

下面則是介紹如何用R包ComplexHeatmap進行組圖，雖然這個R包名帶著Complex，但是并不是說這個R包很復雜，這個Complex應(yīng)該翻譯成復合，也就是說這個R包能在熱圖的基礎(chǔ)上整合很多信息。

先安裝并加載R包。

BiocManager::install("ComplexHeatmap")
library(ComplexHeatmap)

為了手動繪制一個熱圖，要從Seurat對象中提取所需要的表達量矩陣。我提取的是原始的count值，然后用log2(count + 1)的方式進行標準化

mat <- GetAssayData(pbmc, slot = "counts")
mat <- log2(mat + 1)

獲取基因和細胞聚類信息

gene_features <- top10
cluster_info <- sort(pbmc$seurat_annotations)

對表達量矩陣進行排序和篩選

mat <- as.matrix(mat[top10$gene, names(cluster_info)])

用Heatmap繪制熱圖。對于單細胞這種數(shù)據(jù)，一定要設(shè)置如下4個參數(shù)

  • cluster_rows= FALSE: 不作行聚類
  • cluster_columns= FALSE: 不作列聚類
  • show_column_names=FALSE: 不展示列名
  • show_row_names=FALSE: 不展示行名，基因數(shù)目不多時候可以考慮設(shè)置為TRUE
Heatmap(mat,
        cluster_rows = FALSE,
        cluster_columns = FALSE,
        show_column_names = FALSE,
        show_row_names = TRUE)
Heatmap-1

從圖中，我們可以發(fā)現(xiàn)以下幾個問題：

  • 長寬比不合理，當然這和繪圖函數(shù)無關(guān)，可以在保存時修改長寬比
  • 基因名重疊，考慮調(diào)整大小，或者不展示，或者只展示重要的基因
  • 顏色可以調(diào)整
  • 缺少聚類信息

這些問題，我們可以通過在ComplexHeatmap Complete Reference查找對應(yīng)信息來解決。

配色方案

在熱圖中會涉及到兩類配色，一種用來表示表達量的連續(xù)性變化，一種則是展示聚類。有一個神奇的R包就是用于處理配色，他的Github地址為https://github.com/caleblareau/BuenColors。

devtools::install_github("caleblareau/BuenColors")
library("BuenColors")

它提供了一些列預設(shè)的顏色，比方說jdb_color_maps

      HSC       MPP      LMPP       CMP       CLP       MEP       GMP 
"#00441B" "#46A040" "#00AF99" "#FFC179" "#98D9E9" "#F6313E" "#FFA300" 
      pDC      mono     GMP-A     GMP-B     GMP-C       Ery       CD4 
"#C390D4" "#FF5A00" "#AFAFAF" "#7D7D7D" "#4B4B4B" "#8F1336" "#0081C9" 
      CD8        NK         B 
"#001588" "#490C65" "#BA7FD0"

這些顏色就能用于命名單細胞的類群，比如說我選擇了前9個

col <- jdb_color_maps[1:9]
names(col) <- levels(cluster_info)

增加列聚類信息

Heatmap的row_split和column_split參數(shù)可以通過設(shè)置分類變量對熱圖進行分隔。更多對熱圖進行拆分，可以參考Heatmap split

Heatmap(mat,
        cluster_rows = FALSE,
        cluster_columns = FALSE,
        show_column_names = FALSE,
        show_row_names = FALSE,
        column_split = cluster_info)
Heatmap-2

只用文字描述可能不夠好看，最好是帶有顏色的分塊圖，其中里面的顏色和t-SNE或UMAP聚類顏色一致，才能更好的展示信息。

為了增加聚類注釋，我們需要用到HeatmapAnnotation函數(shù)，它對細胞的列進行注釋，而rowAnnotation函數(shù)可以對行進行注釋。這兩個函數(shù)能夠增加各種類型的注釋，包括條形圖，點圖，折線圖，箱線圖，密度圖等等，這些函數(shù)的特征是anno_xxx，例如anno_block就用來繪制區(qū)塊圖。

top_anno <- HeatmapAnnotation(
  cluster = anno_block(gp = gpar(fill = col), # 設(shè)置填充色
                       labels = levels(cluster_info), 
                       labels_gp = gpar(cex = 0.5, col = "white"))) # 設(shè)置字體

其中anno_block中的gp參數(shù)用于設(shè)置各類圖形參數(shù)，labels設(shè)置標簽，labels_gp設(shè)置和標簽相關(guān)的圖形參數(shù)?？梢杂?code>?gp來了解有哪些圖形參數(shù)。

Heatmap(mat,
        cluster_rows = FALSE,
        cluster_columns = FALSE,
        show_column_names = FALSE,
        show_row_names = FALSE,
        column_split = cluster_info,
        top_annotation = top_anno, # 在熱圖上邊增加注釋
        column_title = NULL ) # 不需要列標題
Heatmap-3

突出重要基因

由于基因很多直接展示出來，根本看不清，我們可以強調(diào)幾個標記基因。用到兩個函數(shù)是rowAnnotation和anno_mark

已知不同類群的標記基因如下

Cluster ID Markers Cell Type
0 IL7R, CCR7 Naive CD4+ T
1 IL7R, S100A4 Memory CD4+
2 CD14, LYZ CD14+ Mono
3 MS4A1 B
4 CD8A CD8+ T
5 FCGR3A, MS4A7 FCGR3A+ Mono
6 GNLY, NKG7 NK
7 FCER1A, CST3 DC
8 PPBP Platelet

我們需要給anno_mark提供基因所在行即可。

mark_gene <- c("IL7R","CCR7","IL7R","S100A4","CD14","LYZ","MS4A1","CD8A","FCGR3A","MS4A7","GNLY","NKG7","FCER1A", "CST3","PPBP")
gene_pos <- which(rownames(mat) %in% mark_gene)

row_anno <-  rowAnnotation(mark_gene = anno_mark(at = gene_pos, 
                                     labels = mark_gene))

接著繪制熱圖

Heatmap(mat,
        cluster_rows = FALSE,
        cluster_columns = FALSE,
        show_column_names = FALSE,
        show_row_names = FALSE,
        column_split = cluster_info,
        top_annotation = top_anno,
        right_annotation = row_anno,
        column_title = NULL)
Heatmap-4

關(guān)于如何增加標記注釋，參考mark-annotation

調(diào)增圖例位置

目前的熱圖還有一個問題，也就是表示表達量范圍的圖例太占位置了，有兩種解決方法

  • 參數(shù)設(shè)置show_heatmap_legend=FALSE直接刪掉。
  • 利用heatmap_legend_param參數(shù)更改樣式

我們根據(jù)legends這一節(jié)的內(nèi)容進行一些調(diào)整

Heatmap(mat,
        cluster_rows = FALSE,
        cluster_columns = FALSE,
        show_column_names = FALSE,
        show_row_names = FALSE,
        column_split = cluster_info,
        top_annotation = top_anno,
        right_annotation = row_anno,
        column_title = NULL,
        heatmap_legend_param = list(
          title = "log2(count+1)",
          title_position = "leftcenter-rot"
        ))
heatmap-5

因為ComplextHeatmap是基于Grid圖形系統(tǒng)，因此可以先繪制熱圖，然后再用grid::draw繪制圖例，從而實現(xiàn)將條形圖的位置移動到圖中的任意位置。

先獲取繪制熱圖的對象

p <- Heatmap(mat,
        cluster_rows = FALSE,
        cluster_columns = FALSE,
        show_column_names = FALSE,
        show_row_names = FALSE,
        column_split = cluster_info,
        top_annotation = top_anno,
        right_annotation = row_anno,
        column_title = NULL,
        show_heatmap_legend = FALSE
        )

根據(jù)p@matrix_color_mapping獲取圖例的顏色的設(shè)置，然后用Legend構(gòu)建圖例

col_fun  <- circlize::colorRamp2(c(0, 1, 2 ,3, 4),
                                 c("#0000FFFF", "#9A70FBFF", "#D8C6F3FF", "#FFC8B9FF", "#FF7D5DFF"))
lgd <-  Legend(col_fun = col_fun, 
               title = "log2(count+1)", 
               title_gp = gpar(col="white", cex = 0.75),
               title_position = "leftcenter-rot",
               #direction = "horizontal"
               at = c(0, 1, 4), 
               labels = c("low", "median", "high"),
               labels_gp = gpar(col="white")
               )

繪制圖形

grid.newpage() #新建畫布
draw(p) # 繪制熱圖
draw(lgd, x = unit(0.05, "npc"), 
     y = unit(0.05, "npc"), 
     just = c("left", "bottom")) # 繪制圖形
heatmap-6

ComplexHeatmap繪制熱圖非常強大的工具，大部分我想要的功能它都有，甚至我沒有想到的它也有，這個教程只是展示其中一小部分功能而已，還有很多功能要慢慢探索。

版權(quán)聲明：本博客所有文章除特別聲明外，均采用 知識共享署名-非商業(yè)性使用-禁止演繹 4.0 國際許可協(xié)議 (CC BY-NC-ND 4.0) 進行許可。

掃碼即刻交流
?著作權(quán)歸作者所有,轉(zhuǎn)載或內(nèi)容合作請聯(lián)系作者
【社區(qū)內(nèi)容提示】社區(qū)部分內(nèi)容疑似由AI輔助生成，瀏覽時請結(jié)合常識與多方信息審慎甄別。
平臺聲明：文章內(nèi)容（如有圖片或視頻亦包括在內(nèi)）由作者上傳并發(fā)布，文章內(nèi)容僅代表作者本人觀點，簡書系信息發(fā)布平臺，僅提供信息存儲服務(wù)。

相關(guān)閱讀更多精彩內(nèi)容

友情鏈接更多精彩內(nèi)容