文章期号:20190604
第二章统计进阶,多元统计:典型相关分析
电信相关分析(canonical correlation analysis)利用综合变量对之间的相关关系来反映两组指标之间的整体相关性的多元统计分析方法。
它的基本原理是:为了从总体上把握两组指标之间的相关关系,分别在两组变量中提取有代表性的两个综合变量U1和V1(分别为两个变量组中各变量的线性组合),利用这两个综合变量之间的相关关系来反映两组指标之间的整体相关性
案例分析:我国科研资金投入与科技成果产出
1 ,计算相关系数矩阵:
R相关矩阵拆成四个子块
> X<-read.table(pipe("pbpaste"),header=T)
> R<-cor(X);R
R.D人员全时当量 R.D经费支出 政府资金 企业资金 发表科技论文 专利申请受理数 发明专利
R.D人员全时当量 1.0000000 0.9908535 0.9911356 0.9504280 0.9844693 0.9840037 0.9363118
R.D经费支出 0.9908535 1.0000000 0.9988826 0.9584249 0.9692064 0.9931741 0.9510069
政府资金 0.9911356 0.9988826 1.0000000 0.9503002 0.9739008 0.9869929 0.9363016
企业资金 0.9504280 0.9584249 0.9503002 1.0000000 0.9037032 0.9679390 0.9547370
发表科技论文 0.9844693 0.9692064 0.9739008 0.9037032 1.0000000 0.9516937 0.8849569
专利申请受理数 0.9840037 0.9931741 0.9869929 0.9679390 0.9516937 1.0000000 0.9796435
发明专利 0.9363118 0.9510069 0.9363016 0.9547370 0.8849569 0.9796435 1.0000000
> 其中:R11=R[1:4, 1:4], R12=[1:4,5:7], R21=R[5:7,1:4], R22=R[5:7, 5:7] 注释:R是7*7的矩阵
2,计算特征值和典型相关系数:
> A = solve(R[1:4, 1:4])%*%R[1:4,5:7]%*%solve(R[5:7, 5:7])%*%R[5:7,1:4]
> A
R.D人员全时当量 R.D经费支出 政府资金 企业资金
R.D人员全时当量 0.6179116 0.148934508 0.16561570 -0.09490055
R.D经费支出 0.6752529 1.125554468 0.28245759 3.70189199
政府资金 -0.1990924 -0.267044777 0.56207651 -2.68974200
企业资金 -0.1012152 -0.008291917 -0.01149849 0.03873382
> ev<-eigen(A)$values; ev # 特征值
[1] 9.993448e-01 8.555178e-01 4.894137e-01 7.438494e-15
> sprt(ev)
Error in sprt(ev) : 没有"sprt"这个函数
> sqrt(ev)
[1] 9.996724e-01 9.249421e-01 6.995811e-01 8.624671e-083,典型相关系数检验,确定典型变量
> xy<-scale(X)
> ca<-cancor(xy[,1:4],xy[,5:7])
> ca$cor
[1] 0.9996724 0.9249421 0.6995811
> ca$xcoef #x 的典型载荷
[,1] [,2] [,3] [,4]
R.D人员全时当量 -0.08749404 0.33149835 2.4874802 0.5207858
R.D经费支出 -0.59662755 -6.82768318 2.0161033 -5.5220954
政府资金 0.34091961 6.57868886 -4.0255389 3.7658502
企业资金 0.01028809 -0.06437817 -0.4930343 1.2966516
> ca$ycoef #y 的典型载荷
[,1] [,2] [,3]
发表科技论文 -0.02225256 0.2496113 1.677748
专利申请受理数 -0.41077257 1.0529651 -3.769277
发明专利 0.10170566 -1.3439072 2.217174
>
> ca$xcenter
R.D人员全时当量 R.D经费支出 政府资金 企业资金
2.005340e-16 5.689893e-17 1.033895e-16 -1.106754e-16
> ca$ycenter
发表科技论文 专利申请受理数 发明专利
1.568190e-16 -1.387779e-17 -9.228729e-17
>
> u<-as.matrix(xy[,1:4])%*%ca$xcoef
> v<-as.matrix(xy[,5:7])%*%ca$ycoef
>
> corcoef.test<-function(r,n,p,q,alpha=0.5){
+ m<-length(r);Q<-rep(0,m);lambda<-1
+ for(k in m:1){
+ lambda<-lambda*(1-r[k]^2);
+ Q[k]<-log(lambda)
+ }
+ s<-0;i<-m
+
+ for(k in 1:m){
+ Q[k]<-(n-k+1-1/2*(p+q+3)+s)*Q[k]
+ chi<-1-pchisq(Q[k],(p-k+1)*(q-k+1))
+ if(chi>alpha){
+ i<-k-1;break;
+ }
+ s<-s+1/r[k]^2
+ }
+ i
+ }
> corcoef.test(r=ca$cor,n=10,p=4,q=3)
[1] 0从R计算结果可以看出x,y 前两对[,1] , [,2]典型载荷的典型相关变量可以用表达表示:
第一组:
u1 = -0.08749404*R.D人员全时当量 - 0.59662755*R.D经费支出 + 0.34091961 *政府资金 +
0.01028809 * 企业资金;
v1 = -0.02225256 * 发表科技论文 - 0.41077257 *专利申请受理数 + 0.10170566* 发明专利
大家可以写出第二组典型相关系数表达式。
从第一组可以看出 R.D经费支出 ,政府资金 相关其它变量有较大的载荷,说明经费和政府资金对科研产出有很大的影响力,占据科研投入的主导地位。 同时,专利申请受理数,发明专利 是科研产出的主力军。
前两组典型相关系数效果图:(第一组两个典型相关系数,存在强的线性相关性,也体现来科研投入与科研产出的线性相关性)
> par(mfrow=c(1,2))
> plot(u[,1],v[,1],xlab='u1',ylab='v1');
> abline(0,1)
> plot(u[,2],v[,2],xlab='u2',ylab='v2');
> abline(0,1)前两组典型相关系数效果图