一文搞定
R
– 6·10 –
今天,来和大家讨论下统计学中一个常用的概念——相关系数~
1
相关分析概述
具体来说,相关分析的内容主要包括如下内容:
(1)判断确定现象之间有无关系以及相关关系的具体表现形式;
在进行相关分析时,首先通过理论定性的方法或利用相关表和相关图观察的方法,判断现象之间是否有关系,只有现象之间有关系,进行相关分析才有意义。然后判断现象之间相关关系的形态皮尔森相关系数怎么看,以便在之后的分析中选择相应的分析方法。
(2)确定相关关系的密切程度;
根据变量数据的类型,选择适当的方法,计算出相关系数,确定现象之间相关关系的密切程度,为进一步的分析提供依据。
(3)检验相关关系的显著性;
显著性检验包括:检验相关关系的存在性、检验相关关系强度是否达到一定水平、检验现象相关程度的差异性、估计相关系数的取值。
一般来说,在相关分析中,我们使用相关系数(correlation coefficient)来表明变量之间的相关密切程度,通常用r表示,它的取值范围为[-1,1]。r的正、负号可以反映相关的方向:
当r>0时表示线性正相关;
当r
r的大小可以反映相关的程度,r=0表示两个变量之间不存在线性关系。
注意:
① 变量的数据类型不同,相关系数的计算方法不同;
② 变量的数量之间不同,用到的相关分析的方法也不同;
A)简单相关系数:又叫相关系数或线性相关系数,一般用字母r表示,用来度量两个变量间的线性关系;
B)复相关系数:又叫多重相关系数,是指一个因变量与多个自变量之间的相关关系;不同于偏相关系数,它度量多个变量(自变量)的联合效应和某一个变量(因变量)间的相关,反映了变量集和某个变量间的线性关联程度;
C)偏相关系数:在多要素所构成的系统中,先不考虑其他要素的影响,而单独研究两个要素之间的相互关系的密切程度,这称为偏相关。用以度量偏相关程度的统计量,称为偏相关系数。偏相关系数分布的范围在-1到1之间,偏相关系数的绝对值越大,表示其偏相关程度越大。偏相关系数考察变量Xi与Xj的相关关系,同时排除了其他变量的影响。
D)典型相关系数:是先对原来各组变量进行主成分分析,得到新的线性关系的综合指标,再通过综合指标之间的线性相关系数来研究原各组变量间的相关关系。
相关性分析目的与分析方法选择:
① 研究两个随机变量之间的相关性时:
可选用Pearson相关分析或Spearman秩相关分析;
② 研究一个响应变量与k(k≥2)个自变量之间的相关性(k+1个变量服从k+1维正态分布)时:
可选用复相关分析;
③ 现有k个随机变量,以其中一个为因变量,其余为自变量,将其中k-2个自变量分别固定在各自的特定取值上。研究另外一个自变量与因变量之间的相关性时:
可选用偏相关分析;
④ 研究两组定量变量之间的相关性时:
可选用典型相关分析;
在医学研究中,常常会遇到2种情况:
① 多个自变量与一个因变量的相关关系:
此时可以考虑多元相关分析,计算复相关系数和偏相关系数并进行假设检验;
②两组变量之间的线性相关关系:
如研究患者的各种临床症状与所患各种疾病之间的线性相关性,如果是一组自变量(体质指标)与一组因变量(体能指标)之间的关系皮尔森相关系数怎么看,则考虑典型相关分析;
/第一部分·简单相关/
简单相关系数:又叫相关系数或线性相关系数,一般用字母r表示,用来度量两个变量间的线性关系。
1
Pearson相关系数
概述:皮尔逊相关系数( Pearson correlation coefficient)也称为积差相关(或积矩相关)是英国统计学家皮尔逊于20世纪提出的一种计算直线相关的方法;
适用:两个变量都是连续型数据(定距数据);
条件:
当两个变量的标准差都不为零时,相关系数才有定义,皮尔逊相关系数适用于:
① 两个变量之间是线性关系,都是连续数据;
② 两个变量的总体是正态分布,或接近正态的单峰分布;
③ 两个变量的观测值是成对的,每对观测值之间相互独立;
上述任一条件不满足,就用Spearman相关系数,不能用pearson相关系数。
计算公式:两个变量之间的协方差与标准差的商;
R语言代码实现:
# 皮尔森相关性系数
X1<-c(1, 2, 3, 4, 5, 6)
Y1<-c(0.3, 0.9, 2.7, 2, 3.5, 5)
mean(X1) #平均值
mean(Y1)
var(X1) #方差
var(Y1)
sd(X1) #标准差
sd(Y1)
cov(X1,Y1) #协方差
cor(X1,Y1,method="pearson") #皮尔森相关性系数
cor.test(X1,Y1,method="pearson") # 皮尔森相关性系数的检验
注意:cor.test()每次只能检验一种相关关系,一种计算相关矩阵的方式是使用psych包中的corr.test()函数。
corr.test()函数用法如下:
corr.test(x, y = NULL, use = "pairwise",method="pearson",adjust="holm",
alpha=.05,ci=TRUE,minlength=5)
corr.p(r,n,adjust="holm",alpha=.05,minlength=5,ci=TRUE)
# method="pearson","spearman" and "kendall"
一个示例如下:
install.packages("psych") # 安装psych包
library(psych) # 加载包
states <- state.x77[,1:6] # 导入数据
corr.test(states,use="complete")
corr.test(states,use="pairwise",method="pearson") #计算相关系数以及显著性检验,method可为pearson/spearman/kendall,use表示处理缺失的方式
corr.test(states,use="pairwise",method="pearson",adjust="none") #P值不经过多重比较的调整
print(corr.test(states,use="pairwise",method="pearson"),short=F) #显示相关系数和P值的95%CI
Python代码实现:
# 导入库
import pandas as pd
import numpy as np
#原始数据
X1=pd.Series([1, 2, 3, 4, 5, 6])
Y1=pd.Series([0.3, 0.9, 2.7, 2, 3.5, 5])
X1.mean() # 平均值
Y1.mean() #
X1.var() # 方差
Y1.var() #
X1.std() # 标准差不能为0
Y1.std() # 标准差不能为0
X1.cov(Y1) # 协方差
X1.corr(Y1,method="pearson") # 皮尔森相关性系数
X1.cov(Y1)/(X1.std()*Y1.std()) # 皮尔森相关性系数
import scipy.stats as stats
stats.pearsonr(X1,Y1) # 输出结果第一个值为pearsonr相关系数;第二个为p-value;
X1.corr(Y1,method='spearman') # 斯皮尔曼相关系数
X1.corr(Y1,method='kendall') # 肯德尔秩相关系数
2
Spearman秩相关系数
概述:斯皮尔曼秩相关系数(Spearman’s rank correlation coefficient),是一种等级相关系数,是一个非参数性质(与分布无关)的秩统计参数,由 Spearman在 1904年提出,用来度量两个变量之间联系的强弱;
适用:两个变量均为定序变量或两个变量均为不满足正态分布假设的等间隔定距数据;
条件:其利用两变量的秩次大小作线性相关分析,对原始变量的分布不作要求,属于非参数统计方法,适用范围要广些。对于服从 Pearson 相关系数的数据亦可计算 Spearman 相关系数,但统计效能要低一些;
计算公式:斯皮尔曼相关系数与皮尔森相关系数的计算公式是基本相同的,差异在于斯皮尔曼相关系数将皮尔森相关系数中的变量实际数值替换成它们的秩(排名)。表示的是秩次(排名)的线性相关关系。
R代码实现:
X<-c(11,490,14,43,30,3)
Y<-c(2,75,3,44,7,42)
cor(X,Y,method=”spearman”) # 斯皮尔曼相关系数
Python代码实现:
# 导入库
import pandas as pd
import numpy as np
#原始数据
X1=pd.Series([1, 2, 3, 4, 5, 6])
Y1=pd.Series([0.3, 0.9, 2.7, 2, 3.5, 5])
#处理数据删除Nan
x1=X1.dropna()
y1=Y1.dropna()
n=x1.count()
x1.index=np.arange(n)
y1.index=np.arange(n)
#分部计算斯皮尔曼秩相关系数
d=(x1.sort_values().index-y1.sort_values().index)**2
dd=d.to_series().sum()
p=1-n*dd/(n*(n**2-1))
#s.corr()函数计算斯皮尔曼秩相关系数
r=x1.corr(y1,method='spearman')
print(r,p)
import scipy.stats as stats
stats.spearman(X1,Y1) # 输出结果第一个值为spearman相关系数;第二个为p-value;
# 构造秩次表h函数
def show(x1,y1):
print('原始位置x 原x 秩次x 排序x 原始位置y 原y 秩次y 排序y 秩次差的平方')
for i in range(len(x1)):
xx1=x1.sort_values();yy1=y1.sort_values()
ix=x1.index[i]
ixx=xx1.index[i]
iy=y1.index[i]
iyy=yy1.index[i]
d_2=(ixx-iyy)**2
print(' {:5} {:10} {:5} {:5} {:10} {:10} {:5} {:5} {:10.2f}'.format(
ix,x1[i],ixx,xx1[i],iy,y1[i],iyy,yy1[i],d_2) )
show(x1,y1)
3
Kendall相关系数
概述:肯德尔秩相关系数,即Kendall相关系数,与Spearman秩相关系数类似,也是反映分类变量相关性的指标。对相关的有序变量进行非参数相关检验;取值范围在-1-1 之间,此检验适合于正方形表格。
适用:两个变量均为定序变量或两个变量均为不满足正态分布假设的等间隔定距数据;
计算公式:同样作为非参数相关检验,与Spearman秩相关系数不同,Kendall相关系数是从两个定序变量中,对应数值的秩次是否一致着手考察两个变量之间的相关关系,秩次一致的数值越多,那么相关性越强。
Kendall相关系数有三种计算方式,分别为:Kendall’s tau-a; Kendall’s tau-b; Kendall’s tau-c;
其中:
① Kendall’s tau-a:剔除排名相同(秩次相同)的数值对,导致分母偏大,低估两个定序变量之间的相关关系,应用较少;
②Kendall’s tau-b:在tau-a的基础上于分母上减去了秩次相同的数值对,结果更客观准确;
③Kendall’s tau-c:在tau-b的基础上于分母的位置矫正了两个定序变量水平数不同的问题,只考虑水平数少的变量水平;
R语言代码实现:
X<-c(3,1,2,2,1,3)
Y<-c(1,2,3,2,1,1)
cor(X,Y,method='kendall') # 肯德尔秩相关系数
Python代码实现:
# 导入库
import pandas as pd
import numpy as np
#原始数据
height=pd.Series([180,170,160,150,140,130,120,110])
weight=pd.Series([80,75,90,85,70,60,55,65])
# 处理数据
n=weight.count()
df=pd.DataFrame({'a':height,'b':weight})
df1=df.sort_values('b',ascending=False)
df1.index=np.arange(1,df.a.count()+1)
df2=df1.sort_values('a',ascending=False)
s=pd.Series(df2.index)
# 方法一
P=0
for i in s.index:
P=P+s[s > s[i]].count()
s=s.drop(i)
R=(P-(n*(n-1)/2-P))/(n*(n-1)/2)=(4P/(n*(n-1)))-1
R1 =(P - (n * (n - 1) / 2 - P)) / (n * (n - 1) / 2)
print('R1=',R1)
R2= (4*P / (n * (n-1)))-1
print('R2=',R2)
# 方法二
r1=height.corr(weight,method='kendall')
print('r1=',r1)
r2=weight.corr(height,method='kendall')
print('r2=',r2)
import scipy.stats as stats
stats.kendall(X1,Y1) # 输出结果第一个值为kendall相关系数;第二个为p-value;
4
Gamma相关系数
概述:Gamma相关系数也被称为Goodmanand Kruskal’s gamma相关系数。
适用:两个定序型变量的相关性分析
计算公式:
即,如果两个定序型变量的所有数值对都是一致对,那么Gamma的值为+1,反之则为-1;
注意:Gamma相关系数没有考虑秩次相同的数值对;
5
Somers’d相关系数
概述:针对上方Gamma相关系数没有考虑秩次相同数值对的问题,统计学家又设计了Somers’d相关系数进行矫正。
适用:两个定序型变量水平数不同时;
计算公式:Somers’d相关系数有以下两种:
① 公式一的分母加上了y变量的秩次相同的数值对,表示x变量对y变量的相关强度;
② 公式二的分母加上了x变量的秩次相同的数值对,表示y变量对x变量的相关强度;
注意:与其它相关系数不同,Somers’d相关系数有方向性,这是它与其它表示两个定序型变量相关程度的相关系数有所区别的一个重要特点;
6
其他相关系数
相关系数的选择除了考虑两个变量数据类型不同以外,还要考虑到列联表的尺寸以及是否对称。
以上来自东山草堂君的笔记~
/第二部分·复相关和偏相关/
1
复相关系数
在学习一元线性回归分析时,讨论了与之紧密联系的一元相关分析或简单相关分析。将这个概念扩展到多元,就是多元相关分析或复相关分析。简单相关分析研究两个变量之间的关联性,复相关则是研究多个变量之间的关联性。
概述:复相关系数是指在具有多元相关关系的变量中,用来测定因变量y与一组自变量x1, x2,⋅⋅⋅, xm之间相关程度的指标。
注意:
多元回归分析要求将变量分为自变量和因变量,要求因变量服从正态分布;
多元相关分析不要求将变量分为自变量和因变量,要求所有的变量服从正态分布;
计算公式:
复相关系数可以通过计算多因素线性回归的复确定系数,然后开平方根就可以计算出。
R语言代码实现:
在R语言中,可以先通过lm()函数获得复确定系数,然后算出复相关系数。
2
偏相关系数
偏相关,是指在控制一个或多个定量变量时,另外两个定量变量之间的相互关系。
R语言代码实现:
偏相关系数的计算依赖于ggm包中的pcor()函数获得。
pcor()函数调用格式为:
pcor(u,S)
u:是一个数值向量,前两个数值表示要计算相关系数的变量下标,其余的数值为条件变量(即要排除影响的变量)的下标;
S:为变量的协方差阵;
其检验采用ggm包中的pcor.test()函数。
pcor.test()函数调用格式为:
pcor.test(r,q,n)
r:是由pcor()函数计算得到的偏相关系数;
q:要控制的变量数(以数值表示位置);
n:样本大小;
一个示例:
install.packages("ggm") # 安装ggm包
library(ggm)
states <- state.x77[,1:6] # 导入数据
pcor(c(1,5,2,3,4),cov(states)) #偏相关系数计算,向量中的前两位代表计算相关系数的变量位置,后三位代表需要校正的变量位置
pcor.test(0.448,3,30) #偏相关系数显著性检验,0.448是pcor计算的相关系数,3是控制的变量数目,30是样本量
另一个示例:
##偏相关
library(ggm) # 加载包,ggm 包 pcor()
data(mtcars) # 导入数据
#标准化不影响相关系数计算值,但可以让数据服从均值 0,标准差 1 的等方差结构
mtcars <- scale(mtcars)
#要计算相关系数的两个变量,或指定下标
x1 <- c('mpg', 'cyl')
#要控制的条件变量,或指定下标
x2 <- c('drat', 'wt', 'qsec')
#指定协方差矩阵,计算偏相关
pcor_pearson <- pcor(c(x1, x2), cov(mtcars, method = 'pearson'))
pcor_pearson
pcor_spearman <- pcor(c(x1, x2), cov(mtcars, method = 'spearman'))
pcor_spearman
psych包中的r.test()函数提供了多种实用的显著性检验方法,可以用于检验:
① 某种相关系数的显著性;
② 两个独立相关系数的差异是否显著;
③ 两个基于一个共享变量得到的非独立相关系数的差异是否显著;
④ 两个基于完全不同的变量得到的非独立相关系数的差异是否显著;
一个示例:
r.test(30,0.5,0.7) #两个独立相关系数的差异是否显著,30为样本量,0.5和0.7为相关系数,如果两样本量不一样,则可在括号中增加n2
r.test(30,r12=0.5,r13=0.7,r23=0.2) #两个基于一个共享变量得到的非独立相关系数的差异是否显著
r.test(30,r12=0.5,r34=0.7,r13=0.5,r14=0.6,r23=0.7,r24=0.8) #两个基于完全不同的变量得到的非独立相关系数的差异是否显著
/第三部分·相关性的绘图/
相关系数矩阵,常用相关性系数图来进行展示。相关系数图一般都是三维及以上的数据,但是使用二维图表显示。
其中:
X列、Y列:为都为类别数据,分别对应图表的X轴和Y轴;
Z列:为数值信息,通过颜色饱和度、面积大小等视觉暗示表示;
图1-6-21(a)使用颜色饱和度和颜色色相综合表示Z列数据;
图1-6-21(b)使用方块的面积大小及颜色综合表示Z列数据;
从图中很容易观察到哪两组变量的相关性最好;
R语言中实现相关性热图的方法有很多种:例如heatmap、pheatmap、corrplot、ggcorrplot等,接下来为大家一一演示。
① heatmap()函数:
# 导入数据
data(mtcars)
#标准化不影响相关系数计算值,但可以让数据服从均值 0,标准差 1 的等方差结构
mtcars <- scale(mtcars)
#协方差计算,cov()
cov_pearson <- cov(mtcars, method = 'pearson')
cov_pearson
cov_spearman <- cov(mtcars, method = 'spearman')
cov_spearman
cov_kendall <- cov(mtcars, method = 'kendall')
cov_kendall
#相关系数计算,cor()
cor_pearson <- cor(mtcars, method = 'pearson')
cor_pearson
cor_spearman <- cor(mtcars, method = 'spearman')
cor_spearman
cor_kendall <- cor(mtcars, method = 'kendall')
cor_kendall
# ① heatmap相关图
library(stats)
heatmap(cor_pearson,symm = TRUE)
② pheatmap()函数:
# ② pheatmap相关图
library(pheatmap)
pheatmap(cor_pearson,display_numbers = TRUE,number_color = "#C0C0C0",cluster_row = FALSE, cluster_col = FALSE,color = colorRampPalette(c("navy", "white", "firebrick3"))(50))
③ggcorrplot相关矩阵图:
# ③ ggcorrplot相关矩阵图
library(ggcorrplot)
library(cowplot)
p1 <- ggcorrplot(cor_pearson,method='square')
p2 <- ggcorrplot(cor_pearson,method="circle")
plot_grid(p1,p2,ncol=2,labels=c('A','B'),align=c('v','h'))
p1 <- ggcorrplot(cor_pearson,type='lower',outline.color='white',colors=c('#6D9EC1','white','#E46726'))
p2 <- ggcorrplot(cor_pearson,type='upper',outline.color='white',colors=c('#084594','white','#ef3b2c'))
p3 <- ggcorrplot(cor_pearson,type='lower',lab=TRUE,lab_size = 3)
p.mat <- cor_pmat(mtcars) # 显示P值,是否有统计学标示
p4 <- ggcorrplot(cor_pearson,type='lower',p.mat=p.mat)
plot_grid(p1,p2,p3,p4,ncol=2,labels=LETTERS[1:4],align=c('v','h'))
④corrplot相关图:
# corrplot相关图
# 换用不同method,分别显示数字和颜色(method can be "circle", "square", "ellipse", "number", "shade", "color", "pie")
# method 默认为“circle”
library(corrplot)
corrplot(cor_pearson, method = 'number', number.cex = 0.8, diag = FALSE, tl.cex = 0.8)corrplot(cor_pearson, add = TRUE, type = 'upper', method = 'pie', diag = FALSE, tl.pos = 'n', cl.pos = 'n')
#输出,例如
write.table(cor_pearson, 'cor_pearson.txt', sep = 't', col.names = NA, quote = FALSE)
corrplot(cor_pearson, add = TRUE, type = 'upper', method = 'pie', diag = FALSE, tl.pos = 'n', cl.pos = 'n')
#输出,例如
write.table(cor_pearson, 'cor_pearson.txt', sep = 't', col.names = NA, quote = FALSE)
# 当order = "hclust"时,可使用hclust.method选择层次聚类的方法
# hclust.method可以为“complete”, “ward”, “single”, “average”, “mcquitty”, “median”, “centroid”
corrplot(cor_pearson,order = "hclust",addrect = 3,hclust.method="ward.D2")
#控制对角标签旋转45度
corrplot(cor_pearson, type = "lower", order = "hclust", tl.col = "black", tl.srt = 45)
其他的绘图方法,还有ggcor包,也能实现十分惊艳的效果~
———END———
限 时 特 惠: 本站每日持续更新海量各大内部创业教程,一年会员只需98元,全站资源免费下载 点击网站首页每天更新
站 长 微 信: aiwo51889