利用R求简单相关系数及其显著性检验

1、题目：测定“丰产三号”小麦的每株穗数X1,主茎上每穗结实小穗数x2,百粒重 x3 （单位：g），主茎株高x4 （单位：cm）和每株籽粒产量y （单位：g）的15组观测值如下表，试计算其简单相关系数并作相关系数的显著性检验。ixi1xi2xi3xi4y110233.611315.729203.610614.5310223.711117.5413213.710922.5510223.611015.5610233.510316.978233.31008.6810243.411417910203.410413.71010213.411013.41110233.910420.3128213.510910

2、.2136233.21147.4148213.711311.6159223.610512.3丰产3号 小麦栽培试验的观测值利用R软件自带函数计算简单相关系数及其检验p值:利用R软件自带的计算简单相关系数的函数cor（）得到简单相关系数表:x1x2x3x4yx11.0000-0.13570.5007-0.09390.8973x2-0.13571.0000-0.14890.12340.0461x30.0057-0.14891.0000-0.03580.6890 x4-0.09390.1234-0.03581.0000-0.0065y0.89730.04620.6890-0.00651.0000再调

3、用cor.test（涵数得到相关系数检验的p值表:x2x3x4yx10.6294-0.057280.73925.75e-06*x20.59640.66130.8702x30.89910.004499*x40.9816手工编写相关系数计算函数和检验函数:COR.test = function(X,R);#求F检验的p值，为矩阵形式COR = function(X);#求相关系数，为矩阵形式lxy = function(x,y);#各个观测值之间的离均差乘积得到的简单相关系数表:x1x2x3x4yX11.0000-0.13570.5007-0.09390.8973x2-0.13571.0000-0

4、.14890.12340.0461x30.0057-0.14891.0000-0.03580.6890 x4-0.09390.1234-0.03581.0000-0.0065y0.89730.04620.6890-0.00651.0000得到的检验p值表:x1x2x3x4yx10.000000.629600.057280.739205.75E-06x20.629600.000000.596410.661300.87020 x30.057280.596410.000000.899100.00450 x40.739200.661300.899130.000000.98160y0.000010.87

5、0200.004500.981600.00000与R软件自带的函数得到的结果相比较，可以看出相关系数结果一致，检验p 值偏差可以忽略。各个函数具体代码如下：COR.test = function(X,R)options(digits = 4) #求 F 检验的 p 值，为矩阵形式n = dim(X)2;#得到p矩阵的阶数p = diag(0,n) ; #n 阶零矩阵for (i in 1:n)for (j in 1:n)f = Ri)jA2/(1-RiJjA2)/(dim(X)1-2); #用F检验对相关系数作显著性检验 pi,j=1-pf(f,1,dim(X)1-2); #用 F 检验计算

6、p 值pCOR = function(X)options(digits = 3)# 求相关系数，为矩阵形式n = dim(X)2;#得到相关系数矩阵的阶数Y = diag(n);#产生单位矩阵for (i in 1:n)for( j in 1:n)Yi,j = lxy(X,i,X,j)/sqrt(lxy(X,i,X,i)*lxy(X,j,X,j);# 简单相关系数计算公式Ylxy = function(x,y)#各个观测值之间的离均差乘积n = length(x)sum(x*y) - sum(x)*sum(y)/nf_COR = function(R)#复相关系数，R为简单相关系数矩阵R =

7、solve(R)#简单相关系数矩阵的逆矩阵n = dim(R)1r = numeric(n)#产生 n 阶零向量for (i in 1:n)ri = sqrt(1 - 1/Ri,i)复相关系数的计算公式rp_COR = function(X,R)#偏相关系数，X为输入数据框或者矩阵，R为简单相关系数矩阵R = solve(R)#简单相关系数矩阵的逆矩阵n = dim(R)1n1 = dim(X)1p = dim(X凡2 - 1#自由度为 n1-p-1r = diag(0,n) #产生n阶零矩阵for (i in 1:(n - 1)#对角线以上部分for (j in (i + 1):n)ri,j = -(1)*Ri,j/sqrt(Ri,i*Rj,j) # 偏相关系数计算公式for (i in 1:(n - 1)for (j in (i + 1):n)#对角线以下部分F = ri,j人2/(1