我试图使用 PCA 来选择好的预测变量，用于 arima 模型的 xreg 参数，以便预测下面的 tVar 变量。我只是使用了下面的简化数据集，只包含几个变量，以便简化示例。

我正在尝试理解 princomp 中的公式参数是如何工作的。对于下面的 pc 对象，它是否是在说“使用 xVar1 和 xVar2 来解释 na.omit(dfData[,c("tVar","xVar1","xVar2")]) 中的方差”？

我最终希望做的是创建一个新变量，该变量解释了 tVar 中的大部分方差。这是我可以使用 PCA 做的事情吗？如果可以的话，有人能解释一下如何做或者指向一个示例吗？

代码：

pc <- princomp(~xVar1+xVar2,               data = na.omit(dfData[,c("tVar","xVar1","xVar2")]),                cor=TRUE)

数据：

dput(na.omit(dfData[1:100,c("tVar","xVar1","xVar2")]))structure(list(tVar = c(11, 14, 17, 5, 5, 5.5, 8, 5.5,           6.5, 8.5, 4, 5, 9, 10, 11, 7, 6, 7, 7, 5, 6, 9, 9, 6.5, 9, 3.5,           2, 15, 2.5, 17, 5, 5.5, 7, 6, 3.5, 6, 9.5, 5, 7, 4, 5, 4, 9.5,           3.5, 5, 4, 4, 9, 4.5, 6, 10, 9.5, 15, 9, 5.5, 7.5, 12, 17.5,           19, 7, 14, 17, 3.5, 6, 15, 11, 10.5, 11, 13, 9.5, 9, 7, 4, 6,           15, 5, 18, 5, 6, 19, 19, 6, 7, 7.5, 7.5, 7, 6.5, 9, 10, 5.5,           5, 7.5, 5, 4, 10, 7, 5, 12), xVar1 = c(0L, 0L, 0L, 0L,           0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L,           0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 1L,           0L, 0L, 0L, 0L, 0L, 0L, 0L, 1L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L,           0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 1L, 0L, 0L, 0L, 0L, 0L, 0L,           1L, 0L, 0L, 0L, 1L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L,           0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L),          xVar2  = c(0L,           1L, 0L, 1L, 0L, 0L, 0L, 0L, 0L, 1L, 1L, 1L, 1L, 0L, 1L, 1L, 0L,           2L, 3L, 0L, 0L, 1L, 0L, 0L, 1L, 1L, 0L, 1L, 0L, 0L, 0L, 1L, 0L,           0L, 1L, 0L, 1L, 1L, 0L, 0L, 0L, 0L, 1L, 1L, 1L, 0L, 1L, 0L, 0L,           0L, 0L, 0L, 0L, 1L, 0L, 1L, 1L, 0L, 0L, 0L, 3L, 1L, 0L, 1L, 2L,          0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 1L, 1L, 1L, 0L,           1L, 1L, 0L, 0L, 1L, 0L, 1L, 1L, 0L, 0L, 1L, 0L, 0L, 0L, 0L, 1L,           0L)), .Names = c("tVar", "xVar1", "xVar2"          ), row.names = c(1L, 2L, 3L, 4L, 5L, 6L, 7L, 9L, 10L, 11L, 12L,           13L, 14L, 15L, 16L, 17L, 18L, 19L, 20L, 21L, 22L, 23L, 24L,25L,           26L, 27L, 28L, 29L, 30L, 31L, 32L, 33L, 34L, 35L, 36L,37L,           39L, 40L, 41L, 42L, 43L, 44L, 45L, 46L, 47L, 48L, 49L, 50L,51L,           52L, 54L, 55L, 56L, 57L, 58L, 59L, 60L, 61L, 62L, 63L, 64L, 65L,          66L, 67L, 68L, 69L, 70L, 71L, 72L, 73L, 74L, 75L, 76L, 77L, 78L,           79L, 80L, 81L, 82L, 83L, 84L, 85L, 86L, 87L, 88L, 89L, 90L, 91L,           92L, 93L, 94L, 95L, 96L, 97L, 98L, 99L, 100L),          class  = "data.frame", na.action = structure(c(8L,53L),          .Names = c("8", "53"), class = "omit"))

回答：

(这是一个非常好的帖子！今天有另一个关于 PCA 的帖子也很有趣。虽然那个问题更基础，涉及到 princomp 和 prcomp 之间的区别，但我在 答案 中用 R 代码解释的数学细节可能对任何学习 PCA 的人都有帮助。)

PCA 用于降维（低秩近似），当：

1. 你有很多（假设是 p 个）相关的变量 x1, x2, ..., xp；
2. 你想将它们缩减为少量（假设是 k < p 个）新的、线性独立的变量 z1, z2, ..., zk；
3. 你想使用 z1, z2, ..., zk 而不是 x1, x2, ..., xp 来预测响应变量 y。

基本图景和一点数学

假设你有一个响应变量 y，一个完整的线性回归模型不应丢弃任何变量，其公式应为：

y ~ x1 + x2 + ... + xp

然而，我们可以在 PCA 之后做一个合理的近似模型。设 X 为上述模型矩阵，即通过列组合所有 x1, x2, ... , xp 的观测值的矩阵，则

S <- cor(X)  ## 获取相关系数矩阵SE <- eigen(S)  ## 计算 S 的特征分解root_eigen_value <- sqrt(E$values)  ## 特征值的平方根eigen_vector_mat <- E$vectors  ## 特征向量矩阵X1 <- scale(X) %*% eigen_vector_mat  ## 转换原始矩阵

现在，root_eigen_value（一个长度为 p 的向量）是单调递减的，即对总协方差的贡献是递减的，因此我们可以只选择前 k 个值。相应地，我们可以选择转换后的矩阵 X1 的前 k 列。让我们这样做：

Z <- X1[, 1:k]

现在，我们已经成功地将 p 个变量减少到 k 个变量，每列 Z 就是新的变量 z1, z2, ..., zk。请记住，这些变量不是原始变量的子集；它们是全新的，没有名称。但由于我们只对预测 y 感兴趣，给 z1, z2, ..., zk 取什么名字并不重要。然后我们可以拟合一个近似的线性模型：

y ~ z1 + z2 + ... + zk

使用 princomp()

实际上，事情更简单，因为 princomp() 为我们完成了所有计算。通过调用：

pc <- princomp(~ x1 + x2 + ... + xp, data, cor = TRUE)

我们可以得到我们想要的所有东西。在 pc 中返回的几个值中：

1. pc$sdev 给出 root_eigen_value。如果你执行 plot(pc)，你可以看到一个显示此内容的条形图。如果你的输入数据高度相关，那么你期望看到这个图中近乎指数衰减，只有少数变量主导协方差。（不幸的是，你的玩具数据不会起作用。xVar1 和 xVar2 是二元的，它们已经是线性独立的，因此在 PCA 之后，你会看到它们都给出了相同的贡献。）
2. pc$loadings 给出 eigen_vector_mat；
3. pc$scores 给出 X1。

使用 arima()

变量选择过程很简单。如果你决定通过检查 plot(pc) 从总共 p 个变量中取出前 k 个变量，那么你提取 pc$scores 矩阵的前 k 列。每列形成 z1, z2, ..., zk，并通过参数 reg 传递给 arima()。

回到你关于公式的问题

对于下面的 pc 对象，它是否是在说“使用 xVar1 和 xVar2 来解释 na.omit(dfData[,c(“tVar”,”xVar1″,”xVar2″)]) 中的方差”

经过我的解释，你应该知道答案是“不是”。不要将回归步骤中使用的响应变量 tVar 与 PCA 步骤中使用的预测变量 xVar1、xVars 等混淆。

princomp() 允许三种方式传递参数：

1. 通过公式和数据；
2. 通过模型矩阵；
3. 通过协方差矩阵。

你选择了第一种方式。公式用于告诉 princomp() 从 data 中提取数据，它将在之后计算模型矩阵、协方差矩阵、相关系数矩阵、特征分解，直到我们最终得到 PCA 的结果。

关于你的评论的跟进

所以如果我理解正确的话，PCA 主要用于减少变量数量，我不应该在公式或数据中包含响应变量 tVar。但我想知道为什么 princomp(~xVar1+xVar2, data = na.omit(dfData[,c("tVar","xVar1","xVar2")]), cor=TRUE) 和 princomp(na.omit(dfData[,c("xVar1","xVar2")]), cor=TRUE) 基本上是等价的？

公式告诉如何从数据框中提取矩阵。由于你使用了相同的公式 ~ xVar1 + xVar2，无论你是否在传递给 princomp 的数据框中包含 tVars 都没有区别，因为该列不会被 princomp 触及。

不要在你的 PCA 公式中包含 tVars。正如我所说，回归和 PCA 是不同的问题，不应混淆。

为了清楚起见，PCA 的策略不是创建一个新的变量，它是 xVar1 和 xVar2 的组合，并且解释了 tVar 中的大部分方差，而是创建一个新的变量，它是 xVar1 和 xVar2 的组合，并且解释了 dfData[,c("xVar1","xVar2")] 中的大部分方差？

是的。回归（或在你的设置中是 arima()）用于建立你的响应 tVars 与预测变量 x1, x2, ..., xp 或 z1, z2, ..., zk 之间的关系。一个回归/ARIMA 模型将解释响应的均值和方差，根据预测变量来解释。

PCA 是一个不同的问题。它只选择一个低秩（较少参数）的表示来表示你的原始预测变量 xVar1, xVar2, ...，这样你可以在后续的回归/ARIMA 建模中使用更少的变量。

尽管如此，你可能还需要考虑是否应该为你的问题做 PCA。

1. 你有很多变量，比如 10 个以上吗？在统计建模中，达到数十万个参数是很常见的。如果我们使用所有这些参数，计算可能会变得非常慢。在这种情况下，PCA 很有用，可以减少计算复杂性，同时给出原始协方差的合理表示。
2. 你的变量高度相关吗？如果它们已经是线性独立的，PCA 可能不会删除任何东西。例如，你提供的玩具数据 xVar1 和 xVar2 只是线性独立的，因此无法进行降维。你可以通过 pairs(mydata) 查看数据中的相关性。更好的可视化方法可能是使用 corrplot R 包。参见 这个答案，了解如何使用它来绘制协方差矩阵的示例。