相关性分析及可视化R笔记
2019-11-08 本文已影响0人
想成为你_50a2
近日，因宿舍小赵同学看文章，看到一张关于相关性分析结果可视化的图表，并对其表达含义表示了百思不得其解，于是在我们的男神428宿舍引发了一波生信小热潮。可能对于大佬们来说这是很简单的图了，但是大家还是对其表达内容及绘制方法很感兴趣，因此通过发邮件跟文章作者咨询（老师态度很好，讲的很详细，因为未征求同意，此处就不揭露文章及本人名字了[调皮+1]），详细了解了一波，一语惊醒梦中人。但是梦并不止步，搜集搜集再搜集，通过度娘和简书上相关类似的做法，进行了高层次的对比筛选，找到两种个人感觉最适合展示相关性分析结果的可视化表达方式，并补加了一些友好的注释，整理如下：
相关性分析：是指对两个或多个具备相关性的变量元素进行分析，从而衡量两个变量因素的相关密切程度。相关性的元素之间需要存在一定的联系或者概率才可以进行相关性分析。
如图，就是这件事情的“罪魁祸首”
相关性分析及频率分布图
首先讲解一下其含义：
全图围绕柱状图所在的对角线呈两两对应，两两相关的局势，右上角部分表示显著情况，NS为不显著，即p值>0.05；※为显著，即p值 ≤0.05；**为非常显著，即p值≤0.01。左下角部分表示两两指标间相关性的散点分布图，即2行1图对应GP与PH之间的相关性散点分布，1行2图对应其显著性。对角线处的柱状图对应该指标的频率分布即正相关分布图。
然后，直接上干货
1.相关性分析及频率分布可视化
1.1相关系数及p值的计算
#mtcars为R自带的数据集
mtcars <- data("mtcars")
head(mtcars)
#> head(mtcars)
mpg cyl disp hp drat wt qsec vs am gear carb
Mazda RX4 21.0 6 160 110 3.90 2.620 16.46 0 1 4 4
Mazda RX4 Wag 21.0 6 160 110 3.90 2.875 17.02 0 1 4 4
Datsun 710 22.8 4 108 93 3.85 2.320 18.61 1 1 4 1
Hornet 4 Drive 21.4 6 258 110 3.08 3.215 19.44 1 0 3 1
Hornet Sportabout 18.7 8 360 175 3.15 3.440 17.02 0 0 3 2
Valiant 18.1 6 225 105 2.76 3.460 20.22 1 0 3 1
install.packages("Hmisc")
library(Hmisc)
res <- rcorr(as.matrix(mtcars))
res
#> res
mpg cyl disp hp drat wt qsec vs am gear carb
mpg 1.00 -0.85 -0.85 -0.78 0.68 -0.87 0.42 0.66 0.60 0.48 -0.55
cyl -0.85 1.00 0.90 0.83 -0.70 0.78 -0.59 -0.81 -0.52 -0.49 0.53
disp -0.85 0.90 1.00 0.79 -0.71 0.89 -0.43 -0.71 -0.59 -0.56 0.39
hp -0.78 0.83 0.79 1.00 -0.45 0.66 -0.71 -0.72 -0.24 -0.13 0.75
drat 0.68 -0.70 -0.71 -0.45 1.00 -0.71 0.09 0.44 0.71 0.70 -0.09
wt -0.87 0.78 0.89 0.66 -0.71 1.00 -0.17 -0.55 -0.69 -0.58 0.43
qsec 0.42 -0.59 -0.43 -0.71 0.09 -0.17 1.00 0.74 -0.23 -0.21 -0.66
vs 0.66 -0.81 -0.71 -0.72 0.44 -0.55 0.74 1.00 0.17 0.21 -0.57
am 0.60 -0.52 -0.59 -0.24 0.71 -0.69 -0.23 0.17 1.00 0.79 0.06
gear 0.48 -0.49 -0.56 -0.13 0.70 -0.58 -0.21 0.21 0.79 1.00 0.27
carb -0.55 0.53 0.39 0.75 -0.09 0.43 -0.66 -0.57 0.06 0.27 1.00
n= 32
P
mpg cyl disp hp drat wt qsec vs am gear carb
mpg 0.0000 0.0000 0.0000 0.0000 0.0000 0.0171 0.0000 0.0003 0.0054 0.0011
cyl 0.0000 0.0000 0.0000 0.0000 0.0000 0.0004 0.0000 0.0022 0.0042 0.0019
disp 0.0000 0.0000 0.0000 0.0000 0.0000 0.0131 0.0000 0.0004 0.0010 0.0253
hp 0.0000 0.0000 0.0000 0.0100 0.0000 0.0000 0.0000 0.1798 0.4930 0.0000
drat 0.0000 0.0000 0.0000 0.0100 0.0000 0.6196 0.0117 0.0000 0.0000 0.6212
wt 0.0000 0.0000 0.0000 0.0000 0.0000 0.3389 0.0010 0.0000 0.0005 0.0146
qsec 0.0171 0.0004 0.0131 0.0000 0.6196 0.3389 0.0000 0.2057 0.2425 0.0000
vs 0.0000 0.0000 0.0000 0.0000 0.0117 0.0010 0.0000 0.3570 0.2579 0.0007
am 0.0003 0.0022 0.0004 0.1798 0.0000 0.0000 0.2057 0.3570 0.0000 0.7545
gear 0.0054 0.0042 0.0010 0.4930 0.0000 0.0005 0.2425 0.2579 0.0000 0.1290
carb 0.0011 0.0019 0.0253 0.0000 0.6212 0.0146 0.0000 0.0007 0.7545 0.1290
#输出说明
#r ：第一个矩阵为相关性矩阵
#n : 处理数据的总记录数（行数）
#P : 显著性水平矩阵（越小说明越显著）
1.2 返还相关性矩阵r
signif(res$r,2)
#> signif(res$r,2)
mpg cyl disp hp drat wt qsec vs am gear carb
mpg 1.00 -0.85 -0.85 -0.78 0.680 -0.87 0.420 0.66 0.600 0.48 -0.550
cyl -0.85 1.00 0.90 0.83 -0.700 0.78 -0.590 -0.81 -0.520 -0.49 0.530
disp -0.85 0.90 1.00 0.79 -0.710 0.89 -0.430 -0.71 -0.590 -0.56 0.390
hp -0.78 0.83 0.79 1.00 -0.450 0.66 -0.710 -0.72 -0.240 -0.13 0.750
drat 0.68 -0.70 -0.71 -0.45 1.000 -0.71 0.091 0.44 0.710 0.70 -0.091
wt -0.87 0.78 0.89 0.66 -0.710 1.00 -0.170 -0.55 -0.690 -0.58 0.430
qsec 0.42 -0.59 -0.43 -0.71 0.091 -0.17 1.000 0.74 -0.230 -0.21 -0.660
vs 0.66 -0.81 -0.71 -0.72 0.440 -0.55 0.740 1.00 0.170 0.21 -0.570
am 0.60 -0.52 -0.59 -0.24 0.710 -0.69 -0.230 0.17 1.000 0.79 0.058
gear 0.48 -0.49 -0.56 -0.13 0.700 -0.58 -0.210 0.21 0.790 1.00 0.270
carb -0.55 0.53 0.39 0.75 -0.091 0.43 -0.660 -0.57 0.058 0.27 1.000
1.3 返还P值矩阵
signif(res$P,2)
#> signif(res$P,2)
mpg cyl disp hp drat wt qsec vs am gear carb
mpg NA 6.1e-10 9.4e-10 1.8e-07 1.8e-05 1.3e-10 1.7e-02 3.4e-05 2.9e-04 5.4e-03 1.1e-03
cyl 6.1e-10 NA 1.8e-12 3.5e-09 8.2e-06 1.2e-07 3.7e-04 1.8e-08 2.2e-03 4.2e-03 1.9e-03
disp 9.4e-10 1.8e-12 NA 7.1e-08 5.3e-06 1.2e-11 1.3e-02 5.2e-06 3.7e-04 9.6e-04 2.5e-02
hp 1.8e-07 3.5e-09 7.1e-08 NA 1.0e-02 4.1e-05 5.8e-06 2.9e-06 1.8e-01 4.9e-01 7.8e-07
drat 1.8e-05 8.2e-06 5.3e-06 1.0e-02 NA 4.8e-06 6.2e-01 1.2e-02 4.7e-06 8.4e-06 6.2e-01
wt 1.3e-10 1.2e-07 1.2e-11 4.1e-05 4.8e-06 NA 3.4e-01 9.8e-04 1.1e-05 4.6e-04 1.5e-02
qsec 1.7e-02 3.7e-04 1.3e-02 5.8e-06 6.2e-01 3.4e-01 NA 1.0e-06 2.1e-01 2.4e-01 4.5e-05
vs 3.4e-05 1.8e-08 5.2e-06 2.9e-06 1.2e-02 9.8e-04 1.0e-06 NA 3.6e-01 2.6e-01 6.7e-04
am 2.9e-04 2.2e-03 3.7e-04 1.8e-01 4.7e-06 1.1e-05 2.1e-01 3.6e-01 NA 5.8e-08 7.5e-01
gear 5.4e-03 4.2e-03 9.6e-04 4.9e-01 8.4e-06 4.6e-04 2.4e-01 2.6e-01 5.8e-08 NA 1.3e-01
carb 1.1e-03 1.9e-03 2.5e-02 7.8e-07 6.2e-01 1.5e-02 4.5e-05 6.7e-04 7.5e-01 1.3e-01 NA
1.4 可视化
#install.packages("PerformanceAnalytics")
library(PerformanceAnalytics)
chart.Correlation(mtcars, histogram=TRUE, pch=19)
清晰明了的相关性分析及频率分布图
正常的操作到这里按理说就可以停止了，但是很显然，并不适合精神小伙们的气质，继续阅读，继续总结，，
2 可导出的提取相关性矩阵r及其P值的方法
Function 函数旨在创建一个新的 Function 对象
关于一些矩阵中数据提取的函数，查询了度娘，总结：
diag()可以取出矩阵对角线数据, 修改对角线数据, 生成指定对角线数值的矩阵等；lower.tri 可以取出矩阵对角线下方的数值；upper.tri 可以取出矩阵对角线上方的数值。实际上这两个函数是将对角线两边的位置置为TRUE或FALSE.
2.1 构建一个以"row,column,cor,p"为列名的矩阵提取函数
corMatrix <- function(cor,p){
ut <- upper.tri(cor)
data.frame(row = rownames(cor)[row(cor)[ut]],
column = rownames(cor)[col(cor)[ut]],
cor = (cor)[ut],
p = p[ut])
}
2.2 相关系数及p值的矩阵构建
res <- rcorr(as.matrix(mtcars))
ppp <- corMatrix(res$r,res$P)
head(ppp)
> head(ppp)
row column cor p
1 mpg cyl -0.8521620 6.112688e-10
2 mpg disp -0.8475514 9.380328e-10
3 cyl disp 0.9020329 1.803002e-12
4 mpg hp -0.7761684 1.787835e-07
5 cyl hp 0.8324475 3.477861e-09
6 disp hp 0.7909486 7.142679e-08
2.3 将矩阵导出csv
write.table(ppp,file = "ppp.csv",
row.names = F,sep = ",")
3 Heatmap可视化相关性分析
3.1 计算矩阵相关系数
#cor来自R自带的统计函数stats，计算相关系数
matcar.cor <- cor(mtcars)
matcar.cor
3.2 heatmap可视化
install.packages("corrplot")
library(corrplot)
#round函数返回按照指定的小数位数进行四舍五入运算的结果,2代表两位小数
round(matcar.cor,2)
#相关性分析可视化，order对数值做聚类；
#method控制数值的可视化图形种类，ellipse等；addCoef.col添加对应颜色的数值在图形上
corrplot(matcar.cor,order = "AOE",
method = "color",
addCoef.col = "gray")
清晰明了的Corr Heamap
至此，本次热潮暂歇，但是“幸运”的是，无意间我又看到了一种更加简洁的热图表达方式，但是一直没有时间去感悟，正好可怜的“小白鼠”室友有了一点数据（坏笑+10），，，
最近太颓废了，，又是一个“玩”到两点的夜晚，溜了溜了，，继续挖掘啦，期待下一次的笔记！生信真的太有意思了，，，，
