python中的滚动窗口PCA

Question

我想知道是否有人知道如何实现滚动/移动窗口 PCA，它在添加和删除测量值时重用计算的 PCA。

我的想法是我在很长一段时间内拥有大量数据（测量），我希望从我的数据集的开头开始有一个移动窗口（比如 200 天），每一步，我包括第二天的测量并丢弃最后一次测量，所以我的窗口总是 200 天。但是，我不想每次都简单地重新计算 PCA。

是否有可能制定一种比单独计算每个窗口的 PCA 更有效的算法？提前致谢！

Answer 1

一个完整的答案取决于很多因素。我将介绍我认为最重要的此类因素，希望这些信息足以为您指明正确的方向。

首先，直接回答您的问题，是的，可以制定一种比单独计算每个窗口的 PCA 更有效的算法。

改进 Naive PCA 算法（低维输入）

作为该问题的第一步，我们假设您正在执行没有归一化的简单 PCA 计算（即，您将数据单独放置，计算协方差矩阵，并找到该矩阵的特征值/特征向量）。

当面对我们要计算其 PCA 的输入矩阵 X 时，朴素算法首先计算协方差矩阵 W = X.T @ X。一旦我们计算出某个包含 200 个元素的窗口，我们就可以通过移除它们对协方差的贡献来廉价地从原始数据集中添加或移除元素。

"""
W:   shape (p, p)
row: shape (1, p)
"""

def add_row(W, row):
  return W + (row.T @ row)

def remove_row(W, row):
  return W - (row.T @ row)

您对滑动窗口的描述相当于删除一行并添加一个新行，因此我们可以使用 O(p^2) 计算快速计算一个新的协方差矩阵，而不是 O(n p^2) 一个典型的矩阵乘法需要（对于这个问题，n==200）。

虽然协方差矩阵不是最终答案，但我们仍然需要找到主成分。如果您不自己手动滚动特征解算器，则没有太多工作要做——您每次都要为新的特征值和特征向量付出代价。

但是，如果您正在编写自己的 eigensolver，大多数此类方法都会接受起始输入并迭代直到某个终止条件（通常是最大迭代次数，或者如果错误变得足够低，以您先达到的为准）。换出单个数据点不太可能彻底改变主成分，因此对于典型数据，人们可能期望重新使用现有的特征值/特征向量作为特征求解器的输入将允许您在比开始时少得多的迭代中终止从随机输入，提供额外的加速。

改进无协方差算法（高维输入）

通常（也许总是？），无协方差 PCA 算法具有某种迭代求解器（很像特征求解器），但它们具有计算快捷方式，可以在不显式具体化协方差矩阵的情况下找到特征值/特征向量。

任何单独的此类方法都可能具有额外的技巧，允许您将一些信息从一个窗口保存到下一个窗口，但通常人们希望您可以通过重新使用现有的主成分来减少迭代的总数使用随机输入来启动求解器（很像上面的 eigensolver 案例）。

使用朴素算法的窗口标准化

假设您将每个窗口标准化为每列的平均值为 0（在 PCA 中很常见），那么在修改协方差矩阵时您将需要做一些额外的工作。

首先，我假设您已经有一个滚动机制来跟踪需要从一个窗口应用到下一个窗口的任何差异。如果没有，请考虑以下内容：

"""
We're lazy and don't want to handle a change in sample
size, so only work with row swaps -- good enough for
a sliding window.
old_row: shape (1, p)
new_row: shape (1, p)
"""

def replaced_row_mean_adjustment(old_row, new_row):
  return (new_row - old_row)/200.  # whatever your window size is

对协方差矩阵的影响计算起来还不错，但我还是会在这里放一些代码。

"""
W:      shape (p, p)
center: shape (1, p)
  exactly equal to the mean diff vector we referenced above
X:      shape (200, p)
  exactly equal to the window you're examining after any
  previous mean centering has been applied, but before the
  current centering has happened. Note that we only use
  its row and column sums, so you could get away with
  a rolling computation for those instead, but that's
  a little more code, and I want to leave at least some
  of the problem for you to work on
"""

def update_covariance(W, center, X):
  result = W
  result -= center.T @ np.sum(X, axis=0).reshape(1, -1)
  result -= np.sum(X, axis=1).reshape(-1, 1) @ center
  result += 200 * center.T @ center  # or whatever your window size is
  return result

将标准差重新调整为 1 在 PCA 中也很常见。这也很容易适应。

"""
Updates the covariance matrix assuming you're modifing a window
of data X with shape (200, p) by multiplying each column by
its corresponding element in v. A rolling algorithm to compute
v isn't covered here, but it shouldn't be hard to figure out.

W: shape (p, p)
v: shape (1, p)
"""

def update_covariance(W, v):
  return W * (v.T @ v)  # Note that this is element-wise multiplication of W

使用无协方差算法的窗口归一化

您在此处可用的技巧会因您使用的算法而有很大差异，但我首先尝试的一般策略是使用滚动算法来跟踪平均值和标准差当前窗口的每一列并修改迭代求解器以考虑到这一点（即，给定一个窗口 X，您想在重新缩放的窗口 Y 上迭代 - 将 Y=a*X+b 替换为您选择的迭代算法并象征性地简化以希望产生一个具有少量额外恒定成本的版本。

和以前一样，您需要重新使用找到的任何主成分，而不是为每个窗口使用随机初始化向量。