生物信息R画图

盘一盘画相关性热图的几种方式

2020-10-11  本文已影响0人  凯凯何_Boy

R可视化相关性矩阵的几种方案

R中相关性矩阵的可视化解决方法在已经有很多了,我们在这里总结一些常用的都有哪些:

  1. ggplot2
  2. corrplot
  3. 热图包(pheatmapheatmap包等)
  4. ggally 包中的 ggcorr() 函数
  5. ggcorrplot

准备数据集

我们就使用R中内置数据集mtcars, 计算出矩阵的相关系数corr,后续出图都将是基于这个数据进行绘制。

data("mtcars")
cormat <- round(cor(mtcars), 1) # 计算相关系数,默认是pearson方法

#      mpg  cyl disp   hp drat   wt qsec   vs   am gear carb
#mpg   1.0 -0.9 -0.8 -0.8  0.7 -0.9  0.4  0.7  0.6  0.5 -0.6
#cyl  -0.9  1.0  0.9  0.8 -0.7  0.8 -0.6 -0.8 -0.5 -0.5  0.5
#disp -0.8  0.9  1.0  0.8 -0.7  0.9 -0.4 -0.7 -0.6 -0.6  0.4
#......

方式一: ggplot2

  1. 通过reshape2包将矩阵变换为长数据格式再进行绘制
library(reshape2)
melted_cormat <- melt(cormat)
head(melted_cormat)

##   Var1 Var2 value
## 1  mpg  mpg  1.00
## 2 disp  mpg -0.85
## 3   hp  mpg -0.78
## 4 drat  mpg  0.68
  1. 使用 geom_tile() 函数或者 geom_raster() 函数绘制相关性系数的图,再通过 coord_equal() 函数固定坐标轴的纵横比,简单出图瞅瞅。
library(ggplot2)
ggplot(data = melted_cormat, aes(x=Var1, y=Var2, fill=value))+geom_raster() +coord_equal()
图片.png
  1. 注意,相关性矩阵是有冗余信息的(矩阵的上三角和下三角信息表达一样的意思)。我们可以使用下面函数(参考)将它的一半设为NA,让它只输出上三角或者下三角矩阵的图形(分别利用 upper.trilower.tri 函数)。
# 得到上三角矩阵
  get_lower_tri<-function(cormat){
    cormat[upper.tri(cormat)] <- NA
    return(cormat)
  }
  # 得到下三角矩阵
  get_upper_tri <- function(cormat){
    cormat[lower.tri(cormat)]<- NA
    return(cormat)
  }
  
 # 使用
upper_tri <- get_upper_tri(cormat)
upper_tri

#     mpg  cyl disp   hp drat   wt qsec   vs   am gear carb
#mpg    1 -0.9 -0.9 -0.9  0.7 -0.9  0.5  0.7  0.6  0.5 -0.7
#cyl   NA  1.0  0.9  0.9 -0.7  0.9 -0.6 -0.8 -0.5 -0.6  0.6
#disp  NA   NA  1.0  0.9 -0.7  0.9 -0.5 -0.7 -0.6 -0.6  0.5
#hp    NA   NA   NA  1.0 -0.5  0.8 -0.7 -0.8 -0.4 -0.3  0.7
#drat  NA   NA   NA   NA  1.0 -0.8  0.1  0.4  0.7  0.7 -0.1
#wt    NA   NA   NA   NA   NA  1.0 -0.2 -0.6 -0.7 -0.7  0.5
#qsec  NA   NA   NA   NA   NA   NA  1.0  0.8 -0.2 -0.1 -0.7
#vs    NA   NA   NA   NA   NA   NA   NA  1.0  0.2  0.3 -0.6
#am    NA   NA   NA   NA   NA   NA   NA   NA  1.0  0.8 -0.1
#gear  NA   NA   NA   NA   NA   NA   NA   NA   NA  1.0  0.1
#carb  NA   NA   NA   NA   NA   NA   NA   NA   NA   NA  1.0
  1. 画图
library(reshape2)
melted_cormat <- melt(upper_tri, na.rm = TRUE)
# Heatmap
library(ggplot2)

ggplot(data = melted_cormat, aes(Var2, Var1, fill = value))+
  geom_tile(color = "white")+
  scale_fill_gradient2(low = "green", high = "gold", mid = "white", 
                       midpoint = 0, limit = c(-1,1), space = "Lab", 
                       name="Pearson\nCorrelation") +
  theme_minimal()+ 
  theme(axis.text.x = element_text(angle = 45, vjust = 1, 
                                   size = 12, hjust = 1))+
  coord_fixed() # 这个函数确保x轴上的一个单位长度与y轴上的一个单位长度相同。

图片.png
  1. 根据聚类顺序重新排序绘制相关性图,这样的好处是将相关系数接近的紧挨在一起,更有利于我们观察数据中变化模式。这里我们利用hclust函数方法实现层级聚类。
## 我们写个函数 
reorder_cormat <- function(cormat){
# 使用变量间的相关系数作为距离
dd <- as.dist((1-cormat)/2)
hc <- hclust(dd)
cormat <-cormat[hc$order, hc$order]
}


# 重排序相关矩阵
cormat <- reorder_cormat(cormat)
upper_tri <- get_upper_tri(cormat)
# 变换为长格式数据
melted_cormat <- melt(upper_tri, na.rm = TRUE)
p <- ggplot(data = melted_cormat, aes(Var2, Var1, fill = value))+
  geom_tile(color = "white")+
  scale_fill_gradient2(low = "green", high = "gold", mid = "white", 
                       midpoint = 0, limit = c(-1,1), space = "Lab", 
                       name="Pearson\nCorrelation") +
  theme_minimal()+ 
  theme(axis.text.x = element_text(angle = 45, vjust = 1, 
                                   size = 12, hjust = 1))+
  coord_fixed()
p
图片.png
  1. 将相关性系数增加在图形上
p + 
geom_text(aes(Var2, Var1, label = value), color = "black", size = 4) +
theme(
  axis.title.x = element_blank(),
  axis.title.y = element_blank(),
  panel.grid.major = element_blank(),
  panel.border = element_blank(),
  panel.background = element_blank(),
  axis.ticks = element_blank(),
  legend.justification = c(1, 0),
  legend.position = c(0.6, 0.7),
  legend.direction = "horizontal")+
  guides(fill = guide_colorbar(barwidth = 7, barheight = 1,
                title.position = "top", title.hjust = 0.5))
图片.png

方式二: corrplot

这个包应该算是可视化相关性图最简便的方式了吧,还是推荐大家去看下corrplor包的文档,作者写的算是最详细的。

### R中输入?corrplot看下corrplot中的参数
corrplot(corr, method = c("circle", "square", "ellipse", "number", "shade",
  "color", "pie"), type = c("full", "lower", "upper"), add = FALSE,
  col = NULL, bg = "white", title = "", is.corr = TRUE, diag = TRUE,
  outline = FALSE, mar = c(0, 0, 0, 0), addgrid.col = NULL,
  addCoef.col = NULL, addCoefasPercent = FALSE, order = c("original",
  "AOE", "FPC", "hclust", "alphabet"), hclust.method = c("complete", "ward",
  "ward.D", "ward.D2", "single", "average", "mcquitty", "median", "centroid"),
  addrect = NULL, rect.col = "black", rect.lwd = 2, tl.pos = NULL,
  tl.cex = 1, tl.col = "red", tl.offset = 0.4, tl.srt = 90,
  cl.pos = NULL, cl.lim = NULL, cl.length = NULL, cl.cex = 0.8,
  cl.ratio = 0.15, cl.align.text = "c", cl.offset = 0.5, number.cex = 1,
  number.font = 2, number.digits = NULL, addshade = c("negative",
  "positive", "all"), shade.lwd = 1, shade.col = "white", p.mat = NULL,
  sig.level = 0.05, insig = c("pch", "p-value", "blank", "n", "label_sig"),
  pch = 4, pch.col = "black", pch.cex = 3, plotCI = c("n", "square",
  "circle", "rect"), lowCI.mat = NULL, uppCI.mat = NULL, na.label = "?",
  na.label.col = "black", win.asp = 1, ...)

可以看到包内参数特别之多,但是都很容易理解,我们简单使用一下吧~

method 参数

支持七种作图效果,分别是 "circle" (默认), "square", "ellipse", "number", "pie", "shade" and "color".

type参数控制只输出上三角或者下三角

# 随便用一种
corrplot(cormat,method = 'pie',type = "upper")
图片.png
# 也可以混合搭配,上下三角各用一种效果
corrplot.mixed(cormat, lower = "ellipse", upper = "circle")
图片.png

重排序矩阵进行可视化

内置了四种排序算法:AOE, FPC,hclust,alphabet

# 我们使用最常见的hclust层级聚类方法
corrplot(cormat, order = "hclust",method = 'pie')
图片.png

自定义颜色

addrect = 2 根据层次聚类在图上绘制的矩形的数量,仅在排序方法为hclust时有效。如果为空(默认值),则不添加任何矩形。

col1 <- colorRampPalette(c("#7F0000", "red", "#FF7F00", "yellow", "white",
                           "cyan", "#007FFF", "blue", "#00007F"))
# 这里再推荐几种配色
# col = heat.colors(100)
# col = terrain.colors(100)
# col = cm.colors(100)
# col = brewer.pal(n = 8, name = "PuOr") 利用RColorBrewer包调色
corrplot(cormat, , order = "hclust", addrect = 2, col = col1(100)) # COL1(100)代表渐变色阶划分为100份。
图片.png

计算p值相关矩阵并可视化

利用cor.mtest函数基于原始数据计算相关性,并执行显著性检验。

res1 <- cor.mtest(mtcars, conf.level = .95) # 计算相关显著性P值矩阵
corrplot(cormat, p.mat = res1$p, sig.level = .05) # 设置显著性阈值为0.05
图片.png

除此之外:

设置insig = "p-value" 可将不显著的系数显示出数值

corrplot(cormat, p.mat = res1$p, insig = "p-value")   #  增加sig.level = -1可显示全部系数的P值
图片.png

设置insig = "blank" 可将不显著的设置为空白

corrplot(cormat, p.mat = res1$p, insig = "blank")
图片.png

设置insig = "label_sig" 可将P值以 * 号方式展示出来, 通过设置sig.level = c(.001, .01, .05) 可以将不同范围内的p值用 * 来表示 ,如p值小于0.001的系数设置为三个 * 。

corrplot(cormat, p.mat = res1$p, insig = "label_sig",
         sig.level = c(.001, .01, .05), pch.cex = .9, pch.col = "white")
图片.png

我们还可以自定义显著的系数标签

corrplot(cormat, p.mat = res1$p, insig = "label_sig", pch.col = "grey",
         pch = "显著", pch.cex = .7, order = "AOE")
图片.png

OK, 基本使用如上所述,更多参数自己摸索吧~~

方式三:热图包

热图包有很多种,有空会专门介绍一下R中流行的几个热图包(网上其实也有很多教程了),这里就简单拿俩个常用的举个例子吧~~

pheatmap包

library(RColorBrewer)
coul <- colorRampPalette(brewer.pal(8, "PiYG"))(25)

pheatmap(cormat,cellwidth = NA, cellheight =NA, treeheight_row = 50,
         treeheight_col = 50 ,color =coul,legend=TRUE,
         border_color=NA, fontsize_row=8, fontsize_col=10,
         ,main="Pheatmap")
图片.png

gplot包中的heatmap.2函数

library(gplots)
library(RColorBrewer)
coul <- colorRampPalette(brewer.pal(8, "PiYG"))(25)#换个好看的颜色


heatmap.2(cormat,
          trace="none",#不显示trace
          col=coul,#修改热图颜色
          density.info = "none",#图例取消density
          key.xlab ='Correlation',
          key.title = "",
          cexRow = 1,cexCol = 1,#修改横纵坐标字体
          Rowv = F,Colv = F, #去除聚类
          margins = c(6, 6),
          cellnote = cormat,notecol='black'#添加相关系数的值及修改字体颜色
)

图片.png

方式四:GGally包

“GGally”是“ggplot2”包的一个扩展,它内部添加了几个函数来降低几何对象与转换后的数据组合的复杂性。其中包括相关性图、平行坐标绘图,绘制网络等等的函数,我们这次只介绍下其中的ggcorr函数来绘制相关性图,因为是ggplot2的扩展,所以它支持代码中添加我们熟悉的ggplot语法去设置主题,图例,颜色等等。

参数

R中输入?ggcorr查看参数使用, 也可看该包的Documentation或者Github进行学习。

?ggcorr

ggcorr(data, method = c("pairwise", "pearson"), cor_matrix = NULL,
  nbreaks = NULL, digits = 2, name = "", low = "#3B9AB2",
  mid = "#EEEEEE", high = "#F21A00", midpoint = 0, palette = NULL,
  geom = "tile", min_size = 2, max_size = 6, label = FALSE,
  label_alpha = FALSE, label_color = "black", label_round = 1,
  label_size = 4, limits = c(-1, 1), drop = is.null(limits) ||
  identical(limits, FALSE), layout.exp = 0, legend.position = "right",
  legend.size = 9, ...)  # 参数也不少,下面只会提及几个重要的参数。

简单使用

# Cran
install.packages('GGally')
# Github
library(devtools)
install_github("ggobi/ggally")
library(GGally)

#简单出图
ggcorr(cormat)
图片.png

自定义颜色

通过设置low,mid,high可以设置为三色渐变。

# 图形种类、字体大小、颜色等设置
ggcorr(cormat, ,geom = 'circle',size = 5, 
       low = 'green', mid = 'white', high = 'gold',
       min_size = 5,  max_size = 15,label = TRUE,
       label_alpha = TRUE)  
图片.png

另外,通过nbreaks参数,我们可以将颜色划分为多个范围值。

library(RColorBrewer) # 提供调色板
ggcorr(
  cormat, 
  geom = 'circle',  # 设置图形的展示样式效果,还有"blank","text","tile"三种效果,可以自己试一下
  max_size = 15,    # 控制最大值
  min_size = 5,   # 控制最小值
  size = 3,    
  hjust = 0.75,
  nbreaks = 6,     # 划分为色阶范围。n为整数
  digits = 2,     # 控制小数位数
  palette = 'Set3',  # 如果nbreak参数使用,我们
  label = TRUE,   # 显示标签
  legend.position = "bottom", legend.size = 12)
图片.png

显示值高的系数

参考了这个帖子,了解到我们输入ggcorr(cormat)$data其实可以查看内部相关性数据格式,有时候我们不太关注低相关性的系数,因此我们可以根据coefficient相关系数的那一列数据在相关性图中只显示出绝对值大于多少的系数,并用颜色区分是正相关还是负相关。

head(ggcorr(cormat)$data)

#     x   y coefficient label
#  cyl mpg  -0.9910782  -1.0
# disp mpg  -0.9892051  -1.0
#   hp mpg  -0.9542264  -1.0
# drat mpg   0.9365901   0.9
#   wt mpg  -0.9867463  -1.0

ggcorr(cormat, geom = "blank",
       label = TRUE, hjust = 0.75) +
  geom_point(size = 10, aes(color = coefficient > 0, 
                            alpha = abs(coefficient) > 0.7)) +
  scale_alpha_manual(values = c("TRUE" = 0.5, "FALSE" = 0)) +
  scale_color_manual(name = "Coefficient",labels = c('Negetive','Postive'),values = c("blue","red"))
图片.png

方式五:ggcorrplot

这个包的灵感来自corrplot包, 相比于ggplot2方式该函数该包对不太擅长写函数的同学很友好我们只需要利用包中的参数就可以实现重新排序相关矩阵的方法,而且能在相关图上显示显著性水平。除此之外,该包还包括一个计算相关p值矩阵的函数,总之使用起来还是十分方便的。

安装加载包

# 通过Cran安装
install.packages("ggcorrplot")
# Github安装最新版
if(!require(devtools)) install.packages("devtools")
devtools::install_github("kassambara/ggcorrplot")
# 加载
library(ggcorrplot)

可视化

默认是方形

# method = "square" (默认)
ggcorrplot(cormat)
图片.png

圆圈

ggcorrplot(cormat, method = "circle") 
图片.png

重新排序相关性矩阵

# Reordering the correlation matrix
# --------------------------------
# 使用hclust层级聚类
ggcorrplot(cormat, hc.order = TRUE, outline.col = "white") # 使用order参数
图片.png

只输出上三角或者下三角矩阵的图形

# 得到上三角
ggcorrplot(cormat, hc.order = TRUE, type = "upper",
     outline.col = "white")
# 得到下三角
ggcorrplot(cormat, hc.order = TRUE, type = "lower",
     outline.col = "white")

<img src="https://gitee.com/kai_kai_he/PicGo/raw/master/imgage/image-20201005222129696.png" alt="image-20201005222129696" style="zoom:50%;" />

自定义主题与颜色

ggcorrplot(cormat, hc.order = TRUE, type = "lower",
   outline.col = "white",
   ggtheme = ggplot2::theme_gray,
   colors = c("#6D9EC1", "white", "#E46726"))
图片.png

显示标签

ggcorrplot(cormat, hc.order = TRUE,
    type = "lower", p.mat = p.mat)
图片.png

计算P值显著性矩阵,并在图中显示

p.mat <- cor_pmat(mtcars)
head(p.mat[, 1:4])
##               mpg          cyl         disp           hp
## mpg  0.000000e+00 6.112687e-10 9.380327e-10 1.787835e-07
## cyl  6.112687e-10 0.000000e+00 1.803002e-12 3.477861e-09
## disp 9.380327e-10 1.803002e-12 0.000000e+00 7.142679e-08
## hp   1.787835e-07 3.477861e-09 7.142679e-08 0.000000e+00
## drat 1.776240e-05 8.244636e-06 5.282022e-06 9.988772e-03
## wt   1.293959e-10 1.217567e-07 1.222311e-11 4.145827e-05


ggcorrplot(cormat, hc.order = TRUE,
           type = "lower", p.mat = p.mat)  # 除非没有显著系数

图片.png

除了将不显著的相关系数标记为×号,也可以直接将其留白去除颜色。

# 不显著的系数留白
ggcorrplot(cormat, p.mat = p.mat, hc.order = TRUE,
    type = "lower", insig = "blank")
图片.png
上一篇下一篇

猜你喜欢

热点阅读