巧用热图展示基因分布的总体趋势

Posted 2023-03-14

热图是最常见的基因表达量数据的可视化方式，将每个单元格的表达量按照数值高低映射为不同的颜色，可以直观展示表达量在不同样本间的分布，再综合聚类的结果和基因/样本的注释信息，进一步丰富了展示的信息，一个经典的热图如下

图中提供了两大类的信息，第一大部分也是热图的主体部分，即表达量信息，上图中，每一列表示样本，每一行表示基因，用不同颜色表征表达量的不同数值；第二部分为行或者列的注释信息，对应上图中顶部的样本注释信息，从图例可以看到，有3类注释信息。那么这样的一幅图如何来做？

可用的R包当然很多，我最常用的是pheatmap这个R包，在作图之前，先准备好数据，即表达量数据和样本的注释信息

# 1.基因表达量数据
# 纯文本文件，每一行为基因，每一列为样本
> data <- read.table("data.txt", header = T, sep = "\\t", row.names = 1)
> data[1:5, 1:5]
       GSM3076582 GSM3076584 GSM3076586 GSM3076588 GSM3076590
PNP      7.123107   7.115196   7.103920   7.377837   6.960771
PTGIS    2.483152   4.215764   4.174901   3.658576   3.296362
PTGS2    6.156447   6.539128   6.294466   4.918229   6.269206
NT5C1B   1.912044   3.691517   3.771510   4.546582   3.784036
SIRT3    4.210383   3.960023   3.623266   4.658751   4.347605


# 2.样本的注释信息
# 纯文本文件，每一行为样本，每一列为一种注释信息
> group <- read.table("sample.group.txt", header = T, sep = "\\t", row.names = 1, stringsAsFactors = T)
> head(group)
           group sex    series
GSM3076582   ALS   m GSE112676
GSM3076584   ALS   m GSE112676
GSM3076586   ALS   f GSE112676
GSM3076588   CON   m GSE112676
GSM3076590   CON   m GSE112676
GSM3076592   CON   f GSE112676

准备好文件之后，只需要调用pheatmap即可快速出图，接下来的工作就是熟悉这个方法的各项参数，通过调参来优化我们的可视化结果。我推荐的做法是做加法，首先只使用最少的参数，画出基本图形，然后再添加必要的参数，美化输出，这样可以更好的掌握每个选项的作用，基本用法如下

pheatmap(data)

通过观察输出，我们来确定优化的方向。表达量的基本信息已经有了，对照文献中的图，需要进行以下调整

1.调整颜色梯度，默认的颜色梯度和文献中的不同

2.调整行和列的文字大小，缩小基因名字使其不至于重叠，同时不展示样本名

3.调整聚类结果，不展示列的聚类结果

带着这几个目的，从选项中查找我们需要的选项，并设置合适的值，代码如下

> color <- colorRampPalette(c("navy", "white", "firebrick3"))(50)
> pheatmap(data,
+     color = color,    # 图例颜色
+     cluster_col = FALSE,    # 不显示样本聚类
+     show_colnames = FALSE,  # 不显示样本名称
+     fontsize_row = 10,     # 调整基因名称大小
+ )

调参后的结果如下

这里只展示了调整的方法，具体的颜色设置可以根据你的审美进行更换，其实一幅好看的图表，画图的代码并不是最难的，配色，图片的纵横比才是决定颜值的灵魂。

搞定了表达量之后，再来添加样本的注释信息，代码如下

> pheatmap(data,
+     color = color,
+     cluster_col = FALSE,  
+     show_colnames = FALSE,
+     annotation_col = group, # 添加列的注释
+     fontsize_row = 10
+ )

效果图如下

可以看到，只需要提供列注释的数据框，函数会自动进行映射，当然我们可以通过一下参数手动调整颜色

> ann_colors = list(
+     series = c(GSE112676 = "blue",  GSE112680 = "green"),
+     group  = c(ALS = "brown", CON = "green", MIM = "pink"),
+     sex    = c(f = "blue", m = "red")
+ )
>
> pheatmap(data,
+     color = color,
+     cluster_col = FALSE,  
+     show_colnames = FALSE,
+     annotation_col = group,
+     annotation_colors = ann_colors,
+     fontsize_row = 10
+ )
>

虽然直接从网站粘贴别人的代码，可能可以实现你想要的效果，但是更多的也只是人云亦云，只用通过基本用法加筛选参数的两步走战略，才能更好的帮助我们掌握每个选项的作用。