# 5.3.2 Seurat V3 | 如何改造Seurat包的DoHeatmap函数？

> 分析过单细胞数据的小伙伴应该都使用过Seurat包，其中有个函数叫`DoHeatmap`，具体操作可以看： [单细胞转录组学习笔记-17-用Seurat包分析文章数据](https://www.jianshu.com/p/f6f54ce92e24)

## 前言

走完Seurat流程，会得到分群结果`FindClusters()`，并找到marker基因`FindAllMarkers()`，然后想要对每群的前10个marker基因进行热图可视化

```r
rm(list = ls()) 
options(warn=-1) 
options(stringsAsFactors = F)

install.packages('Seurat')
library(Seurat)
library(stringr)   
library(dplyr)  

load('sce_out_for_heatmap_all.Rdata')
top10 <- sce.markers %>% group_by(cluster) %>% top_n(10, avg_logFC)
DoHeatmap(sce,top10$gene,size=3)
```

![plot\_zoom\_png](https://jieandze1314-1255603621.cos.ap-guangzhou.myqcloud.com/blog/2019-12-04-071009.png)

但是这个图在后期调整时会遇到很多障碍，因此最好用pheatmap重新画一下

## 如何用Pheatmap画这个结果？

其实，**画一个热图最重要的是表达矩阵和分组信息**

### **第1步：得到表达矩阵**

那么现在有了sce对象，从中提取表达矩阵也不难

```r
# 提取原始表达矩阵
cts <- GetAssayData(sce, slot = "counts")
> cts[1:4,1:4]
4 x 4 sparse Matrix of class "dgCMatrix"
               A1_01 A1_10 A1_11 A1_12
00R-AC107638.2      .      2      .      2
0610005C13Rik       .      .      .      .
0610007P14Rik       1     49      3    328
0610009B22Rik       .     12      .     78
# 然后对这个矩阵取log
cts <- log10(cts + 1)
```

### **第2步：得到小的top10表达矩阵**

当然不能使用整个表达矩阵进行处理，可以直接使用Seurat得到的`top10`结果，它帮我们得到了每个cluster的marker基因，而我们**只需要取出这个小表达矩阵即可**

Seurat得到的top10计算结果是这样，那么**矩阵的行就按基因名取：**

![](https://jieandze1314-1255603621.cos.ap-guangzhou.myqcloud.com/blog/2019-12-04-072247.png)

因为这个热图结果是按照cluster进行排序展示的，因此我们也要**将小表达矩阵的列按cluster从小到大排序：**

```r
# 原来的cluster分组信息存储在 sce$seurat_clusters 中
> head(sce$seurat_clusters)
A1_01 A1_10 A1_11 A1_12 A1_13 A1_14 
     4      4      4      4      4      4 
Levels: 0 1 2 3 4 5 6 7 8
# 可见并不是从第0个cluster开始的

# 排序之后
new_cluster <- sort(sce$seurat_clusters)
> head(new_cluster)
A10_09 A10_16 A10_18 A10_33 A10_36 A10_42 
      0       0       0       0       0       0 
Levels: 0 1 2 3 4 5 6 7 8
```

👌有了行和列的规定，我们就能很轻松地提取出整个小的表达矩阵：

```r
cts <- as.matrix(cts[top10$gene, names(new_cluster)])
```

### **第3步：做一个列的注释**

因为要做出来DoHeatmap的顶部0-8 cluster的展示，需要使用pheatmap的一个参数：`annotation_col`

这个参数接收一个数据框作为输入。因为是对列进行注释，所以这个数据框的行是矩阵的列名，而它的列在这里对应的就是cluster分群信息

```r
ac=data.frame(cluster=new_cluster)
rownames(ac)=colnames(mat)
> head(ac)
        cluster
A10_09       0
A10_16       0
A10_18       0
A10_33       0
A10_36       0
A10_42       0
```

最后，就可以画图了：

```r
library(pheatmap)
pheatmap(cts,show_colnames =F,show_rownames = T,
         cluster_rows = F,
         cluster_cols = F,
         annotation_col=ac)
```

![plot\_zoom\_png-2](https://jieandze1314-1255603621.cos.ap-guangzhou.myqcloud.com/blog/2019-12-04-073859.png)

当然，这个pheatmap有很多参数可以调整，看它的参数就知道：

比如想加上每个cluster的分隔，需要用到参数`gaps_col =` ，不过需要提供每个cluster最后一个细胞的序号

【当然这个是后话了，最重要的是知道如何进行数据的转换】

![](https://jieandze1314-1255603621.cos.ap-guangzhou.myqcloud.com/blog/2019-12-04-073958.png)

### 除了pheatmap，当然还能用其他的包

画热图的包有很多，其中一个比较常用的就是`ComplexHeatmap`

例如：

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

# 列标题的颜色框
color <- rainbow(9)
names(color) <- levels(new_cluster)
top_color <- HeatmapAnnotation(
  cluster = anno_block(gp = gpar(fill = color), 
                       labels = levels(new_cluster), 
                       labels_gp = gpar(cex = 0.5, col = "white"))) 

Heatmap(mat,
        cluster_rows = FALSE,
        cluster_columns = FALSE,
        show_column_names = FALSE,
        show_row_names = TRUE,
        column_split = new_cluster,
        heatmap_legend_param = list(
          title = "log10(count+1)",
          title_position = "leftcenter-rot"
        ),
        top_annotation = top_anno,
        column_title = NULL)
```

也可以实现列的注释，并且将每个cluster的列都进行分隔，和`DoHeatmap`的结果更像

![plot\_zoom\_png-3](https://jieandze1314-1255603621.cos.ap-guangzhou.myqcloud.com/blog/2019-12-04-075333.png)


---

# Agent Instructions: Querying This Documentation

If you need additional information that is not directly available in this page, you can query the documentation dynamically by asking a question.

Perform an HTTP GET request on the current page URL with the `ask` query parameter:

```
GET https://jieandze1314.osca.top/05/5.3-seurat-de-shi-yong/5.3.3-seurat-v3-ru-he-gai-zao-seurat-bao-de-doheatmap-han-shu.md?ask=<question>
```

The question should be specific, self-contained, and written in natural language.
The response will contain a direct answer to the question and relevant excerpts and sources from the documentation.

Use this mechanism when the answer is not explicitly present in the current page, you need clarification or additional context, or you want to retrieve related documentation sections.
