《多元统计分析-基于R（第3版）》课件 第9-11章-主成分分析、因子分析、对应分析

第9章主成分分析2024/7/312

主成分分析（principalcomponentanalysis，PCA）是将具有相关关系的多个变量有效地转化为少数几个综合变量来处理，从而简化相关统计分析的一种多元统计方法．主成分分析也称主分量分析，是由Pearson于1901年首先提出，到1933年由Hotelling加以完善后发展起来的．目前，在涉及到高维数据分析处理的诸多领域主成分分析都有广泛的应用．本章介绍主成分的基本思想、总体主成分和样本主成分的定义、性质、主成分的计算和解释、主成分回归等．2024/7/313

基本思想：用个数较少，但是保留了原始变量的大部分信息几个不相关的综合变量(即主成分)来代替原来较多的变量，从而可以简化数据，对原来复杂的数据关系进行简明有效的统计分析．

主成分分析的本质就是“有效降维”，既要减少变量个数，又不能损失太多信息．换句话说，就是“降噪”、“挤水分”或者说“冗余消除”，将高维数据有效地转化为低维数据来处理．揭示变量之间的内在联系，进而分析解决实际问题.

第9章

主要内容9.1

主成分分析的基本思想9.2

总体主成分

主成分的含义、计算、主要性质、主成分个数的确定、变量的标准化．9.3

样本主成分

样本主成分的性质和计算；主成分分析的步骤和相关R函数．9.4案例：主成分回归分析2024/7/314

9.1

主成分分析的基本思想2024/7/315

实际统计分析中，经常处理多变量、大维数数据分析问题，分析过程较复杂，难度较大.全部变量中可能存在信息的重叠．为去除这些信息重叠，希望用个数较少，但是保留了原始变量的大部分信息几个不相关的综合变量(即主成分)来代替原来较多的变量，达到降维的目的，从而进行简明有效的统计分析.主成分分析中的信息，就是指变量的变异性，常用标准差或方差来表示它.9.2总体主成分

9.2.1

主成分的含义

以二维正态分布样本点来直观说明：如图，设有n个样品点大致分布在平面上一个椭圆内：2024/7/316

图9－1二维情形主成分示意图二维情形主成分的解释

样本点之间的差异是由

x1和

x2的变化引起的，两者变动的相差不大，但如果用新坐标

y1和

y2来代替，易见，这些样本点的差异主要体现在

y1轴上，n个点在

y1轴方向上的方差达到最大，即在此方向上包含了有关n个样品的最多的信息．

将这些点投影到

y1轴方向能使信息的损失最小，如果

y1轴方向的差异占了全部样本点差异的绝大部分，那么将

y2忽略是合理的，这样就把两个变量简化为一个，显然这里的

y1轴代表了数据变化最大的方向，称之为第一主成分．y2称为第二主成分，并要求已经包含在

y1中的信息不出现在

y2中,即2024/7/317

二维情形主成分的解释(续)

注意两个主成分

y1和

y2都是

x1和

x2的线性组合：其中P为旋转变换矩阵，它是正交矩阵．

极端情形1：椭圆变成圆，第一主成分

y1只体现了原始二维样品点约一半的信息，若此时将

y2忽略，则将损失约50％的信息，这显然是不可取的．

极端情形2：椭圆扁平到了极限，变成

y1轴上一条线段，第一主成分

y1几乎包含有二维样品点的全部信息，仅用

y1代替原始数据几乎不会有任何的信息损失，此时降维效果是非常理想的.2024/7/318

利用R程序来模拟这一过程（先下载安装mvtnorm）2024/7/319

>library(mvtnorm)

#先加载多元正态及t分布程序包mvtnorm>set.seed(8)

#设置随机数种子>sigma<-matrix(c(1,0.9,0.9,1),ncol=2)#设定协方差矩阵,相关系数为0.9>mnorm<-rmvnorm(n=200,mean=c(0,0),sigma=sigma)>plot(mnorm)#产生200个二维正态分布随机数并画散点图（见图9-2）>abline(a=0,b=1);abline(a=0,b=-1)#画坐标轴旋转45度后的二条直线>eig<-eigen(cor(mnorm));eig#求特征值和特征向量>vector1<-eig$vectors[,1];vector2<-eig$vectors[,2]>y1<-scale(mnorm)%*%vector1;y2<-scale(mnorm)%*%vector2#将数据标准化>plot(y1,y2,ylim=c(-2,2));abline(h=0,v=0)

#见图9-3>cbind(var(y1),var(y2),cor(y1,y2))

[,1][,2][,3][1,]1.885428

0.11457184.418324e-16#可见y1方差很大为1.885，y2方差相对很小为0.115，且y1和y2不相关

#可以认为原来二维数据的变化都体现在y1一个维度上了2024/7/3110

图9－2二维正态分布模拟数据的主成分示意图2024/7/3111

图9－3坐标轴旋转以后的散点图一般,设总体的p个主成分为：

第i

个主成分yi的方差为y1

是

X的一切线性组合中方差最大者；y2

是与y1

不相关是

X的一切线性组合中方差最大者；······2024/7/3112

9.2.2主成分的计算2024/7/3113

9.2.3主成分的主要性质2024/7/3114

设总体X的p个主成分所成向量为：性质1

：性质2：性质3：载荷与载荷矩阵：2024/7/3115

9.2.4主成分个数的确定2024/7/3116

第i

个主成分yi的方差贡献率为：通常取使得累积贡献率满足的最小的k为主成分个数．有的文献取累积贡献率首次超过85%的k．

9.2.5变量的标准化及意义从总体协方差矩阵Σ出发做主成分分析倾向于反映方差大的变量的信息，会出现“大数吃小数”的现象．为了均等地对待每一个原始变量，常常将各原始变量作标准化处理：标准化后的总体

的协方差矩阵恰好是原总体X的相关系数矩阵

ρ．综上讨论，既可从Σ出发，也可以从ρ出发做主成分分析，考虑到现实经济意义，后者用得更多.2024/7/3117

9.3样本主成分实际问题中Σ和ρ往往是未知的,需要用样本的协方差矩阵

S和样本的相关系数矩阵

R来估计：2024/7/3118

9.3.1样本主成分的性质和计算设S的p个特征值为，对应的单位正交特征向量为，则样本的第

i个主成分为性质1

：性质2

：性质3

：性质4：2024/7/3119

9.3.2主成分分析的步骤

实际问题中更常用的是从样本相关系数矩阵R

出发进行主成分分析,方法是用

R

替换

S,其余操作不变,其步骤可归纳为：将原始样本标准化求样本的相关系数矩阵

R

求R的特征值以及对应的单位正交特征向量按主成分累积贡献率超过80%确定主成分的个数k,并写出主成分表达式为(5)对分析结果做统计意义和实际意义两方面的解释.2024/7/3120

主成分分析特别说明：2024/7/3121

9.3.2

(续)主成分相关的R函数1.

princomp函数(这是主成分分析最常用的函数)princomp(x,cor=FALSE,score=TRUE,…)2.summary函数(提取主成分的信息)summary(object,loadings=FALSE,…)3.loadings函数(显示主成分的载荷阵)loadings(object)4.

predict函数(预测主成分的值)predict(object,newdata,…)5.

screeplot函数(画出主成分的碎石图)screeplot(object,type=c(“barplot”,“lines”,…)2024/7/3122

2024/7/3123表9-1给出了52名学生的数学

(x1)、物理

(x2)、化学

(x3)、语文

(x4)、历史

(x5)和英语

(x6)成绩，对其进行主成分分析．例9.1学生六门课成绩数据的主成分分析>setwd("C:/data")#设定工作路径>d9.1<-read.csv(“exam9.1.csv”,header=T)#读入数据>R=round(cor(d9.1),3);R#样本相关系数阵保留三位小数

x1x2x3x4x5x6x11.0000.6470.696-0.561-0.456-0.439x2

0.6471.0000.573-0.503-0.351-0.458x30.6960.5731.000

-0.380-0.274-0.244x4-0.561-0.503-0.3801.0000.8130.835x5-0.456-0.351-0.2740.8131.0000.819x6-0.439-0.458-0.244

0.8350.8191.000解

先读取数据，求样本相关系数矩阵，R程序如下：表6-017个地质勘探点样品的标准化数据2024/7/3124

学号x1x2x3x4x5x6165617284817927777766470553676349656757478847562716456671675265576831007941675078694975163558678453586656……………………………………4599100995363604678685275746647729073768079486964606874804952626510096100507072567482745172747588918652687470878783表

9－152名学生六门课程成绩数据2024/7/3125

易见，文科三门课程语文(x4)、历史(x5)和英语(x6)相关性较强；理科三门课程数学(x1)、物理(x2)和化学(x3)相关性也较强．可以进一步作主成分分析，求样本相关矩阵的特征值和主成分载荷．

由下面的R程序运行结果可知主成分的标准差，即相关系数矩阵的六个特征值开方各为：同时前两个主成分的累积贡献率为0.618+0.210=0.829，已经超过80%，所以取两个主成分就可以了．2024/7/3126>PCA9.1=princomp(d9.1,cor=T)

#用样本相关系数阵做主成分分析>PCA9.1Call:princomp(x=d9.1,cor=T)Standarddeviations:Comp.1Comp.2Comp.3Comp.4Comp.5Comp.61.926

1.1240.6640.5200.4120.3836variablesand52observations.>summary(PCA9.1,loadings=T)

#列出主成分分析结果Importanceofcomponents:

Comp.1Comp.2Comp.3Comp.4Comp.5Comp.6Standarddeviation1.926

1.124

0.6640.5200.4120.383ProportionofVariance

0.6180.2100.0730.0450.0280.024CumulativeProportion

0.6180.8290.9020.9470.9761.000Loadings:

Comp.1Comp.2Comp.3Comp.4Comp.5Comp.6x1-0.412-0.3760.2160.788-0.145x2-0.381-0.357-0.806

-0.1180.212

0.141x3-0.332-0.563

0.467-0.588x40.461-0.279

0.599-0.590x50.421-0.415-0.250-0.738

-0.205x60.430-0.4070.1460.1340.2220.749第一主成分和第二主成分分别为：参见教材分析，可将将它们分别理解为“课程差异主成分”

和“课程均衡主成分”.2024/7/31

27>round(predict(PCA9.1),3)

#作预测,计算主成分得分并解释>screeplot(PCA9.1,type=“lines”)

#画线型碎石图(见图9-4)>load=loadings(PCA9.1)

#提取主成分载荷矩阵为load>plot(load[,1:2],xlim=c(-0.6,0.6),ylim=c(-0.6,0.6))#作散点图>rnames=c(“数学”,“物理”,“化学”,“语文”,“历史”,“英语”)#命名>text(load[,1],load[,2],labels=rnames,adj=c(-0.3,1.5))#用中文为散点图标注>abline(h=0,v=0,lty=3)

#用虚线划分四个象限(见图9-5)2024/7/31

28图9－452名学生六门课程成绩的主成分线型碎石图2024/7/31

29图9－5前两个主成分的载荷散点图2024/7/31

30图9－652名学生成绩数据的双坐标散点图>biplot(PCA9.1,scale=0.5)

#绘制52个样本点关于前两个主成分的散点图2024/7/3131

由于第一主成分是文理课程差异因子，理科课程在第一主成分上的载荷绝对值大且取负值，文科课程在第一主成分上的载荷绝对值大且取正值，因此图中Comp.1轴方向靠左的样本点，如6，7和45号样本点，对应理科成绩好、文科成绩差的学生；相对的Comp.1轴方向靠右的样本点，如30和49号样本点，对应文科成绩好、理科成绩差的学生．又第二主成分表示课程均衡因子，在图中Comp.2轴方向靠下的样本点，如26，33号样本点，对应各科成绩都较好学生，相对的Comp.2轴方向靠上的样本点，如3，5和8号样本点，对应各科成绩都较差的学生，而居中的样本点，如42，24和39号样本点，对应各科成绩都属于中等且差异不大的学生．这样就可以对52名学生按对应样本点所在的位置进行大致分类．2024/7/3132

在某沉积盆地一坳陷区的17个取样点经勘探测定了六个地质变量：x1为有机碳（%）；x2为生油层埋深（米）；x3油层孔隙度（%）；x4为储层厚度（米）；x5为地下水含碘量（p.p.m）；x6为地下水矿化度(克/升)，见表9-2．要求根据这些数据进行主成分分析.例9.2石油勘探样品数据分析(数据exam9.2)表6-017个地质勘探点样品的标准化数据2024/7/3133

点号x1x2x3x4x5x61-0.9142-0.7119-0.9293-0.4385-0.57100.73612-0.3095-0.5206-1.3309-0.2764-0.57100.57143-1.0654-0.71190.2756-0.7626-1.09570.90074-1.3073-0.95111.25740.3718-1.09571.394650.1743-0.47270.3203-0.9895-0.0463-0.25186-0.8235-0.59230.40951.3441-0.83330.406870.90002.1583-0.1260-0.85981.7901-1.89838-0.0071-0.3532-1.4201-1.0219-0.0463-0.581191.20231.6799-0.7508-0.60052.3148-1.2397100.1743-0.3532-0.97391.3441-0.04630.2421112.26061.44070.72192.64050.7407-1.075012-1.4282-0.95110.0079-0.7950-1.09571.065313-0.3397-0.52062.1499-0.1144-0.57100.4068140.7790-0.23361.19700.69590.21610.9104150.41620.72321.0789-0.30880.47840.745716-0.6118-0.71190.36490.0477-0.57101.5593170.90001.08200.1418-0.27641.0031-0.5811表

9－217个地质勘探点样品的标准化数据2024/7/3134>setwd("C:/data")

#设定工作路径>d9.2<-read.csv("exam9.2.csv",header=T)#将exam9.2数据读入到d9.2中>R=round(cor(d9.2),3);R

#求样本相关系数矩阵

x1x2x3x4

x5x6x11.0000.8400.003

0.3470.839-0.747x20.8401.000-0.051

0.077

0.939-0.839x30.003-0.0511.0000.259-0.1640.285x40.3470.0770.259

1.000-0.0370.022x50.8390.939-0.164-0.0371.000-0.827x6-0.747-0.8390.2850.022-0.8271.000

易见，x2与x5相关性最强，其绝对值在0.9~0.95，x1与x2，x1与x5，x2与x6，x5与x6的相关性较强，其绝对值在0.8~0.9，说明六个变量之间确实存在较强的相关关系，应当进行“降维”处理，可以作主成分分析．2024/7/3135>options(digits=3)

#设置小数点位数为3>PCA9.2=princomp(d9.2,cor=T,scores=T);PCA9.2#作主成分分析Call:princomp(x=d9.2,cor=T,scores=T)Standarddeviations:Comp.1Comp.2Comp.3Comp.4Comp.5Comp.6

1.8851.1700.8600.4300.3400.197>summary(PCA9.2,loadings=T)

#列出主成分分析结果Importanceofcomponents:

Comp.1Comp.2Comp.3Comp.4Comp.5Comp.6Standarddeviation

1.8851.1700.8600.43010.33990.19653ProportionofVariance0.5920.2280.1230.03080.01930.00644CumulativeProportion0.592

0.8200.9430.97430.99361.00000Loadings:

Comp.1Comp.2Comp.3Comp.4Comp.5Comp.6x10.4850.2390.2910.7350.274x20.510-0.166-0.587

0.600x30.646-0.728-0.181x40.7020.640-0.254-0.153x50.509-0.1540.409-0.187-0.713x6-0.484

0.1590.837-0.1180.155前两个主成分的累积贡献率为0.592+0.228=0.82，已经超过80%，所以只需取两个主成分．第一主成分和第二主成分各为（为简明起见，样本主成分表达式中的所有“*”省略，以下同）：2024/7/3136

四个变量x1

(有机碳)，x2(生油层埋深)，x5(地下水含碘量)和x6(地下水矿化度)在主成分z1上载荷较大，故第一主成分z1可解释为“生油条件”主成分；第二主成分z2与x3(油层孔隙度)和x4(储层厚度)这两个变量关系特别密切，可解释为“储油条件”主成分．这样的分析结果与石油地质理论是相符合的．2024/7/3137>screeplot(PCA6.1,type=“lines”)#画碎石图，用直方图类型（见图9-7）图9-717个石油地质勘测点样本数据的主成分碎石图2024/7/3138用主成分载荷矩阵前两列数据作主成分载荷散点图(见图9-8)，R程序如下：．load=loadings(PCA9.2)#提取主成分载荷矩阵plot(load[,1:2],xlim=c(-0.5,1),ylim=c(-0.2,0.8))

#作散点图rnames=c(“x1有机碳”,“x2生油层埋深”,“x3油层孔隙度”,“x4储层厚度","x5地下水含碘量","x6地下水矿化度")#见图9－8text(load[,1],load[,2],labels=rnames,cex=0.8,adj=c(-0.1,0.6))

#用中文为散点标号abline(h=0,v=0,lty=3)

#用虚线划分象限六个变量在主成分z1和z2坐标面上的载荷散点图表明了两个主成分z1和z2具有明显的“生油”和“储油”倾向特征．2024/7/3139

图9-8两个主成分的载荷散点图2024/7/3140>A=round(PCA9.2$scores,3)

#计算主成分得分，取3位小数>B=round(apply(A[,1:2],1,crossprod),2)

#按行加总前2个主成分上的载荷平方>cbind(A,"综合得分“=B,"排名“=rank(B))

#按列合并主成分得分、综合得分和排名

Comp.1Comp.2Comp.3Comp.4Comp.5Comp.6综合得分

排名

[1,]

-1.333

-1.1500.6350.136-0.150

-0.015

3.10

9

[2,]

-0.798

-1.1751.0730.252

0.1530.261

2.024

[3,]

-1.901

-0.540

-0.429

-0.187-0.047

0.287

3.9112

[4,]

-2.421

0.956

-0.424

-0.124-0.402-0.125

6.7814

[5,]

-0.037

-0.640

-0.669

-0.357

0.749

-0.144

0.412

[6,]

-1.2091.0340.836

-0.626-0.472-0.190

2.536

[7,]3.560

-1.001

-0.862

-0.504-0.495

0.093

13.6816

[8,]0.265

-1.930

0.606

-0.359

0.514

0.0033.8011

[9,]3.474

-1.110

-0.2160.441-0.255-0.37513.30

15[10,]0.0420.2641.8140.143-0.052-0.2290.071[11,]3.0332.6361.183

-0.323

0.189

0.30416.1517[12,]

-2.269

-0.814

-0.245

-0.119-0.2010.0955.8113[13,]

-1.0471.250

-1.458

-0.573

0.258-0.144

2.667[14,]

-0.0401.500

-0.2940.678

0.464-0.181

2.255[15,]0.4510.589

-1.0980.725-0.161

0.258

0.553[16,]

-1.6920.285

-0.0300.668-0.095

0.000

2.948[17,]1.922

-0.152

-0.4200.129

0.003

0.102

3.7210

11号样本点综合排名最高，为17分；7号和9号排名次之，分别为16分和15分；之后样本点得分排名从高到低依次为4、12、3、8和17号．

利用函数biplot来绘制它们在z1和z2构成的坐标面

z1Oz2上的散点图，并且加入六个变量在同一坐标面

z1Oz2上的载荷散点图，得到所谓的“双坐标”散点图（见图9-9）．

借助该图可以对17个勘测样本点进行大致分类：11号样本点独居右上，它在

“生油”主成分z1和“储油”主成分z2上得分均高，应该首先重点关注．7、9号样本点相邻且最靠右，且在z1得分很高，可合为一类，次重点考虑；此外，在z1和z2上至少有一个得分较高的3、4、8、12和17号样本点也应该重点考察．这与上面的综合得分和排名一致．2024/7/3141

2024/7/3142

图9－917个石油地质勘测点样本数据的双坐标散点图>biplot(PCA9.2,scale=0.5)#绘制17个样本点和6个变量对z1和z2的散点图9.4案例:主成分回归分析

案例9.1(数据文件为case9.1)表9-3给出了2019年全国31个地区相关数据．它们分别为：货运量x1(万吨)，货物周转量x2(亿吨公里)，GDPx3(亿元)，人均GDPx4(元)，城镇居民人均可支配收入y(元)．根据这些数据做线性回归分析和主成分回归分析，并比较它们的异同.2024/7/3143

2024/7/3144

表9-32019年全国部分地区货运量、GDP、人均可支配收入等数据城市x1x2x3x4y北京22808108935371.316422073848.5天津50093266214104.39037146118.9河北2424451356335104.54634835737.7山西192192546617026.74572433262.4内蒙古182702458717212.56785240782.5辽宁178253892124909.55719139777.2吉林43193180311726.84347532299.2黑龙江50475161513612.73618330944.6………………………………贵州83402123516769.34643334404.2云南122727155223223.84794436237.7西藏40251541697.84890237410.0陕西154749348225793.26664936098.2甘肃6361024968718.33299532323.4青海149453982966.04898133830.3宁夏425116513748.55421734328.5新疆844231948135977解(1)先做线性回归分析，R程序及结果如下：>setwd("C:/data")

#设定工作路径>c9.1<-read.csv("case9.1.csv",header=T)

#将数据读入到c9.1中>options(digits=3)

#取三位有效数字>lmc9.1<-lm(y~1+x1+x2+x3+x4,data=c9.1)>summary(lmc9.1)从输出结果（见教材）可以看出，回归方程是非常显著的，R2为0.923，模型拟合效果很好，但x1、x2和x3的回归系数没有通过显著性检验(在0.05的显著性水平下)．回归方程为：然后作逐步回归，R程序及结果如下：>summary(step(lmc9.1))

回归方程和回归系数均显著，R2为0.921，逐步回归方程为：2024/7/3145

(2)再作主成分回归分析,先求样本相关系数阵：>R=round(cor(c9.1[,2:6]),3);R

#求样本相关系数矩阵发现x4与

y高度相关，

x1

,

x2

,

x3相关性较强，可用主成分降维>c9.1pr<-princomp(~x1+x2+x3+x4,data=c9.1,cor=T)#使用公式法>summary(c9.1pr,loadings=T)

前两个主成分累积贡献率已达88%，故选择前两个主成分>pre<-predict(c9.1pr)

#计算主成分得分>c9.1$z1<-pre[,1];c9.1$z2<-pre[,2]>lmpr<-lm(y~z1+z2,data=c9.1)

#做y关于主成分z1和z2的回归>summary(lmpr)输出结果显示：y关于两个主成分z1和z2的回归方程和三个回归系数均是非常显著的，R2为0.887，主成分回归方程为：2024/7/3146

(3)还利用主成分与原来自变量间的关系将主成分还原为原来的自变量，参见参考文献[4]．R程序及结果如下：>beta<-coef(lmpr);A<-loadings(c9.1pr)[,1:2]>x.bar<-c9.1pr$center;x.sd<-c9.1pr$scale>coef<-A%*%beta[2:3]/x.sd>beta0<-beta[1]-x.bar%*%coef>c(beta0,coef)[1]2.36e+04-2.86e-023.98e-014.64e-022.52e-01

由输出结果知主成分z1和z2还原为原始变量后的回归方程为：可将它和最初得到的回归方程进行比较．前者是从主成分回归方程（方程和回归系数均显著）变形而来，更合理，预测效果也更好．2024/7/3147

2024/7/31487/31/2024主编：费宇49

7/31/2024主编：费宇50第10章因子分析因子分析（factoranalysis）最早起源于KarlPearson和ChalesSpearman等人关于智力的定义和测量工作，因子分析的基本目的是，只要可能，就用少数几个潜在的不能观察的随机变量（称为因子）去描述许多个随机变量之间的协方差关系．从这点上看，因子分析与主成分分析有相似之处，但因子分析中的因子是不可观察的，也不必是相互正交的变量．因子分析可以视为主成分分析的一种推广，它的基本思想是：根据相关性大小把变量分组，使得组内的变量相关性较高，但不同组的变量相关性较低，则每组变量可以代表一个基本结构，称为因子，它反映已经观测到的相关性．因子分析可以用来研究变量之间的相关关系，称为R型因子分析；也可以用来研究样品之间的相关关系，称为Q型因子分析．二者虽然形式上有所不同，但数学处理上是一样的，所以本章只介绍R型因子分析．7/31/2024主编：费宇51第10章因子分析10.1正交因子模型10.2因子模型的估计10.3因子正交旋转10.4因子得分10.5因子分析小结10.6案例分析7/31/2024主编：费宇5210.1正交因子模型1.模型定义设p维随机向量的期望为,方差-协方差矩阵为Σ,

假定X线性地依赖于少数几个不可观测的随机变量f1,…,fm(m<p)和p个附加的方差源

ε1,…,εp,

一般称f1,…,fm为公因子，称ε1,…,εp为特殊因子，或误差.

1.模型定义7/31/2024主编：费宇53那么，因子模型为1.模型定义引入矩阵符号，记

那么因子模型(10.1)可以写为7/31/2024主编：费宇541.

模型定义其中aij称为第i个变量在第j个因子上的载荷,矩阵A称为载荷矩阵.我们假定

如果模型(10.2)满足假定(10.3)，则称该模型为正交因子模型，如果F的各个分量相关，即

不是单位阵，则相应的模型称为斜交因子模型，本书只讨论正交因子模型.7/31/2024主编：费宇551.

模型定义从正交因子模型容易求得X的协方差7/31/2024主编：费宇561.模型定义同样，容易求得由(10.4)可得该式说明xi的方差由两部分构成：m个公因子和一个特殊因子，其中表示第j个公因子对xi的方差贡献，而Фi是第i个特殊因子对xi的方差贡献，称之为特殊度.记

hi2=ai12+ai22+…+aim2，它表示m个公因子对变量xi的方差贡献总和，称之为第i个共同度，它是载荷矩阵A的第i行元素平方和.7/31/2024主编：费宇571.模型定义由(10.5)可得上式说明aij表示变量xi与公因子fj的协方差.另一方面，我们也可以考虑某个公因子fj对各个变量x1,…,xp的影响，采用来度量这个影响的大小，bj2是载荷矩阵A第j列元素的平方和，称之为公因子fj对各p个变量的方差贡献，bj2越大，表示fj对各p个变量的影响越大，它可以作为公因子fj重要性的一个度量.7/31/2024主编：费宇581.模型定义需要指出的是，当

m>1时，因子模型是不唯一的,设T为m×m正交矩阵，即TTT=TTT=I，模型(10.2)可改写为

式中，A*=AT,F*=TTF7/31/2024主编：费宇591.模型定义注意到即F*也满足(10.3)，显然因子F与F*有相同的统计性质，但相应的载荷矩阵A与A*是不相同的，但它们产生相同的方差-协方差矩阵Σ

，即7/31/2024主编：费宇601.模型定义一方面，因为F*=TTF，即F*是由F经正交变换得到,而A*=AT，即A*=(a*ij)是由A=(aij)经正交变换得到,另一方面，由(10.11)易知，变量xi的共同度为即正交变换不改变公因子的共同度.7/31/2024主编：费宇6110.2

因子模型的估计建立因子模型首先要估计载荷矩阵及特殊方差，常用的估计方法有主成分法、主因子法和极大似然法.1.

主成份法设Σ的特征值为λ1,

λ2,…,λp(λ1≥λ2≥…≥λp≥0)，e1,e2,…,ep为对应的标准正交化特征向量，那么Σ可以写为7/31/2024主编：费宇621.

主成份法这个分解是公因子个数为p，特殊因子方差为0的因子模型的方差-协方差矩阵结构形式，即虽然上式给出的Σ因子分析表达式是精确的，但实际应用中没有价值，因为因子分析的目的是要寻找少数m(m<p)个公因子解释原来p个变量的协方差结构，所以，采用主成分分析的思想，如果Σ的最后p-m个特征值很小，在(10.13)中略去λm+1em+1eTm+1+…+λpepeTp对Σ的贡献，7/31/2024主编：费宇631.

主成份法于是得这里假定了(10.2)中的特殊因子是可以在Σ的分解中忽略的，如果特殊因子不能忽略，那么它们的方差可以取Σ-AAT的对角元，7/31/2024主编：费宇641.

主成份法此时有其中7/31/2024主编：费宇651.

主成份法实际应用中Σ是未知的，通常用它的估计，即样本协方差矩阵S来代替，考虑到变量的量纲差别，往往需要将数据标准化，这样求得的样本方差-协方差矩阵就是原来数据的相关系数矩阵R，所以可以从R出发来估计因子载荷矩阵和特殊因子的方差.7/31/2024主编：费宇661.

主成份法设R的特征值为 ,

为对应的标准正交化特征向量，设

m<p，则由R出发因子模型的载荷矩阵的估计为特殊因子的方差фi的估计为7/31/2024主编：费宇671.

主成份法这时，共同度hi2的估计为变量xi与公因子fj协方差的估计为，公因子fj对各个变量的贡献bj2的估计为7/31/2024主编：费宇681.

主成份法那么，如何确定公因子数目m呢？可以仿照主成分分析的思想，比如寻找m使得来确定公因子数m.7/31/2024主编：费宇692.

主因子法假定原始向量X的各分量已作了标准化变换.如果其满足正交因子模型，则有7/31/2024主编：费宇70则称为X的约相关矩阵.其中，

R为X的相关矩阵.令2.

主因子法7/31/2024主编：费宇71

中的对角线元素是

，而不是1，非对角线元素和R中是完全一样的，并且是一个非负定矩阵.设是特殊方差的一个合适的初始估计，则约相关矩阵可估计为：2.

主因子法7/31/2024主编：费宇72其中，，是的初始估计.又设的前个特征值依次为，相应的正交单位特征向量为，则A的主因子解为：2.

主因子法7/31/2024主编：费宇73由此我们可以重新估计特殊方差，的最终估计为：如果我们希望求得拟合程度更高的解，则可以采用迭代的方法，即利用式(10.26)中的再作为特殊方差的初始估计，重复上述步骤，直至解稳定为止．2.

主因子法7/31/2024主编：费宇74特殊(或共性)方差的常用初始估计方法有：(1)取

，其中是的第个对角线元素，此时共性方差的估计为

，它是

和其他个变量间样本复相关系数的平方，该初始估计方法最为常用．(2)取

，此时.(3)取

，此时

，得到的是一个主成分解.3.

极大似然法7/31/2024主编：费宇75设公共因子，特殊因子，且相互独立，则必然有原始向量．由样本计算得到的似然函数是和的函数．由于，故似然函数可更清楚地表示为．记的极大似然估计为，即有可以证明，而和满足方程组：3.

极大似然法7/31/2024主编：费宇76式中，由于A的解是不唯一的，故为了得到唯一解，可附加计算上方便的唯一性条件：是对角矩阵.3.

极大似然法7/31/2024主编：费宇77方程组（10.28）的和一般可用迭代方法解得.对极大似然解，当因子数增加时，原来因子的估计载荷及对x的贡献将发生变化，这与主成分解及主因子解不同．10.3

因子正交旋转在第10.1节我们已经看到，满足方差结构Σ

=AAT+Ф的因子模型并不惟一，模型的公因子与载荷矩阵不惟一.如果F是模型的公因子，A是相应的载荷矩阵，而T是m×m正交矩阵，则F*=TTF也是公因子，相应的载荷矩阵为A*=AT，A*也满足Σ=A*A*T+Ф这说明，公因子和因子载荷矩阵作正交变换后，并不改变共同度，我们称因子载荷的正交变换和伴随的因子正交变换为因子正交旋转.7/31/2024主编：费宇7810.3

因子正交旋转设

是用某种方法（比如主成分法）得到的因子载荷矩阵的估计，T为

m×m正交阵,则是旋转载荷矩阵.问题是：为什么要进行因子旋转？其目的是什么？7/31/2024主编：费宇7910.3

因子正交旋转如果初始载荷不易解释时，就需要对载荷作旋转，以便得到一个更简单的结构.最理想的情况是这样的载荷结构，每个变量仅在一个因子上有较大的载荷，而在其余因子上的载荷比较小，至多是中等大小，这样公因子fi的具体含义可由载荷较大的变量根据具体问题加以解释.如何进行因子旋转寻找一个简单结构的载荷矩阵，这里不作详细介绍.7/31/2024主编：费宇8010.4

因子得分在因子分析中，虽然我们关心模型中载荷矩阵的估计和对公因子的解释，但对于公因子的估计，即因子得分，有时也是需要的.但是因子得分的计算并不同于通常意义下的参数估计，而是对不可观测的因子fj取值的估计，下面介绍用加权最小二乘法估计因子得分.7/31/2024主编：费宇811.

加权最小二乘法给定因子模型X=μ+AF+ε,假定均值向量μ,载荷矩阵A和特殊方差阵Ф已知，把特殊因子ε看作误差，因为Var(εi)=фi(i=1,2,…,p)未必相等,所以我们用加权最小二乘法估计公因子F.首先将因子模型

(10.2)改写为7/31/2024主编：费宇821.

加权最小二乘法两边左乘Ф-1/2得记X*=Ф-1/2(X-μ),A*=Ф-1/2A,ε*=Ф-1/2ε,则上式可以写成注意到E(ε*)=Ф-1/2E(ε)=0,

cov(ε*)=E(ε*ε*T)=Ф-1/2E(εεT)Ф-1/2=I7/31/2024主编：费宇831.

加权最小二乘法所以(10.32)是经典的回归模型，由最小二乘法知F的估计为

实际中,A,

Ф和μ都是未知的,通常用它们的某种估计来代替,比如我们采用正交旋转后的载荷矩阵A的估计,和样本均值 ,分别代替A,Ф和μ7/31/2024主编：费宇841.

加权最小二乘法于是可得对应于xj的因子得分7/31/2024主编：费宇852.

回归法在正交因子模型中，假设服从(m+p)元正态分布，用回归预测方法可将估计为：7/31/2024主编：费宇862.

回归法在实际应用中，可用，和分别代替上式中的，和来得到因子得分．样品的因子得分7/31/2024主编：费宇873.

综合因子得分7/31/2024主编：费宇88个因子任意若干个取相反符号，特别是全部取相反符号仍然满足因子分析模型，所以仍然可以作为因子。3.

综合因子得分以各因子的方差贡献率为权重，由各因子的线性组合得到综合评价指标函数：式中，7/31/2024主编：费宇893.

综合因子得分7/31/2024主编：费宇90那么这样的因子得分函数将会有种不同的组合。所以这样的因子得分实际上是不好解释的，此外，使用不同的因子旋转会得到不同的因子，从而综合评价函数也就不同，哪一个才是对的呢？还有，综合起来表示的是什么因子呢？所以，因子综合得分是没有合理的解释的。例10.1数据文件为eg9.1前面第9章例9.1表9-1给出了52名学生的数学（x1）、物理（x2）、化学（x3）、语文（x4）、历史（x5）和英语（x6）成绩，试进行学生成绩的因子分析.解:采用R软件对样本数据进行因子分析，首先计算样本数据的相关系数矩阵，观察各变量之间的相关性.

R程序及结果如下：7/31/2024主编：费宇91例10.1数据文件为eg9.1#假设已经读取了52名学生成绩数据cor(X)#计算样本数据的相关系数矩阵x1x2x3x4x5x6x11.000.650.70-0.56-0.46-0.44x20.651.000.57-0.50-0.35-0.46x30.700.571.00-0.38-0.27-0.24x4-0.56-0.50-0.381.000.810.83x5-0.46-0.35-0.270.811.000.82x6-0.44-0.46-0.240.830.821.007/31/2024主编：费宇92例10.1数据文件为eg9.1从样本数据各变量的相关系数上可以看出,x4、x5和x6之间存在较强的相关性.为了消除各变量之间的相关性，下面分别采用R软件中基于极大似然法的因子分析函数factanal()和基于主成分法的因子分析函数factpc()对数据进行因子分析提取因子.

R程序及结果如下：7/31/2024主编：费宇93例10.1数据文件为eg9.1#极大似然法做因子分析factanal(X,factors=2,rotation="none")Call:factanal(x=X,factors=2,rotation="none")Uniquenesses:x1x2x3x4x5x60.230.460.330.150.210.15Loadings:Factor1Factor2x1-0.680.56x2-0.600.43x3-0.490.66x40.920.10x50.860.24x60.880.277/31/2024主编：费宇94例10.1数据文件为eg9.1

Factor1Factor2SSloadings3.401.07ProportionVar0.570.18CumulativeVar0.570.74

Testofthehypothesisthat2factorsaresufficient.Thechisquarestatisticis3.6on4degreesoffreedom.Thep-valueis0.46#主成分法做因子分析library(mvstats)#加载mvstats包fac=factpc(X,2)fac7/31/2024主编：费宇95例10.1数据文件为eg9.1$VarsVarsVars.PropVars.CumFactor13.7100.618361.83Factor21.2620.210482.87$loadingsFactor1Factor2X1-0.79370.4224x2-0.73420.4008x3-0.63970.6322x40.88830.3129x50.81010.4661x60.82850.45677/31/2024主编：费宇96例10.1

数据文件为eg9.1从上述极大似然法和主成分法得出的因子分析结果上可以看出，极大似然法前两个因子累计贡献率只有74%，而主成分法累计贡献率达到了82.87%，说明主成分法效果比极大似然分析法效果好，其原因在于，极大似然法做因子分析要求数据样本要服从多元正态分布，但在实际中大多数数据都很难满足多元正态要求。接下来为了更好地解释因子的含义，我们基于主成分法采用方差最大化作因子正交旋转。R程序及结果如下：7/31/2024主编：费