R的典型相关分析感觉比较复杂,还是SAS简单
数据来自SAS:
Weight Waist Pulse Chins Situps Jumps
191 36 50 5 162 60
189 37 52 2 110 60
193 38 58 12 101 101
162 35 62 12 105 37
189 35 46 13 155 58
182 36 56 4 101 42
211 38 56 8 101 38
167 34 60 6 125 40
176 31 74 15 200 40
154 33 56 17 251 250
169 34 50 17 120 38
166 33 52 13 210 115
154 34 64 14 215 105
247 46 50 1 50 50
193 36 46 6 70 31
202 37 62 12 210 120
176 37 54 4 60 25
157 32 52 11 230 80
156 33 54 15 225 73
138 33 68 2 110 43
WilksL F df1 df2 p1
0.3503905 2.04823353 9 34.22293 0.063530942
0.9547227 0.17578229 4 30.00000 0.949120253
0.9947336 0.08470926 1 16.00000 0.77475327
数据来自SAS:
Weight Waist Pulse Chins Situps Jumps
191 36 50 5 162 60
189 37 52 2 110 60
193 38 58 12 101 101
162 35 62 12 105 37
189 35 46 13 155 58
182 36 56 4 101 42
211 38 56 8 101 38
167 34 60 6 125 40
176 31 74 15 200 40
154 33 56 17 251 250
169 34 50 17 120 38
166 33 52 13 210 115
154 34 64 14 215 105
247 46 50 1 50 50
193 36 46 6 70 31
202 37 62 12 210 120
176 37 54 4 60 25
157 32 52 11 230 80
156 33 54 15 225 73
138 33 68 2 110 43
<br />
fit <- read.table("D:/data/fit.txt", header=TRUE, sep="\t", na.strings="0.0", dec=".", strip.white=TRUE)<br />
attach(fit)<br />
#空间统计包<br />
library(fields)<br />
# t() 转置函数<br />
#stats() fields包中的函数,显示数据的一些描述性统计量<br />
t(stats(fit))<br />
<br />
# 将变量分成两类<br />
wwp<-fit[,1:3]<br />
csj<-fit[,4:6]<br />
#cancor()在stats包中,典型相关分析<br />
cancor(wwp,csj)<br />
#CCA 典型相关分析包<br />
library(CCA)<br />
# 计算两类数据的相关矩阵 ,输出结果也包括但数据内变量相关系数<br />
matcor(wwp,csj)<br />
# cc 增强典型相关分析<br />
cc1 <- cc(wwp,csj)<br />
# 显示相关性<br />
cc1[1]<br />
# raw canonical coefficients <br />
cc1[3:4]<br />
# 进一步分析,comput为cc的补充<br />
cc2<-comput(wwp,csj, cc1)<br />
#显示典型变量及其相关系数<br />
cc2[3:6]<br />
<br />
# 标准维数检验<br />
<br />
ev<-cc1$cor^2<br />
ev2<-1-ev<br />
n<-dim(wwp)[1]<br />
p<-length(wwp)<br />
q<-length(csj)<br />
m<-n -3/2 - (p+q)/2<br />
w<-cbind(NULL) # initialize wilks lambda<br />
<br />
for (i in 1:3){<br />
w<-cbind(w,prod(ev2[i:3]))<br />
}<br />
<br />
d1<-cbind(NULL)<br />
d2<-cbind(NULL)<br />
f<-cbind(NULL) # initialize f<br />
for (i in 1:3){<br />
s<-sqrt((p^2*q^2-4)/(p^2+q^2-5))<br />
si<-1/s<br />
df1<-p*q<br />
d1<-cbind(d1,df1)<br />
df2<-m*s-p*q/2+1<br />
d2<-cbind(d2,df2)<br />
r<-(1-w[i]^si)/w[i]^si<br />
f<-cbind(f,r*df2/df1)<br />
p<-p-1<br />
q<-q-1<br />
}<br />
pv<-pf(f,d1,d2,lower.tail=FALSE) #F分布<br />
dmat<-cbind(t(w),t(f),t(d1),t(d2),t(pv))<br />
colnames(dmat)<-c("WilksL","F","df1","df2","p")<br />
rownames(dmat)<-c(seq(1:length(w)))<br />
dmat<br />
# wwp标准化典型相关系数<br />
<br />
sd<-sd(wwp)<br />
s1<-diag(sd) # diagonal matrix of wwp sd's<br />
s1 %*% cc1$xcoef # %*% 矩阵乘法<br />
sd<-sd(csj)<br />
s2<-diag(sd) # diagonal matrix of acad sd's<br />
s2 %*% cc1$ycoef<br />
<br />
WilksL F df1 df2 p1
0.3503905 2.04823353 9 34.22293 0.063530942
0.9547227 0.17578229 4 30.00000 0.949120253
0.9947336 0.08470926 1 16.00000 0.77475327