查看: 2249|回复: 15

[R语言] 如何提取热图聚类簇中的基因?

  [复制链接]

管理员

Rank: 15Rank: 15Rank: 15Rank: 15Rank: 15

主题
739
注册时间
2020.6.16
在线时间
412 小时

发表于 2021.12.13 10:08:41 | 显示全部楼层 |阅读模式
上周在微信推文《如何为热图添加指定的基因标记?》为大家介绍了如何为热图指定显示特定的基因标签,如下图,内容比较受欢迎。


刚解决热图指定标签问题,可是又有小伙伴在生信交流QQ群(群号:133146835)提问如何提取热图聚类簇中的基因。废话不多说,下面就直接为大家介绍两种常用的提取热图基因簇的方法。

1.pheatmap包绘制的热图

[AppleScript] 纯文本查看 复制代码
#载入pheatmap包;
library(pheatmap)
#设置工作目录;
setwd("C:/Users/86136/Desktop/获取热图不同簇的基因")

#读入测试数据,这里有两个,先读入包含50个基因的表达量文件用于测试;
df1 <- read.csv("mat_50genes4.csv",header = T,row.names = 1)
#查看数据前6行;
head(df1)



[AppleScript] 纯文本查看 复制代码
#自定义渐变色;
mycol = colorRampPalette(c("purple", "white","orange"))(100)
#创建分组颜色条所需的数据框;
Group = factor(c(rep("TESA",3), rep("TESB",3), rep("TESC",3),rep("TESD",3)))
anndf <- data.frame(Group)
rownames(anndf) <- colnames(df1)
#自定义分组颜色条的颜色;
anncol = list(Group=c(TESA="#c77cff",TESB="#FF9999",TESC="#99CC00",TESD="#FF9900"))

#绘制热图;
#对热图中的文字字号,聚类树高度,渐变颜色条做了自定义;
#同时,根据聚类结果在行和列方向都添加了“gap”;
htmap<- pheatmap(df1, scale = "row",
border="white",
cellwidth = 18,
cellheight = 6,
color = mycol,
cutree_rows=4,
cutree_cols=4,
fontsize = 7,
fontsize_col = 6,
fontsize_row = 6,
treeheight_row=30,
treeheight_col=30,
annotation_col=anndf,
annotation_colors=anncol,
annotation_legend=F,
legend_breaks=c(-2.5,-1.5,0,1.5,2.5),
legend_labels=c("-2.5","-1.5","0","1.5","2.5"),
angle_col = 90)



[AppleScript] 纯文本查看 复制代码
#按照热图的聚类顺序获取全部基因列表;
order <- htmap$tree_row$order
genelist <- row.names(df1)[order]
#查看前8个;
genelist[1:8]



[AppleScript] 纯文本查看 复制代码
#提取热图的行方向(基因)的聚类树;
cluster <- htmap$tree_row
#绘制聚类树;
plot(cluster,hang = -1,cex=0.6,axes=FALSE,ann=FALSE)



[AppleScript] 纯文本查看 复制代码
#对聚类树进行分簇;
cut <- cutree(cluster,4)
cut
#然后,获取第4簇的genes,注意第4簇位于热图的上方;
(genes<- names(cut)[cut==4])
#提取相应的表达量数据;
(genesdf <- df1[genes,])



[AppleScript] 纯文本查看 复制代码
#导出到本地;
write.csv(genesdf,file = "genesdf_4.csv")


2.使用ComplexHeatmap包绘制的热图

[AppleScript] 纯文本查看 复制代码
#载入ComplexHeatmap和circlize包;
library(ComplexHeatmap)
library(circlize)
#读入测试数据,这里有两个,先读入包含50个基因的表达量文件用于测试;
df1 <- read.csv("mat_50genes4.csv",header = T,row.names = 1)
#对数据进行归一化;
#由于scale函数默认对列进行归一化,因此这里做了两次转置;
df_scaled <- t(scale(t(df1)))
#查看归一化后的数据前6行;
head(df_scaled)



[AppleScript] 纯文本查看 复制代码
#计算距离矩阵,dist默认对矩阵的行进行计算;
d <- dist(df_scaled,method = "euclidean")
#然后进行分层聚类(Hierarchical cluster);
cluster2 <- hclust(d, method = "complete")
#绘制聚类树;
plot(cluster2,hang = -1,cex=0.6,axes=FALSE,ann=FALSE)



[AppleScript] 纯文本查看 复制代码
#同样的方法,对聚类树进行分簇;
cut <- cutree(cluster2,4)
cut
#然后,获取第4簇的genes,注意第4簇位于热图的上方;
(genes<- names(cut)[cut==4])
#提取相应的表达量数据;
(genesdf <- df1[genes,])



[AppleScript] 纯文本查看 复制代码
#导出到本地;
write.csv(genesdf,file = "genesdf_4.csv")

#接下来获得数据的范围;
range(df_scaled)
#[1] -1.532186 2.878927
#然后,根据数据范围建立自定义颜色映射关系;
#green-red;
col_fun = colorRamp2(c(-1.6,0.6,2.8), c("greenyellow","white", "red"))
#准备分组颜色条数据;
class = anno_block(gp = gpar(fill = c("#c77cff","#FF9999","#99CC00","#FF9900"),
col="white"),height = unit(5, "mm"),
labels = c("TESA", "TESB", "TESD","TESC"),
labels_gp = gpar(col = "white", fontsize = 8,fontface="bold"))
group= HeatmapAnnotation(group=class)

#绘制聚类热图,注意,Heatmap函数其实可以直接指定现有的聚类树的;
Heatmap(df_scaled,
rect_gp = gpar(col="white"),
row_names_gp = gpar(fontsize = 6),
column_names_gp = gpar(fontsize = 8),
col = col_fun,
column_split = 4,
row_split = 4,
column_title = NULL,
row_title = NULL,
cluster_rows = TRUE,
show_row_dend = TRUE,
top_annotation =group,
name = "Exp")



[AppleScript] 纯文本查看 复制代码
#当基因数较多时可以隐藏行方向的聚类树;
Heatmap(df_scaled,
rect_gp = gpar(col="white"),
row_names_gp = gpar(fontsize = 6),
column_names_gp = gpar(fontsize = 8),
col = col_fun,
column_split = 4,
row_split = 4,
column_title = NULL,
row_title = NULL,
cluster_rows = TRUE,
show_row_dend = FALSE,
top_annotation =group,
name = "Exp")



可见, ComplexHeatmap包和pheatmap包完全一样,默认都是使用hclust()函数进行聚类。因此,二者提取基因簇的方法也基本相同。好啦,今天的内容就分享到这里啦~


本文作者:基迪奥-莫北

本帖子中包含更多资源

您需要 登录 才可以下载或查看,没有帐号?立即注册

x
新的一天加油!
回复

使用道具 举报

迅猛龙

Rank: 8Rank: 8

主题
0
注册时间
2017.6.6
在线时间
182 小时

灌水之王


发表于 2021.12.13 14:24:13 | 显示全部楼层
新的一天加油!
回复

使用道具 举报

迅猛龙

Rank: 8Rank: 8

主题
0
注册时间
2021.4.21
在线时间
33 小时

发表于 2021.12.13 18:04:01 | 显示全部楼层
666666666
回复

使用道具 举报

迅猛龙

Rank: 8Rank: 8

主题
3
注册时间
2021.6.22
在线时间
57 小时

发表于 2021.12.13 21:46:26 | 显示全部楼层
学会了
回复

使用道具 举报

功夫熊猫

Rank: 10Rank: 10Rank: 10

主题
0
注册时间
2020.11.21
在线时间
109 小时

发表于 2021.12.13 23:11:16 | 显示全部楼层
新的一天加油!
回复

使用道具 举报

中华鲟

Rank: 5Rank: 5

主题
1
注册时间
2016.8.25
在线时间
98 小时

发表于 2021.12.13 23:30:59 | 显示全部楼层
新的一天加油!
回复

使用道具 举报

功夫熊猫

Rank: 10Rank: 10Rank: 10

主题
3
注册时间
2017.9.8
在线时间
83 小时

发表于 2021.12.14 08:10:21 | 显示全部楼层
新的一天加油!
回复

使用道具 举报

迅猛龙

Rank: 8Rank: 8

主题
0
注册时间
2017.6.6
在线时间
182 小时

灌水之王


发表于 2021.12.14 08:24:23 | 显示全部楼层
新的一天加油!
回复

使用道具 举报

迅猛龙

Rank: 8Rank: 8

主题
31
注册时间
2016.1.8
在线时间
621 小时

发表于 2021.12.14 08:36:59 | 显示全部楼层
如此好贴,必须点赞
新的一天加油!
回复 支持 反对

使用道具 举报

迅猛龙

Rank: 8Rank: 8

主题
0
注册时间
2021.1.28
在线时间
32 小时

发表于 2021.12.14 09:33:57 | 显示全部楼层
强!!!!!!
回复

使用道具 举报

钵水母

Rank: 3Rank: 3

主题
6
注册时间
2017.4.8
在线时间
36 小时

发表于 2021.12.14 09:48:58 | 显示全部楼层
如此好贴,必须点赞
回复 支持 反对

使用道具 举报

帝王蝶

Rank: 4

主题
0
注册时间
2019.3.26
在线时间
7 小时

发表于 2021.12.14 14:34:57 | 显示全部楼层
6666666666
新的一天加油!
回复 支持 反对

使用道具 举报

功夫熊猫

Rank: 10Rank: 10Rank: 10

主题
0
注册时间
2020.11.21
在线时间
109 小时

发表于 2021.12.14 16:29:17 | 显示全部楼层
新的一天加油!
回复

使用道具 举报

功夫熊猫

Rank: 10Rank: 10Rank: 10

主题
3
注册时间
2017.9.8
在线时间
83 小时

发表于 2021.12.15 08:10:10 | 显示全部楼层
新的一天加油!
回复

使用道具 举报

功夫熊猫

Rank: 10Rank: 10Rank: 10

主题
3
注册时间
2017.9.8
在线时间
83 小时

发表于 2021.12.16 08:04:41 | 显示全部楼层
新的一天加油!
回复

使用道具 举报

功夫熊猫

Rank: 10Rank: 10Rank: 10

主题
3
注册时间
2017.9.8
在线时间
83 小时

发表于 2021.12.17 08:15:20 | 显示全部楼层
新的一天加油!
回复

使用道具 举报

您需要登录后才可以回帖 登录 | 立即注册

本版积分规则

快速回复 返回顶部 返回列表