《R语言实战》第二版碎碎念之五:主成分分析与因子分析

主成分分析(PCA)是一种数据降维技巧,能将大量相关变量转化为一组很少的不相关变量,这些无关变量称为主成分,并能尽可能的保留原始数据集的信息。

探索性因子分析(EFA)是一系列用来发现一组变量的潜在结构的方法。它通过寻找一组更小的、潜在的或隐藏的结构来解释已观测到的、显式的变量间的关系。经验表明,因子分析需要5-10倍于变量数的样本数。

R中的基础包princomp()函数提供PCA分析,factanal()函数提供EFA分析。psych包的函数更丰富。

principal() 含多种可选的方差旋转方法的主成分分析

fa() 可用主轴、最小残差、加权最小平方或最大似然法估计的因子分析

fa.parallel() 含平行分析的碎石图

factor.plot() 绘制因子分析或主成分分析的结果

fa.diagram() 绘制因子分析或主成分的载荷矩阵

scree() 因子分析和主成分分析的碎石图

PCA和EFA都根据观测变量间的相关性来推导结果。用户可以输入原始数据矩阵或者相关系数矩阵到principal()和fa()函数中。若输入初始数据,相关系数矩阵将会被自动计算,在计算前请确保数据中没有缺失值。

(1)主成分分析

①首先判断主成分的个数: 根据先验经验和理论知识判断主成分数,根据要解释变量方差的积累值的阈值来判断需要的主成分数, 通过检查变量间k×k的相关系数矩阵来判断保留的主成分数。利用fa.parallel()函数:

library(psych)

fa.parallel(df, fa="pc", n.iter=100,show.legend=FALSE, main="Scree plot with parallel analysis")

代码生成图形展示基于观测特征值的碎石检验(由线段和x符号组成) 、根据100个( n.iter=100)随机数据矩阵推导出来的特征值均值(虚线),以及大于1的特征值准则(y=1的水平线)。

②提取主成分:principal()函数可以根据原始数据矩阵或者相关系数矩阵做主成分分析。格式为:principal(r, nfactors=, rotate=, scores=),其中:

 r是相关系数矩阵或原始数据矩阵;

 nfactors设定主成分数(默认为1);

 rotate指定旋转的方法(默认最大方差旋转(varimax));

 scores设定是否需要计算主成分得分(默认不需要)。

library(psych)

pc <- principal(USJudgeRatings[,-1], nfactors=1)

pc

Principal Components Analysis

Call: principal(r = USJudgeRatings[, -1], nfactors=1)

Standardized loadings based upon correlation matrix

PC1 h2 u2

INTG 0.92 0.84 0.157

DMNR 0.91 0.83 0.166

DILG 0.97 0.94 0.061

CFMG 0.96 0.93 0.072

。。。

      PC1

SS loadings 10.13

Proportion Var 0.92

[……已删除额外输出…… ]

PC1栏包含了成分载荷,指观测变量与主成分的相关系数。如果提取不止一个主成分,那么还将会有PC2、 PC3等栏。成分载荷(component loadings)可用来解释主成分的含义。此处可以看到,第一主成分(PC1)与每个变量都高度相关,也就是说,它是一个可用来进行一般性评价的维度。h2栏指成分公因子方差,即主成分对每个变量的方差解释度。 u2栏指成分唯一性,即方差无法被主成分解释的比例(1–h2)。例如,体能(PHYS) 80%的方差都可用第一主成分来解释,20%不能。相比而言, PHYS是用第一主成分表示性最差的变量。SS loadings行包含了与主成分相关联的特征值,指的是与特定主成分相关联的标准化后的方差值(本例中,第一主成分的值为10) 。最后, Proportion Var行表示的是每个主成分对整个数据集的解释程度。此处可以看到,第一主成分解释了11个变量92%。

③主成分旋转

旋转是一系列将成分载荷阵变得更容易解释的数学方法,它们尽可能地对成分去噪。旋转方法有两种:使选择的成分保持不相关(正交旋转),和让它们变得相关(斜交旋转)。旋转方法也会依据去噪定义的不同而不同。最流行的正交旋转是方差极大旋转,它试图对载荷阵的列进行去噪,使得每个成分只由一组有限的变量来解释(即载荷阵每列只有少数几个很大的载荷,其他都是很小的载荷)。

rc <- principal(Harman23.cor$cov, nfactors=2, rotate="varimax")

rc

Principal Components Analysis

Call: principal(r = Harman23.cor$cov, nfactors = 2, rotate = "varimax")

Standardized loadings based upon correlation matrix

RC1 RC2 h2 u2

height 0.90 0.25 0.88 0.123

arm.span 0.93 0.19 0.90 0.097

forearm 0.92 0.16 0.87 0.128

lower.leg 0.90 0.22 0.86 0.139

。。。。。

                          RC1 RC2

SS loadings 3.52 2.92

Proportion Var 0.44 0.37

Cumulative Var 0.44 0.81

[……已删除额外输出…… ]

列的名字都从PC变成了RC,以表示成分被旋转。

④获得主成分得分

library(psych)

pc <- principal(USJudgeRatings[,-1], nfactors=1, score=TRUE)

head(pc$scores)

PC1

AARONSON,L.H. -0.1857981

ALEXANDER,J.M. 0.7469865

ARMENTANO,A.J. 0.0704772

BERDON,R.I. 1.1358765

。。。。

当scores = TRUE时,主成分得分存储在principal()函数返回对象的scores元素中。如果有需要,你还可以获得律师与法官的接触频数与法官评分间的相关系数:cor(USJudgeRatings$CONT, pc$score)

PC1

[1,] -0.008815895

显然,律师与法官的熟稔度与律师的评分毫无关联。

(2)探索性因子分析

EFA的目标是通过发掘隐藏在数据下的一组较少的、更为基本的无法观测的变量,来解释一组可观测变量的相关性。这些虚拟的、无法观测的变量称作因子。(每个因子被认为可解释多个观测变量间共有的方差,因此准确来说,它们应该称作公共因子。)

用fa.parallel()函数可判断需提取的因子数。如果使用PCA方法,你可能会选择一个成分(碎石检验和平行分析)或者两个成分(特征值大于1)。当摇摆不定时,高估因子数通常比低估因子数的结果好,因为高估因子数一般较少曲解“真实”情况。

提取公共因子。fa(r, nfactors=, n.obs=, rotate=, scores=, fm=),其中:

 r是相关系数矩阵或者原始数据矩阵;

 nfactors设定提取的因子数(默认为1);

 n.obs是观测数(输入相关系数矩阵时需要填写);

 rotate设定旋转的方法(默认互变异数最小法);

 scores设定是否计算因子得分(默认不计算);

 fm设定因子化方法(默认极小残差法)。

与PCA不同,提取公共因子的方法很多,包括最大似然法(ml) 、主轴迭代法(pa) 、加权最小二乘法(wls)、广义加权最小二乘法(gls)和最小残差法(minres)。统计学家青睐使用最大似然法,因为它有良好的统计性质。不过有时候最大似然法不会收敛,此时使用主轴迭代法效果会很好。

因子正交旋转fa.varimax <- fa(correlations, nfactors=2, rotate="varimax", fm="pa")

斜交旋转fa.promax <- fa(correlations, nfactors=2, rotate="promax", fm="pa")

对于正交旋转,因子分析的重点在于因子结构矩阵(变量与因子的相关系数),而对于斜交旋转,因子分析会考虑三个矩阵:因子结构矩阵、因子模式矩阵和因子关联矩阵。因子模式矩阵即标准化的回归系数矩阵。它列出了因子预测变量的权重。 因子关联矩阵即因子相关系数矩阵

发表评论

电子邮件地址不会被公开。 必填项已用*标注