查看原文
其他

关于 scale 函数和 pheatmap 的图例问题

JunJunLab 老俊俊的生信笔记 2022-08-15

还没关注?快来点这里!?



问题


最近看见几个小伙伴在群里询问 pheatmap 绘制热图的相关问题:

  • 数据离散度太大,热图趋势看不出来?
  • 这个图例范围怎么调整?
  • scale 了呀,怎么热图还是对比不明显?


看法


1、数据离散度大,可以尝试将数据取 log2log10 来缩小数据之间的差异,比如像用的 FPKM/RPKM/TPM 等表达量数据,然后再做个 scale(Z-score)处理应该就可以了。

当然如果用的 log2FoldChange 差异倍数的数据画,就不用这么做了。

2、pheatmap 图例是根据数据的范围来调整的,如果想自己定义显示范围和标签也是可以的。legend_breaksbreakslegend_labels 三个参数互相配合调整即可。

3、如果你的数据已经 scale 了,但还是看不出 对比度 或者 表达趋势 ,说明可能有 异常值 存在,可以 设置图例范围 达到目的,或者 去除异常值



scale 函数详解

pheatmap 函数里有个 scale 参数,可以指定是按行 scale 还是按列 scale,即:scale = 'row'/scale = 'col'。一般情况我们是想比较同一个基因在不同样本中的相对表达高低,而不关心该基因在这个样本里相对其它基因的表达水平,指定按行 scale 即可。

R base 包里也有个 scale 函数,基本上原理是一样的:

scale(x, # 数据,数值型
      center = TRUE# 是否中心化,默认TRUE
      scale = TRUE # 是否标准化,默认TRUE
      )

所谓的 center 参数意思是中心化的意思,即每个数据减去数据的 均值 ,达到一个在 0 左右分布的数据。

scale 参数理解为是否标准化,有两种情况:

  • 1、当 center = TRUE 时,scaleTRUE,计算公式为:中心化的数据除以数据的标准差(sd)。
  • 2、当 center = FALSE 时,scaleTRUE,计算公式为:数据除以均方根(sqrt(sum(x^2)/(n-1)))。

我们具体实践看看,先模拟 center = TRUE,scale = TRUE 情况:

# center = TRUE,scale = TRUE
# 数据
d <- c(1:5)

# 减去均值
dm <- d - mean(d)
# 除以标准差
dm / sd(d)
[1] -1.2649111 -0.6324555  0.0000000  0.6324555  1.2649111
# scale
# 同时中心化和标准化
dcs <- scale(d,center = T,scale = T)
dcs
           [,1]
[1,] -1.2649111
[2,] -0.6324555
[3,]  0.0000000
[4,]  0.6324555
[5,]  1.2649111
attr(,"scaled:center")
[13
attr(,"scaled:scale")
[11.581139

可以看到我们手动算的和函数算的结果是一致的。

再看看 center = FALSE,scale = TRUE 情况:

# 数据除以均方根
st <- sqrt(sum(d^2)/(5-1))
d / st
[10.2696799 0.5393599 0.8090398 1.0787198 1.3483997
# 标准化
ds <- scale(d,center = F,scale = T)
ds
          [,1]
[1,] 0.2696799
[2,] 0.5393599
[3,] 0.8090398
[4,] 1.0787198
[5,] 1.3483997
attr(,"scaled:scale")
[13.708099

可以看到结果也是一样的。如果我们想在不中心化的数据里(center = F)之间用标准差(sd)来做标准化呢?你可以自己手动算:

# 手动算
d / sd(d)
[10.6324555 1.2649111 1.8973666 2.5298221 3.1622777

画个图看看这些区别:

# 绘图展示

# 中心化
dc <- scale(d,center = T,scale = F)
# 标准化
ds <- scale(d,center = F,scale = T)
# 同时中心化和标准化
dcs <- scale(d,center = T,scale = T)

# 整合数据
da <- data.frame(name = c(rep('d',5),rep('dc',5),rep('ds',5),rep('dcs',5)),
                 value = c(d,dc,ds,dcs),
                 y = c(rep(1,5),rep(2,5),rep(3,5),rep(4,5)))
# 绘图
ggplot(da,aes(x = value,y = y,fill = name)) +
  geom_point(size = 7,shape = 21) +
  theme_bw(base_size = 14) +
  theme(legend.position = c(0.9,0.8)) +
  scale_fill_brewer(palette = 'Set1') +
  scale_y_continuous(breaks = seq(0,5,1)) +
  scale_x_continuous(breaks = seq(-2,5,1)) +
  annotate('text',x = 3,y = 1.3,label = 'Row data',size = 6) +
  annotate('text',x = 0,y = 2.3,label = 'center = T',size = 6) +
  annotate('text',x = 0.8,y = 3.3,label = 'scale = T',size = 6) +
  annotate('text',x = 0,y = 4.3,label = 'scale = T & center = T',size = 6)

绘制热图展示:

这里有个小 tips ,就是使用 ggplotify 这个包把 pheatmap 对象变成 ggplot 对象用来拼图:

# 加载R包
library(ggplotify)
library(patchwork)

# Create test matrix
test = matrix(rnorm(200), 2010)
test[1:10, seq(1102)] = test[1:10, seq(1102)] + 3
test[11:20, seq(2102)] = test[11:20, seq(2102)] + 2
test[15:20, seq(2102)] = test[15:20, seq(2102)] + 4
colnames(test) = paste("Test"1:10, sep = "")
rownames(test) = paste("Gene"1:20, sep = "")

# 绘制热图
p1 <- as.ggplot(pheatmap(test,main = 'none scale'))
p2 <- as.ggplot(pheatmap(test, scale = "row",main = 'scale = row'))
p3 <- as.ggplot(pheatmap(test, scale = "column",main = 'scale = column'))

# 拼图
p1 + p2 + p3 + plot_annotation(tag_levels = LETTERS[1:3])
pheatmap 图例调整

上图第一个图有个较大值 8,通过调整图例范围来显示较小的范围的值的对比变化。

显示特定范围:

# 原始数据
p4 <- as.ggplot(pheatmap(test,
         border_color = 'black',
         main = 'raw',
         color = colorRampPalette(colors = c("navy""white""firebrick3"))(50)))

# 显示-4到4
p5 <- as.ggplot(pheatmap(test,
                         border_color = 'black',
                         main = '-4-4',
                         color = colorRampPalette(colors = c("navy""white""firebrick3"))(length(seq(-4,4,by = 0.1))),
                         breaks = seq(-4,4,by = 0.1),
                         legend_breaks = seq(-4,4,1)))
# 显示-8到8
p6 <- as.ggplot(pheatmap(test,
                   border_color = 'black',
                   main = '-8-8',
                   color = colorRampPalette(colors = c("navy""white""firebrick3"))(length(seq(-8,8,by = 0.1))),
                   breaks = seq(-8,8,by = 0.1),
                   legend_breaks = seq(-8,8,2)
                   ))

# 拼图
p4 + p5 + p6 + plot_annotation(tag_levels = LETTERS[1:3])

以上是对于原始数据的展示,对于 scale 后的数据同样可以:

# scale = 'row'
p7 <- as.ggplot(pheatmap(test,
                         border_color = 'black',
                         scale = 'row',
                         main = 'scale = row',
                         color = colorRampPalette(colors = c("navy""white""firebrick3"))(50)))
# -1-1
p8 <- as.ggplot(pheatmap(test,
                         border_color = 'black',
                         scale = 'row',
                         main = '-1-1',
                         breaks = seq(-1,1,by = 0.1),
                         legend_breaks = seq(-1,1,0.5),
                         color = colorRampPalette(colors = c("navy""white""firebrick3"))(length(seq(-1,1,by = 0.1)))))

# 拼图
p7 + p8 + plot_annotation(tag_levels = LETTERS[1:2])


欢迎加入生信交流群。加我微信我也拉你进 微信群聊 老俊俊生信交流群 哦,代码已上传至QQ群文件夹。

群二维码:


老俊俊微信:




知识星球:



所以今天你学习了吗?

欢迎小伙伴留言评论!

今天的分享就到这里了,敬请期待下一篇!

最后欢迎大家分享转发,您的点赞是对我的鼓励肯定

如果觉得对您帮助很大,赏杯快乐水喝喝吧!



 往期回顾 




Circular RNAs 的生物发生、功能和挑战

@你需要提高一下 R 技能了(plyr 包)

ggplot 图例(你想要的都在这了!)

把 corrplot 颜色条改成文献里那样?

ggcor 的环形热图

你看过 NCBI 的基因组和注释文件吗?

ComplexHeatmap 之 Legends 续(二)

ComplexHeatmap 之 Legends 续(一)

ComplexHeatmap 之 Legends

ComplexHeatmap 之 Heatmap List 续(二)

◀...

您可能也对以下帖子感兴趣

文章有问题?点此查看未经处理的缓存