在 R 中使用 princomp() 和 principal () 执行 PCA 时的不同结果答案

【问题标题】：Different results when performing PCA in R with princomp() and principal ()在 R 中使用 princomp() 和 principal () 执行 PCA 时的不同结果
【发布时间】：2016-04-16 00:54:39
【问题描述】：

我尝试使用 princomp() 和 principal() 在 R 中使用数据集 USArressts 进行 PCA。但是，对于加载/旋转和分数，我得到了两个不同的结果。

首先，我将原始数据框居中并归一化，以便比较输出结果。

library(psych)

trans_func <- function(x){
  x <- (x-mean(x))/sd(x)
  return(x)
}

A <- USArrests
USArrests <- apply(USArrests, 2, trans_func)

princompPCA <- princomp(USArrests, cor = TRUE)
principalPCA <- principal(USArrests, nfactors=4 , scores=TRUE, rotate = "none",scale=TRUE)

然后我使用以下命令获得了加载和分数的结果：

princompPCA$loadings
principalPCA$loadings

您能帮我解释一下为什么会有差异吗？我们如何解释这些结果？

【问题讨论】：

在 ?principle 文档中，“与 princomp 不同，它只返回最佳 nfactors 的子集。”
进一步说，“回归权重是从相关矩阵乘以组件负载的倒数得出的。结果是组件得分是标准得分（均值 = 0，sd = 1）标准化的输入。与 princomp 的分数的比较显示了这种差异。默认情况下，princomp 不会标准化数据矩阵，组件本身也不会标准化。"

标签： r pca psych

【解决方案1】：

?principal的帮助文档最后：

“特征向量由特征值的sqrt 重新缩放，以产生在因子分析中更典型的分量载荷。”

所以principal 返回缩放后的载荷。其实principal产生了一个通过主成分法估计的因子模型。

【讨论】：

假设您的数据矩阵被命名为“d”。然后你可以计算 d 的协方差矩阵并进行特征分解。使用特征向量的第一列乘以相应特征值的 sqrt，您将得到相同的输出 'principal(d, rotate="none",covar=T)'
所以你是说它实际上是在计算旋转因子而不是 PCA 的？
使用这些主成分有什么缺点吗？
@clancy 它是一个 PCA。因子模型可以通过 PCA 方法进行估计。它将返回唯一的因子加载。但是，如果您希望因子负载中的稀疏结果（大多数负载接近于零，有利于解释），那么您可以应用旋转方法来旋转轴。 PCA 是因子模型的简化版本。如果你不关心具体的因素。使用 PCA 估计因子模型的一个缺点是，当您向模型中添加更多因子时，因子载荷的估计不会改变。
@Nip 我现在刚刚检查了这些功能。我的回答没有太大区别。但是，我的建议是，如果您想做 PCA，最好使用“princomp”。它返回与特征分解完全相同的结果。但是，如果你想估计一个因子模型，那么你最好使用'principal'，因为它提供了更多的功能。我将在下面的答案中进行调查。

【解决方案2】：

4 年后，我想为这个问题提供更准确的答案。我以虹膜数据为例。

data = iris[, 1:4]

首先，通过特征分解做PCA

eigen_res = eigen(cov(data))
l = eigen_res$values
q = eigen_res$vectors

那么最大特征值对应的特征向量就是因子载荷

q[,1]

我们可以将其视为参考或正确答案。现在我们通过不同的 r 函数检查结果。首先，通过函数'princomp'

res1 = princomp(data)
res1$loadings[,1]
# compare with 
q[,1]

没问题，这个函数实际上只是返回与 'eigen' 相同的结果。现在转到“校长”

library(psych)
res2 = principal(data, nfactors=4, rotate="none")
# the loadings of the first PC is
res2$loadings[,1]
# compare it with the results by eigendecomposition
sqrt(l[1])*q[,1] # re-scale the eigen vector by sqrt of eigen value

您可能会发现它们仍然不同。问题是'principal'函数默认对相关矩阵进行特征分解。注意：PCA 在重新调整变量时不是不变的。如果将代码修改为

res2 = principal(data, nfactors=4, rotate="none", cor="cov")
# the loadings of the first PC is
res2$loadings[,1]
# compare it with the results by eigendecomposition
sqrt(l[1])*q[,1] # re-scale the eigen vector by sqrt of eigen value

现在，您将得到与“eigen”和“princomp”相同的结果。

总结：

如果你想做PCA，你最好应用'princomp'函数。
PCA 是因子模型的特例或因子模型的简化版本。它只是等价于特征分解。
我们可以应用 PCA 来获得因子模型的近似值。它不关心具体因素，即因素模型中的 epsilons。因此，如果您更改模型中的因子数量，您将获得相同的载荷估计。它不同于最大似然估计。
如果你在估计一个因子模型，你最好使用'principal'函数，因为它提供了更多的功能，比如旋转，通过不同的方法计算分数等等。
重新调整 PCA 模型的载荷不会对结果产生太大影响。由于您仍将数据投影到相同的最佳方向，即最大化生成 PC 的变化。

【讨论】：

感谢您的回复。在计算分数时，我还试图找出这两个函数之间的区别。我最终得出以下结论：除了符号之外，这两个分数都是相同的。我在比较 psych 的分数和 princomp 的分数时看到的差异基本上在于数据或分数本身的中心化和缩放。最后，我可以安全地使用p$scores 作为其他分析的输入，对吧？我将在下面编写我用作 aswer 的代码。
@Nip 欢迎您。是的，您可以将分数视为新变量（从原始数据中提取）并对它们进行进一步的数据分析，如 t 检验、ANOVA、回归、分类和聚类分析。也可以看作是一种特征提取方法。

【解决方案3】：

ev <- eigen(R) # R is a correlation matrix of DATA
ev$vectors %*% diag(ev$values) %*% t(ev$vectors)

pc <- princomp(scale(DATA, center = F, scale = T),cor=TRUE) 
p <-principal(DATA, rotate="none")  

#eigen values
ev$values^0.5
pc$sdev
p$values^0.5

#eigen vectors - loadings
ev$vectors
pc$loadings
p$weights %*% diag(p$values^0.5)

pc$loading %*% diag(pc$sdev)
p$loadings 

#weights
ee <- diag(0,2)
for (j in 1:2) {
 for (i in 1:2) {
  ee[i,j] <- ev$vectors[i,j]/p$values[j]^0.5
 }
};ee 

#scores
s <- as.matrix(scale(DATA, center = T, scale = T)) %*% ev$vectors
scale(s)
p$scores
scale(pc$scores)

【讨论】：