主成分分析
主成分分析（PCA）是一种数据降维技巧，它能将大量相关变量转化为一组很少的不相关变量，这些无关变量称为主成分。例如，使用PCA可将30个相关（很可能冗余）的环境变量转化为5个无关的成分变量，并且尽可能地保留原始数据集的信息。
1 R中的主成分和因子分析
R的基础安装包提供了PCA和EFA的函数，分别为princomp()和factanal()。
下表列出了psych包中相关度最高的函数：
每种方法都需要一些步骤（和决策）才能获得最终结果。最常见的步骤如下：
(1) 数据预处理。PCA和EFA都根据观测变量间的相关性来推导结果。用户可以输入原始数据矩阵或者相关系数矩阵到principal()和fa()函数中。若输入初始数据，相关系数矩阵将会
被自动计算，在计算前请确保数据中没有缺失值。
(2) 选择因子模型。判断是PCA（数据降维）还是EFA（发现潜在结构）更符合你的研究目
标。如果选择EFA方法，你还需要选择一种估计因子模型的方法（如最大似然估计）。
(3) 判断要选择的主成分/因子数目。
(4) 选择主成分/因子。
(5) 旋转主成分/因子。
(6) 解释结果。
(7) 计算主成分或因子得分。
2 主成分分析
PCA的目标是用一组较少的不相关变量代替大量相关变量，同时尽可能保留初始变量的信息，这些推导所得的变量称为主成分，它们是观测变量的线性组合。如第一主成分为：
它是k个观测变量的加权组合，对初始变量集的方差解释性最大。第二主成分也是初始变量的线性组合，对方差的解释性排第二，同时与第一主成分正交（不相关）。后面每一个主成分都最大化它对方差的解释程度，同时与之前所有的主成分都正交。
数据集USJudgeRatings包含了律师对美国高等法院法官的评分。数据框包含43个观测，12个变量。下表列出了所有的变量：
从实用的角度来看，你是否能够用较少的变量来总结这11个变量（从INTG到RTEN）评估的信息呢？如果可以，需要多少个？如何对它们进行定义呢？因为我们的目标是简化数据，所以可使用PCA。数据保持初始得分的格式，没有缺失值。因此，下一步便是判断需要多少个主成分。
2.1 判断主成分的个数
以下是一些可用来判断PCA中需要多少个主成分的准则：
根据先验经验和理论知识判断主成分数；
根据要解释变量方差的积累值的阈值来判断需要的主成分数；
通过检查变量间k×k的相关系数矩阵来判断保留的主成分数。
最常见的是基于特征值的方法。每个主成分都与相关系数矩阵的特征值相关联，第一主成分
与最大的特征值相关联，第二主成分与第二大的特征值相关联，依此类推。Kaiser-Harris准则建议保留特征值大于1的主成分，特征值小于1的成分所解释的方差比包含在单个变量中的方差更少。Cattell碎石检验则绘制了特征值与主成分数的图形。这类图形可以清晰地展示图形弯曲状况，在图形变化最大处之上的主成分都可保留。最后，你还可以进行模拟，依据与初始矩阵相同大小的随机数据矩阵来判断要提取的特征值。若基于真实数据的某个特征值大于一组随机数据矩阵相应的平均特征值，那么该主成分可以保留，该方法称作平行分析。
利用fa.parallel()函数，你可以同时对三种特征值判别准则进行评价。对于11种评分（删去了CONT变量），代码如下：
library(psych)
fa.parallel(USJudgeRatings[,-1], fa="pc", n.iter=100,
show.legend=FALSE, main="Scree plot with parallel analysis")
#显示主成分（fa="pc"），n.iter要进行的模拟分析的数量
结果分析：碎石图（直线与x符号）、特征值大于1准则（水平线）和100次模拟的平行分析（虚线）都表明保留一个主成分即可。下一步是使用principal()函数挑选出相应的主成分。
2.2 提取主成分
principal()函数可以根据原始数据矩阵或者相关系数矩阵做主成分分析。格式为：
principal(r, nfactors=, rotate=, scores=)
其中：
q r是相关系数矩阵或原始数据矩阵；
q nfactors设定主成分数（默认为1）；
q rotate指定旋转的方法（默认最大方差旋转（varimax））；
q scores设定是否需要计算主成分得分（默认不需要）。
（1）例题1
使用如下代码可获取第一主成分：
library(psych)
pc <- principal(USJudgeRatings[,-1], nfactors=1) #nfactors设定主成分数（默认为1）
pc
结果分析：此处，你输入的是没有CONT变量的原始数据，并指定获取一个未旋转的主成分。由于PCA只对相关系数矩阵进行分析，在获取主成分前，原始数据将会被自动转换为相关系数矩阵。PC1栏包含了成分载荷，指观测变量与主成分的相关系数。如果提取不止一个主成分，那么还将会有PC2、PC3等栏。成分载荷（component loadings）可用来解释主成分的含义。此处可以看到，第一主成分（PC1）与每个变量都高度相关，也就是说，它是一个可用来进行一般性评价的维度。
h2栏指成分公因子方差，即主成分对每个变量的方差解释度。u2栏指成分唯一性，即方差
无法被主成分解释的比例（1–h2）。例如，体能（PHYS）80%的方差都可用第一主成分来解释，20%不能。相比而言，PHYS是用第一主成分表示性最差的变量。
SS loadings行包含了与主成分相关联的特征值，指的是与特定主成分相关联的标准化后的
方差值（本例中，第一主成分的值为10）。最后，Proportion Var行表示的是每个主成分对整个数据集的解释程度。此处可以看到，第一主成分解释了11个变量92%的方差。
（2）例题2
看看第二个例子，它的结果不止一个主成分。Harman23.cor数据集包含305个女孩的8个身体测量指标。本例中，数据集由变量的相关系数组成，而不是原始数据集，见下表：
我们希望用较少的变量替换这些原始身体指标。如下代码可判断要提取的主成分数。此处，你需要填入相关系数矩阵（Harman23.cor对象中的cov部分），并设定样本大小（n.obs）：
library(psych)
fa.parallel(Harman23.cor$cov, n.obs=302, fa="pc", n.iter=100,
show.legend=FALSE, main="Scree plot with parallel analysis")
结果分析：与第一个例子类似，图形中的Kaiser-Harris准则（虚线）、碎石检验（直线与x符号）和平行分析（y=1）都建议选择两个主成分。但是三个准备并不总是相同，你可能需要依据实际情况提取不同数目的主成分，选择最优解决方案。
身体测量指标的主成分分析：
library(psych)
pc <- principal(Harman23.cor$cov, nfactors=2, rotate="none")
pc
结果分析：从代码清单中的PC1和PC2栏可以看到，第一主成分解释了身体测量指标 58%的方差，而第二主成分解释了22%，(Proportion Var行表示的是每个主成分对整个数据集的解释程度)。两者总共解释了81%的方差。对于高度变量，两者则共解释了其88%的方差。荷阵解释了成分和因子的含义。第一主成分与每个身体测量指标都正相关，看起来似乎是一个一般性的衡量因子；第二主成分与前四个变量（height、arm.span、forearm和lower.leg）负相关，与后四个变量（weight、bitro.diameter、chest.girth和chest.width）正相关，因此它看起来似乎是一个长度－容量因子。
2.3 主成分旋转
旋转是一系列将成分载荷阵变得更容易解释的数学方法，它们尽可能地对成分去噪。旋转方
法有两种：使选择的成分保持不相关（正交旋转），和让它们变得相关（斜交旋转）。旋转方法也会依据去噪定义的不同而不同。最流行的正交旋转是方差极大旋转，它试图对载荷阵的列进行去噪，使得每个成分只由一组有限的变量来解释（即载荷阵每列只有少数几个很大的载荷，其他都是很小的载荷）。对身体测量数据使用方差极大旋转，你可以得到如下所示的结果。
方差极大旋转的主成分分析：
rc <- principal(Harman23.cor$cov, nfactors=2, rotate="varimax")
rc
结果分析：列的名字都从PC变成了RC，以表示成分被旋转。观察RC1栏的载荷，你可以发现第一主成分主要由前四个变量来解释（长度变量）。RC2栏的载荷表示第二主成分主要由变量5到变量8来解释（容量变量）。注意两个主成分仍不相关，对变量的解释性不变，这是因为变量的群组没有发生变化。另外，两个主成分旋转后的累积方差解释性没有变化（81%），变的只是各个主成分对方差的解释度（成分1从58%变为44%，成分2从22%变为37%）。各成分的方差解释度趋同，准确来说，此时应该称它们为成分而不是主成分（因为单个主成分方差最大化性质没有保留）。
我们的最终目标是用一组较少的变量替换一组较多的相关变量，因此，你还需要获取每个观测在成分上的得分。
2.4 获取主成分得分
在美国法官评分例子中，我们根据原始数据中的11个评分变量提取了一个主成分。利用principal()函数，你很容易获得每个调查对象在该主成分上的得分：
library(psych)
pc <- principal(USJudgeRatings[,-1], nfactors=1, score=TRUE)
head(pc$scores)
结果分析：当scores = TRUE时，主成分得分存储在principal()函数返回对象的scores元素中。
如果有需要，你还可以获得律师与法官的接触频数与法官评分间的相关系数：
cor(USJudgeRatings$CONT, pc$score)
结果分析：显然，律师与法官的熟稔度与律师的评分毫无关联。
在身体测量数据中，你有各个身体测量指标间的相关系数，但是没有305个女孩的个体测量值。按照如下代码，你也可得到得分系数：
library(psych)
rc <- principal(Harman23.cor$cov, nfactors=2, rotate="varimax")
round(unclass(rc$weights), 2)
利用如下公式可得到主成分得分：
PC1 = 0.28*height + 0.30*arm.span + 0.30*forearm + 0.29*lower.leg -
0.06*weight - 0.08*bitro.diameter - 0.10*chest.girth - 0.04*chest.width
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）。注意，体重在PC1上的系数约为0.3或0，对于PC2也是一样。从实际角度考虑，你可以进一步简化方法，将第一主成分看作前四个变量标准化得分的均值；类似地，将第二主成分看作后四个变量标准化得分的均值，这正是通常在实际中采用的方法。