在 python 中滚动 window PCA

Rolling window PCA in python

我想知道是否有人知道如何实现 rolling/moving window PCA 以在添加和删除测量时重用计算出的 PCA。

我的想法是我在很长一段时间内拥有大量数据(测量),我想从我开始时开始移动 window(比如 200 天)数据集和每一步,我都包括了第二天的测量值并丢弃了最后的测量值,所以我的 window 总是 200 天。但是,我不想每次都简单地重新计算 PCA。

是否可以制作一种比简单地为每个 window 独立计算 PCA 更有效的算法?提前致谢!

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

首先,直接回答你的问题,是的,可以制作一个比简单地为每个 window 独立计算 PCA 更有效的算法。

改进 Naive PCA 算法(低维输入)

作为该问题的第一步,我们假设您正在执行没有归一化的朴素 PCA 计算(即,您将单独处理数据,计算协方差矩阵,并找到该矩阵的 eigenvalues/eigenvectors).

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

"""
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)

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

虽然协方差矩阵不是最终答案,我们仍然需要找到主成分。如果您不是自己手动运行特征求解器,则无需做很多事情——您每次都需要为新的特征值和特征向量付出代价。

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

改进无协方差算法(高维输入)

通常(也许总是?),无协方差的 PCA 算法有某种迭代求解器(很像特征求解器),但它们有计算捷径,允许在不显式具体化协方差矩阵的情况下找到 eigenvalues/eigenvectors。

任何此类方法都可能有额外的技巧,允许您将一些信息从一个 window 保存到下一个,但通常人们会期望您可以通过重新使用现有的主成分而不是使用随机输入来启动求解器(很像上面的特征求解器案例)。

Window 使用朴素算法的归一化

假设您将每个 window 归一化为每列的平均值为 0(在 PCA 中很常见),那么在修改协方差矩阵时您将需要做一些额外的工作。

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

"""
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

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

根据您使用的算法,此处可用的技巧会有很大差异,但我首先尝试的一般策略是使用滚动算法来跟踪均值和标准当前 window 的每一列的偏差并修改迭代求解器以将其考虑在内(即,给定一个 window X 你想在重新缩放的 window Y 上迭代 - 替换Y=a*X+b 到您选择的迭代算法中并进行符号化简化,以希望产生一个具有少量额外常数成本的版本。

和以前一样,您需要重新使用找到的任何主成分,而不是对每个 window.

使用随机初始化向量