我一直在比较 Python 和 R 的几种 PCA 实现的性能,并注意到一个有趣的行为:
虽然在 Python 中计算稀疏矩阵的 PCA 似乎是不可能的(唯一的方法是 scikit-learn's TruncatedSVD,但它不支持相当于 PCA 协方差解所需的均值中心化。
他们的论点是,这会破坏矩阵的稀疏性。由于类似的原因,Facebook 的 PCA 算法或 scikit learn 中的 PCA/randomPCA 方法等其他实现不支持稀疏矩阵。
虽然所有这些对我来说都很有意义,但有几个 R 包,例如 irlba 、rsvd 等能够处理稀疏矩阵(例如使用 rsparsematrix 生成的矩阵),甚至允许特定的 center=True 参数。
我的问题是,R 如何在内部处理这个问题,因为它似乎比类似的 Python 实现要高效得多。 R 是否仍然通过进行绝对缩放来保持稀疏性(理论上这会伪造结果,但至少保持稀疏性)? 或者是否有任何方法可以显式存储零值的平均值,并且仅存储一次(而不是单独存储每个值)?
要取消保留: R 如何在内部存储以均值为中心的矩阵而不增加 RAM 使用量。 希望足够简洁......
最佳答案
这里的关键是部分 SVD ( restarted Lanczos bidiagonalization C code ) 的底层实现不存储矩阵。相反,您记录矩阵线性运算的结果 应用于从先前迭代获得的一小组向量。
我不会解释 C 代码中使用的具体方法,这是相当先进的(参见 paper 的描述),我将用一个更简单的算法来解释它,该算法捕获了如何保持效率的关键思想来自稀疏性:幂方法(或泛化为多个特征值的子空间迭代方法)。该算法通过迭代应用线性算子返回矩阵 A 的最大特征值,然后归一化(或在子空间迭代的情况下正交化一小组向量)
你在每次迭代中所做的事情是
v=A*v
v=v/norm(v)
矩阵乘法步骤是关键的一步,所以让我们看看当我们对中心 A 尝试相同的操作时会发生什么。中心 A 的矩阵公式(以 center
作为均值向量列值和 ones
作为 1 的向量)是:
A_center=A-ones*transpose(center)
因此,如果我们将迭代算法应用于这个新矩阵,我们将得到
v=A*v-dotproduct(center,v)*ones
由于 A 是稀疏的,我们可以在 (A,v) 上使用稀疏矩阵向量积,并且 -dotproduct(center,v)*ones
只需要从中减去 center 和 v 的点积所得向量在 A
维度上呈线性。
关于python - R 稀疏矩阵的内部处理,我们在Stack Overflow上找到一个类似的问题: https://stackoverflow.com/questions/50853136/