本文深入解析主成分分析(PCA)与因子分析(EFA)的原理与应用，涵盖数据降维、潜在结构探索、主成分提取与旋转、因子模型选择及得分计算等关键步骤。通过实际案例演示，帮助读者掌握复杂数据集的简化与解读技巧。
#主成分分析和因子分析
#主成分分析，PCA
#一种数据降维技巧，它能将大量相关变量转换为一组很少的不相关变量，这些无关变量称为主成分
#探索性因子分析，EFA
#一系列用来发现一组变量的潜在结构的方法。
#它通过寻找一组更小的，潜在的或隐藏的结构来解释已观测到的，显示的变量间的关系
#主成分，因子模型的区别
#主成分（PC1和PC2）是观测变量的线性组合
#形成线性组合的权重都是通过最大化各主成分所解释的方差来获得，同时还要保证各主成分间不相关
#因子(F1和F2)被当做观测变量的结构基础或“原因”，而不是它们的线性组合
#代表观测变量方差的误差（e1到e5)无法用因子来解释
#图中的圆圈表示因子和误差无法直接观测，但是可通过变量间的的相互关系推导得到
#分析中常见步骤
#第一步：数据预处理，注意确保数据中没有确实值
#第二步：选择因子模型，判断PCA(数据降维)，还是EFA(发现潜在结构)更符合你的研究目标
#第三步：判断要选择的主成分、因子数目
#第四步：选择主成分，因子
#第五步：旋转主成分，因子
#第六步：解释结果
#第七步：计算主成分，因子得分
#主成分分析
#PCA的目标是用一组较少的不相关变量代替大量相关变量，同时尽可能保留初始变量的信息，
#这些推导所得的变量称为主成分
#如第一主成分：pc1=ax1+bx2+cx3+.....
#它是k个观测变量的加权组合，对初始变量及的方差解释性最大
#第二主成分也是初始变量的线性组合，对方差的解释性排第二，同时与第一主成分正交（不相关）
#后面每一个主成分都最大化它对花茶的解释程度，同时与之前所有的主成分都正交
#范例一
#第一步，准备数据源
library(psych) #
USJudgeRatings[,1] #仅显示数据框的第一列数据
USJudgeRatings[,-1] #去除第一列数据，显示其他列数据
#第二步，判断主成分个数
#原则一: 根据先验经验和理论知识判断主成分数；
#原则二：根据要解释变量方差的积累值的阈值来判断需要的主成分数
#原则三：通过检查变量间k*k的相关系数矩阵来判断保留的主成分数
#方法一
#基于特征值
#每个主成分都与相关系数矩阵的特征值相关联
#第一主成分与最大特征值相关联
#第二主成分与第二大特征值相关联,依次类推
#方法二
#Kaiser-Harris准则建议保留特征值大于1的主成分，
#特征值小于的1的成分所解释的方差比包含在单个变量中的方差更少
#方法三
#Cattell碎石检验原则绘制了特征值与主成分数的图形
#这类图形可以清晰地展示图形弯曲状况，在图形变化最大处之上的主成分都可保留
#方法四
#平行分析
#进行模拟，依据与初始矩阵相同大小的随机数据矩阵来判断要提取的特征值
#若基于真实数据的某个特征值大于一组随机数据矩阵响应的平均特征值，那么该组成分可以保留
#fa.parallel() 含平行分析的碎石图
#fa
#show the eigen values for a principal components (fa="pc") 显示主组件的特征值
#or a principal axis factor analysis (fa="fa") 主轴因子分析
#or both principal components and principal factors (fa="both")
#n.iter Number of simulated analyses to perform 模拟次数
fa.parallel(USJudgeRatings[,-1],fa="pc",n.iter=100,
show.legend=FALSE,main="Scree plot with parallel analysis")
#结果解释
#Parallel analysis suggests that the number of factors = NA
#and the number of components = 1
#由线段和x符号组成，基于观测特征值的碎石检验
#根据100个随机数据矩阵推导出来的特征值均值（虚线）
#缺一根线：y=1的水平线，对应大于1的特征值准则
#三种准则表明选择一个主成分即可保留数据集的大部分信息
#提取主成分
#范例一
#principal(r,nfactors=,rotate=,scores=)
#r,相关系数矩阵或原始数据矩阵
#nfactors,设定当前所需的主成分数，默认为1
#rotate指定旋转的方法（默认最大方差旋转varimax）
#scores设定是否需要计算主成分得分，默认不需要
pc <- principal(USJudgeRatings[,-1],nfactors = 1)
pc
#PC1栏
#包含了成分载荷，指观测变量与主成分的相关系数
#PC1与每一个变量都高度相关，也就是说，它是一个可用来进行一般性评价的维度
#如果提取不止一个主成分，那么还将会有PC2,PC3等栏
#成分载荷(component loadings)，可用来解释主成分的含义
#h2栏
#公因子方差，即主成分对每个变量的方差解释度
#u2栏
#成分唯一性，即方差无法被主成分解释的比例
#例如
# PC1 h2 u2 com
#PHYS 0.89 0.80 0.2013 1
#PHYS体能，能被PC1解释的方差为h2(80%),不能被解释的方差为u2(20%)
#PHYS与PC1的相关系数，在PC1列中最小
#所以PHYS是用第一主成分表示性最差的变量
#SS loadings 10.13
#包含了与主成分相关联的特征值，指的是与特定主成分相关联的标准化后的方差值，约等于10
#Proportion Var 0.92
#表示了每个主成分对整个数据集的解释程度
#第一主成分解释了11个变量92%的方差
#范例二
library(psych)
Harman23.cor #数据源
class(Harman23.cor) #返回"list"
Harman23.cor$cov #列表中的cov对象，做为数据源
#n.obs,指定样本大小
fa.parallel(Harman23.cor$cov,n.obs=302,fa="pc",
n.iter=100,show.legend=FALSE,main="Scree plot with parallel analysis")
colors() #返回所有可用颜色变量名
abline(h=1,col="orange") #画特征值为1的参考线
#直线和x符号，对应碎石图
#100个随机数据矩阵推导出来的特征值均值虚线，对应虚线
#结论：
#图形中Kaiser-Harris准则，碎石检验和平行建议都建议选择两个主成分
#范例二
library(psych)
#nfactors=2，提取两个主成分
#rotate="none"，
pc <- principal(Harman23.cor$cov,nfactors=2,rotate="none")
pc
#结果解释
# PC1 PC2
#Proportion Var 0.58 0.22
#PC1第一主成分对整个数据集的解释程度为58%
#PC2第二主成分对整个数据集的解释程度为22%
#两者总共解释了81%的方差(88%=58%+22%)
#PC1列，即主成分1与每个变量的相关系数
#PC2列，即主成分2与每个变量的相关系数
# PC1 PC2 h2 u2
#height 0.86 -0.37 0.88 0.123
#height与PC1的相关系数为0.86
#height与pc2的相关系数为-0.37
#第一主成分，第二主成分，总共对变量height的解释度为88%(h2)，未解释度为12%
#看向PC1列的相关系数，显示第一主成分与每一个身体测量指标都正相关，
#看起来似乎是一个一般性的衡量因子
#看向PC2列的相关系数，显示第二主成分与前四个变量(height,arm.span,forearm,lower.leg)负相关，
#与后四个变量(weight,bitro.diameter,chest.girth,chest.width)正相关
#因此，它似乎是是一个长度-容量因子
#主成分旋转
#旋转是一系列将成分载荷阵变得更容易解释的数学方法，它们尽可能的对成分去噪
#两种旋转方法：
#正交旋转，使选择的成分保持不相关
#最流行的正交旋转是方差极大旋转，它试图对载荷阵的列进行去噪，使得每个成分只由
#一组有限的变量来解释，即载荷阵每列只有少数几个很大的载荷，其他都是很小的载荷
#斜交旋转，让它们变的相关
#主成分旋转
rc <- principal(Harman23.cor$cov,nfactors=2,rotate="varimax")
rc
#结果解释
#列名从PC，转变成RC，以表示成分被旋转
#注意两个主成分之间仍然不相关
#Proportion Var 0.44 0.37
#对变量的解释性不变：44%+37%=81%
#这是因为变量的群组没有发生变化
#变的只是各个主成分对整体方差的解释度(成分1从58%变成44%，成分2从22%变成37%)
#观察RC1列的载荷，可以发现第一主成分主要由前四个变量来解释，
#包含(height,arm.span,forearm,lower.leg)，统称为长度变量
#观察RC2列的载荷，可以发现第二主成分主要由后四个变量来解释，
#包含(weight,bitro.diameter,chest.dirth,chest.width)，统称为容量变量
#范例一
#获取主成分得分
library(psych)
#获得每一个调查对象在当前第一主成分上的得分
dim(USJudgeRatings[,-1]) #返回43 11，43行，11列
pc <- principal(USJudgeRatings[,-1],nfactors=1,score=TRUE)
dim(pc$scores) #返回43 1,43行，1列
head(pc$scores) #查看前6行调查对象的得分情况
##当主成分分析基于相关系数矩阵时，原始数据便不可用了，
#也不可能获取每个观测的主成分得分
#但是你可以得到用来计算主成分得分的系数
#范例二
library(psych)
rc <- principal(Harman23.cor$cov,nfactors=2,rotate="varimax") #提取主成分
x <- rc$weights #从主成分模型中提取第一主成分载荷列，第二主成分载荷列
class(x) #返回"matrix"
#unclass(）可以消除对象的类,目前没懂到底啥用处
unclass(x) #仍旧为"matirx"类型
round(x,2) #对x中所有小数，仅保留两位
#结果解释
#利用RC1列载荷(相关系数)*变量名，生成第一主成分
#PC1=0.28*height+0.30*arm.span+0.30*forearm+0.28*lower.leg
#-0.06*weight-0.08*bitro.diameter-0.10*chest.girth-0.04*chest.width
#利用RC2列载荷(相关系数)*变量名，生成第二主成分
#PC2=-0.05*height-0.08*arm.span-0.09*forearm-0.06*lower.leg+
#0.33*weight+0.32*bitro.diameter+0.34*chest.girth+0.27*chest.width
#注意，两个等式都假定身体测量指标都已标准化，mean=0,sd=1
#探索性因子分析(EFA)
#目标是通过发掘隐藏在数据下的一组较少的，更为基本的无法观测的变量来解释一个可观测变量的相关性。
#这些虚拟的，无法观测的变量称为因子
#每个因子被认为可解释多个观测变量间共有的方差，因此准确的来说，它们应该称做公共因子
#范例一
options(digits=2) #全局设定，小数保留两位小数
covariances <- ability.cov$cov #读取协方差矩阵
correlations <- cov2cor(covariances) #将协方差矩阵转换成相关系数矩阵
#判断需提取的公共因字数
fa.parallel(correlations,n.obs=112,fa="both",n.iter=100,
main="Scree plots with parallel analysis")
#PCA
#PC Actual Data PC碎石图(特征值与主成分数的图形)
#PC Simulated Data PC100模拟的平行分析
#结论：PCA建议提取一个或者两个成分
#EFA
#对于EFA,Kaiser-Harris准则的特征值大于0，而不是大于1
#FA Actual Data FA碎石图(特征值与主成分数的图形)
#FA Simulated Data FA100模拟的平行分析
#结论：EFA建议提取两个因子
#当摇摆不定时，高估因子数通常比低估因子数的结果好，因为高估因字数一般较少
#曲解“真实”情况
#提取公共因子，因子未旋转
#fa(r,nfactors=,n.obs=,rotate=,scores=,fm=)
#r 相关系数，或者原始数据矩阵
#nfactors 设定需要提取的因子数量，默认为1
#n.obs 观测值，输入相关系数矩阵时必须填写
#rotate 设定旋转的方法，默认互变异数最小法
#scores 设定是否计算因子得分(默认不计算)
#fm 设定因子化方法，默认极小残差法
#pa 主轴迭代法
#ml 最大似然法
#wls 加权最小二乘法
#gls 广义加权最小二乘法
#minres 最小残差法
fa <- fa(correlations,nfactors=2,rotate="none",fm="pa")
fa
#结果解释
# PA1 PA2
#Proportion Var 0.46 0.14
#两个因子共解释了：0.46+0.14=60%的方差
#因子旋转
#提取公共因子，因子旋转
#fa(r,nfactors=,n.obs=,rotate=,scores=,fm=)
#r 相关系数，或者原始数据矩阵
#nfactors 设定需要提取的因子数量，默认为1
#n.obs 观测值，输入相关系数矩阵时必须填写
#rotate 设定旋转的方法，默认互变异数最小法
#scores 设定是否计算因子得分(默认不计算)
#fm 设定因子化方法，默认极小残差法
#pa 主轴迭代法
#ml 最大似然法
#wls 加权最小二乘法
#gls 广义加权最小二乘法
#minres 最小残差法
#范例一 正交旋转因子
fa.varimax <- fa(correlations,nfactors=2,rotate="varimax",fm="pa")
fa.varimax
class(fa.varimax) #返回"psych" "fa"
fa.varimax$Phi #返回Null,此项值作为正交旋转，和斜交旋转的唯一区别
#结果解释：
#正交旋转将人为的强制两个因子不相关
#看向PA1载荷列，reading,vocab在第一因子上载荷较大
#看向PA2载荷列，picture,blocks,maze在第二因子上载荷较大
#general 在两个因子上载荷较为平均，这表明存在一个语言智力因子和一个非语言智力因子
#范例二 斜交旋转因子
install.packages("GPArotation")
library("GPArotation")
fa.promax <- fa(correlations,nfactors=2,rotate="promax",fm="pa")
class(fa.promax) #返回"psych" "fa"
fa.promax
#结果解释：
#比对正交旋转，斜交旋转
#对于正交旋转，因子分析的重点在于因子结构矩阵
#而对于斜交矩阵，因子分析会考虑三个矩阵:因子结构矩阵，因子模式矩阵和因子关联矩阵
#因子结构矩阵，变量与因子的相关系数
#正交旋转，直接得到因子结构矩阵
#斜交旋转，需要根据，换算出因子结构矩阵
#因子模式矩阵，即标准化的回归系数矩阵，它列出了因子预测变量的权重。
#PA1列，PA2列组成因子模式矩阵
#正交旋转，没有
#斜交旋转，有
#With factor correlations of
#因子关联矩阵，即因子相关系数矩阵
#正交矩阵没有
#斜交矩阵有
#PA1与PA2的相关系数为0.55，相关性很大
#如果因子间的关联性很低，你可能需要重新使用正交旋转来简化问题
#因子结构矩阵，或称为因子载荷阵
#范例
#class(fa.promax)[2] #返回"fa"
#fa.promax$Phi #返回因子关联矩阵
#class(fa.promax$loading) #原始为loadings对象
#class(unclass(fa.promax$loading)) #unclass之后，变身为matirx对象
#unclass(fa.promax$loading) %*% fa.promax$Phi 因子结构矩阵(因子载荷矩阵)=因子模式矩阵*因子关联矩阵
fsm <- function(oblique){
if(class(oblique)[2]=="fa" & is.null(oblique$Phi)){
warning("Object doesn't look like oblique EFA")
}else{
P <- unclass(oblique$loading) #将因子模式对象读成矩阵
F <- P %*% oblique$Phi #因子结构矩阵(因子载荷矩阵)=因子模式矩阵*因子关联矩阵
colnames(F) <- c("PA1","PA2")
return(F)
}
}
fsm(fa.promax) #传入斜交旋转提取因子后的对象oblique
#结论，虽然斜交旋转更为复杂，但是模型更符合真实数据，因为PA1,与PA2的相关系数0.55，相关性比较高
#绘制正交或者斜交结果的图形
factor.plot(fa.promax,labels=rownames(fa.promax$loadings))
#图形说明
#PA1,第一载荷因子，作为横坐标
#vocab,reading在第一因子上载荷较大
#general在两个因子上较为平均
#PA2,第二载荷因子，作为纵坐标
#blocks,picture,maze在第二载荷因子上载荷较大
#general在两个因子上较为平均
#按照因子的维度，连线所有重要的载荷，多了一条PA1到general
fa.diagram(fa.promax,simple=FALSE)
#当simple=TRUE,那么仅将显示每个因子下最大的载荷，以及因子间的相关系数
fa.diagram(fa.promax,simple=TRUE)
#因子得分
fa.promax$weights
#与可精确计算的主成分得分不同，因子得分只是估计得到的
扫一扫
3368
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
