r - PCA:princomp()如何工作，我可以用它为ARIMA提取变量吗？

Question

我正在尝试使用PCA选择好的预测变量，以在xreg模型的arima参数中使用，以尝试预测下面的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,38L, 
          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"))

Answer 1

（这是一篇非常好的文章！今天有关于PCA的另一篇文章很有趣。尽管这个问题更基本，但关于the difference between princomp and prcomp，但是我在the answer中使用R代码编写的数学细节可能对任何人都有益学习PCA。）

在以下情况下，将PCA用于缩小尺寸（低阶近似）：


您有很多（例如p）相关变量x1, x2, ..., xp；
您希望将它们缩小为少量（例如k < p）新的线性独立变量z1, z2, ..., zk；
您想使用z1, z2, ..., zk而不是x1, x2, ..., xp来预测响应变量y。




基本情况和一些数学

假设您有一个响应变量y，不删除任何变量的全线性回归应采用以下公式：

y ~ x1 + x2 + ... + xp

但是，在PCA之后，我们可以做一个合理的近似模型。令 X为上面的模型矩阵，即通过按列组合所有 x1, x2, ... , xp观察值的矩阵，然后

S <- cor(X)  ## get correlation matrix S
E <- eigen(S)  ## compute eigen decomposition of S
root_eigen_value <- sqrt(E$values)  ## square root of eigen values
eigen_vector_mat <- E$vectors  ## matrix of eigen vectors
X1 <- scale(X) %*% eigen_vector_mat  ## transform original matrix

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

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的一些返回值中：


pc$sdev给出 root_eigen_value。如果执行 plot(pc)，则可以看到显示此内容的条形图。如果您的输入数据高度相关，那么您期望在该图中看到接近指数的衰减，只有少数变量支配协方差。 （不幸的是，您的玩具数据不起作用。 xVar1和 xVar2是二进制的，并且它们已经线性独立，因此在PCA之后，您将看到它们都发挥了相同的作用。）
pc$loadings给出 eigen_vector_mat;
pc$scores给出 X1。




使用 arima()

变量选择过程很简单。如果决定通过检查 k从总共 p个变量中取出第一个 plot(pc)变量，则可以提取 k矩阵的第一个 pc$scores列。每列构成 z1, z2, ..., zk，然后通过参数 arima()将它们传递给 reg。



回到关于配方的问题


  对于下面的pc对象，是说“使用xVar1和xVar2来解释na.omit（dfData [，c（“ tVar”，“ xVar1”，“ xVar2”）]）的差异”


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

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


通过公式和数据；
通过模型矩阵
通过协方差矩阵


您选择了第一种方式。该公式用于告诉 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之间的关系。回归/主题模型将根据预测变量解释响应的均值和方差。

PCA是另一个问题。它仅选择原始预测变量 xVar1, xVar2, ...的低秩（较少参数）表示形式，因此您可以在以后的回归/ ARIMA建模中使用较少的变量。

不过，您可能仍需要考虑是否应该对问题进行PCA。


您是否有很多变量，例如10+？在统计建模中，通常达到数十万个参数。如果我们全部使用它们，计算会变得很慢。在这种情况下，PCA很有用，可降低计算复杂度，同时合理地表示原始协方差。
您的变量是否高度相关？如果它们彼此线性独立，PCA可能不会掉落任何东西。例如，您提供的玩具数据 xVar1和 xVar2只是线性独立的，因此无法减小尺寸。您可以按 pairs(mydata)查看数据中的相关性。更好的可视化可能是使用 corrplot R包。有关如何使用它绘制协方差矩阵的示例，请参见 this answer。