继续来看pheatmap那些有趣的事情

刘小泽写于2020.5.18
上次看了热图如何去掉聚类树的同时保留聚类的顺序 ,有朋友留言给出了其他的解决方法,非常感谢!文中会给出答案
另外,这次还是带着数据,去看看pheatmap其他的有用操作

上数据!

test = matrix(rnorm(200), 20, 10)
test[1:10, seq(1, 10, 2)] = test[1:10, seq(1, 10, 2)] + 3
test[11:20, seq(2, 10, 2)] = test[11:20, seq(2, 10, 2)] + 2
test[15:20, seq(2, 10, 2)] = test[15:20, seq(2, 10, 2)] + 4
# 20行(基因)10列(样本)
colnames(test) = paste("Test", 1:10, sep = "")
rownames(test) = paste("Gene", 1:20, sep = "")

# 和上次不同的是,我们这次设置一个极值(为了后面演示scale)
test[10,10]=100
> rownames(test)[10];colnames(test)[10]
[1] "Gene10"
[1] "Test10"
上数据

画图

第一个原始图

可以看到Gene10 Test10这里的100大大超过了一般水平,导致画的图受它的影响极大,变得看不出整体的变化趋势

pheatmap(test)
第一个原始图

第二个:增加z-score参数的图

设置一个参数scale = 'row' 就会对行进行z-score标准化

什么是z-score? 看个公式就明白了(表达量-均值)/标准差,目的就是不让数据受极值的影响太大。另外z-score一般还会设置一个最大最小值(例如-2 ~ 2-1 ~ 1
Z score = (expression G - mean expression G1..Gn)/SD G1..Gn where G is any gene in the pathway and G1..Gn represent the aggregate measure of all of the genes.
Cite:Cheadle C, Vawter MP, Freed WJ, Becker KG. Analysis of microarray data using Z score transformation. J Mol Diagn. 2003;5(2):73‐81. doi:10.1016/S1525-1578(10)60455-2

scale做了什么
pheatmap(test,scale = 'row')
增加z-score参数的图

第三个图:手动zscore

如果想要理解某个参数做了什么,就自己模拟这个参数吧

我们这里也来模拟一下这个参数做了什么【并且用热心朋友推荐的方法去掉聚类树】

cal_z_score <- function(x){
    (x - mean(x)) / sd(x)
}
diy1 <- t(apply(test, 1, cal_z_score))
# 很巧妙,设置高度为0,就不显示了
pheatmap(diy1,treeheight_row = 0,treeheight_col = 0)
手动zscore

第四个图:聚类是怎么聚的呢?

使用参数的话,很简单:
pheatmap(diy1,clustering_method = 'complete',cluster_rows = T, cluster_cols = T,
         clustering_distance_rows = "euclidean",
         clustering_distance_cols = "euclidean")

这个method是直接使用的hlust计算的,全部聚类方法包括:

"ward.D", "ward.D2", "single", "complete", "average" (= UPGMA), "mcquitty" (= WPGMA), "median" (= WPGMC) or "centroid" (= UPGMC).

距离使用dist计算,默认使用欧式距离,全部距离包括:

"euclidean", "maximum", "manhattan", "canberra", "binary", "minkowski"

提取聚类结果:

p1 <- pheatmap(diy1, silent = TRUE) # 只输出结果,不画图
> names(p1)
[1] "tree_row" "tree_col" "kmeans"   "gtable"  
# 做个丑陋的聚类图看看
library(dendextend)
p1$tree_row %>%
    as.dendrogram() %>%
    plot(horiz = TRUE)
# 聚类基因名
> rownames(diy1)[p1$tree_row[["order"]]]
 [1] "Gene1"  "Gene3"  "Gene2"  "Gene8"  "Gene4"  "Gene6"  "Gene9"  "Gene5"  "Gene7"  "Gene10" "Gene13" "Gene17"
[13] "Gene14" "Gene11" "Gene16" "Gene20" "Gene12" "Gene19" "Gene15" "Gene18"
自己聚类是这样的:
my_hclust_gene <- hclust(dist(diy1), method = "complete")

as.dendrogram(my_hclust_gene) %>%
    plot(horiz = TRUE)
自己聚类是这样的

第五个图:把树切分

pheatmap(diy1,cutree_rows = 2,cutree_cols = 2)
把树切分

第六张图:自己给切分后的基因上色

# cutree_rows = 2就做了下面这行代码的事情
gene_anno <- cutree(tree = as.dendrogram(my_hclust_gene), k = 2)

gene_anno <- data.frame(cluster = ifelse(gene_anno == 1, "cluster1", "cluster2"))

> table(gene_anno)
gene_anno
cluster 1 cluster 2 
        9        11 

pheatmap(diy1, annotation_row = gene_anno)
自己给切分后的基因上色
修改颜色

需要设置一个list,如果有列的注释信息,也是添加到这个my_colour的list中

my_colour = list(
    cluster = c(cluster1 = "#e89829", cluster2 = "#cc4ee0")
)

pheatmap(diy1,
         annotation_colors = my_colour,
         annotation_row = gene_anno,
         treeheight_row = 0)
修改颜色

欢迎关注我们的公众号~_~  
我们是两个农转生信的小硕,打造生信星球,想让它成为一个不拽术语、通俗易懂的生信知识平台。需要帮助或提出意见请后台留言或发送邮件到jieandze1314@gmail.com

Welcome to our bioinfoplanet!

最后编辑于
©著作权归作者所有,转载或内容合作请联系作者
  • 序言:七十年代末,一起剥皮案震惊了整个滨河市,随后出现的几起案子,更是在滨河造成了极大的恐慌,老刑警刘岩,带你破解...
    沈念sama阅读 151,688评论 1 330
  • 序言:滨河连续发生了三起死亡事件,死亡现场离奇诡异,居然都是意外死亡,警方通过查阅死者的电脑和手机,发现死者居然都...
    沈念sama阅读 64,559评论 1 273
  • 文/潘晓璐 我一进店门,熙熙楼的掌柜王于贵愁眉苦脸地迎上来,“玉大人,你说我怎么就摊上这事。” “怎么了?”我有些...
    开封第一讲书人阅读 101,749评论 0 226
  • 文/不坏的土叔 我叫张陵,是天一观的道长。 经常有香客问我,道长,这世上最难降的妖魔是什么? 我笑而不...
    开封第一讲书人阅读 42,581评论 0 191
  • 正文 为了忘掉前任,我火速办了婚礼,结果婚礼上,老公的妹妹穿的比我还像新娘。我一直安慰自己,他们只是感情好,可当我...
    茶点故事阅读 50,741评论 3 271
  • 文/花漫 我一把揭开白布。 她就那样静静地躺着,像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上,一...
    开封第一讲书人阅读 39,684评论 1 192
  • 那天,我揣着相机与录音,去河边找鬼。 笑死,一个胖子当着我的面吹牛,可吹牛的内容都是我干的。 我是一名探鬼主播,决...
    沈念sama阅读 31,122评论 2 292
  • 文/苍兰香墨 我猛地睁开眼,长吁一口气:“原来是场噩梦啊……” “哼!你这毒妇竟也来了?” 一声冷哼从身侧响起,我...
    开封第一讲书人阅读 29,847评论 0 182
  • 序言:老挝万荣一对情侣失踪,失踪者是张志新(化名)和其女友刘颖,没想到半个月后,有当地人在树林里发现了一具尸体,经...
    沈念sama阅读 33,441评论 0 228
  • 正文 独居荒郊野岭守林人离奇死亡,尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 29,939评论 2 232
  • 正文 我和宋清朗相恋三年,在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 31,333评论 1 242
  • 序言:一个原本活蹦乱跳的男人离奇死亡,死状恐怖,灵堂内的尸体忽然破棺而出,到底是诈尸还是另有隐情,我是刑警宁泽,带...
    沈念sama阅读 27,783评论 2 236
  • 正文 年R本政府宣布,位于F岛的核电站,受9级特大地震影响,放射性物质发生泄漏。R本人自食恶果不足惜,却给世界环境...
    茶点故事阅读 32,275评论 3 220
  • 文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹,春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 25,830评论 0 8
  • 文/苍兰香墨 我抬头看了看天上的太阳。三九已至,却和暖如春,着一层夹袄步出监牢的瞬间,已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 26,444评论 0 180
  • 我被黑心中介骗来泰国打工, 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留,地道东北人。 一个月前我还...
    沈念sama阅读 34,553评论 2 249
  • 正文 我出身青楼,却偏偏与公主长得像,于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子,可洞房花烛夜当晚...
    茶点故事阅读 34,618评论 2 249