PCA:princomp() 是如何工作的,我可以用它来获取 ARIMA 的变量吗?

PCA: How does princomp() work and can I use it to pick up variables for ARIMA?

我正在尝试使用 PCA 来选择好的预测变量,以在 arima 模型的 xreg 参数中使用,以尝试预测下面的 tVar 变量。我只是使用下面的简化数据集和几个变量来简化示例。

我想了解 princomp 中的公式参数是如何工作的。对于下面的 pc 对象,它是说 "use xVar1 and xVar2 to explain the variance in 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"))

(这个问题很好post!今天还有一个关于PCA的问题post很有趣。虽然这个问题比较基础,但是关于, but the mathematical details with R code I make in 可能对任何学习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)  ## 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向量)单调递减,即对总协方差的贡献正在递减,因此我们只能select第一个k 值。因此,我们可以 select 变换矩阵的前 kX1。让我们做:

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),你可以看到一个条形图显示这一点。如果您的输入数据高度相关,那么您应该会在该图中看到接近指数的衰减,只有少数变量主导协方差。 (不幸的是,您的玩具数据无法正常工作。xVar1xVar2 是二元的,并且它们已经是线性独立的,因此在 PCA 之后,您会看到它们都相等贡献.)
  2. pc$loadings 给出 eigen_vector_mat;
  3. pc$scores 给出 X1.

使用arima()

变量select离子过程简单。如果您决定通过检查 plot(pc) 从总共 p 个变量中取出前 k 个变量,那么您将提取 pc$scores 的前 k 列] 矩阵。每列形成 z1, z2, ..., zk,并通过参数 reg.

将它们传递给 arima()

回到你关于公式的问题

For the pc object below, is it saying "use xVar1 and xVar2 to explain the variance in na.omit(dfData[,c("tVar","xVar1","xVar2")])"

经过我的解释,你应该知道答案是"No"。不要将回归步骤中使用的响应变量 tVar 与 PCA 步骤中使用的预测变量 xVar1xVars、...混合使用。

princomp()允许三种方式传入参数:

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

您选择了第一种方式。公式是用来告诉princomp()data中提取数据,然后计算模型矩阵,协方差矩阵,相关矩阵,特征分解,最终得到PCA的结果。


跟进您的评论

So if I understand correctly, PCA is primarily for reducing the number of variables, and I shouldn't include the response variable tVar in the formula or data. But I was wondering why princomp(~xVar1+xVar2, data = na.omit(dfData[,c("tVar","xVar1","xVar2")]), cor=TRUE) and princomp(na.omit(dfData[,c("xVar1","xVar2")]), cor=TRUE) are basically equivalent?

公式告诉我们如何从数据框中提取矩阵。由于您使用相同的公式 ~ xVar1 + xVar2,是否在数据框中包含 tVars 以传递给 princomp 没有区别,因为该列不会被 princomp.[=92= 触及]

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

To be clear, the strategy with PCA isn't to create a new variable which is a combination of xVar1 and xVar2 and explains most the variance in tVar, but rather to create a new variable which is a combination of xVar1 and xVar2 and explains most the variance of dfData[,c("xVar1","xVar2")]?

是的。回归(或您设置中的 arima())用于设置您的响应 tVars 和预测变量 x1, x2, ..., xpz1, z2, ..., zk 之间的关系。 regression/arima 模型将根据预测变量解释响应的均值和方差。

PCA 是一个不同的问题。它只是 select 原始预测变量 xVar1, xVar2, ... 的低等级(较少参数)表示,因此您可以在以后的回归/ARIMA 建模中使用较少的变量。

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

  1. 你有很多变数吗,比如10+?在统计建模中,达到数十万个参数是很常见的。如果我们使用所有这些,计算会变得非常慢。 PCA 在这种情况下很有用,可以降低计算复杂性,同时给出原始协方差的合理表示。
  2. 你的变量是否高度相关?如果它们很容易彼此线性独立,则 PCA 可能不会丢弃任何东西。比如你给的玩具数据xVar1xVar2只是线性无关的,所以降维是不可能的。您可以通过 pairs(mydata) 查看数据中的相关性。更好的可视化可能是使用 corrplot R 包。有关如何使用它绘制协方差矩阵的示例,请参阅