用 numpy 计算 k 个最大特征值和相应特征向量的最快方法

Posted 2023-03-12

【中文标题】用 numpy 计算 k 个最大特征值和相应特征向量的最快方法

【英文标题】：Fastest way to compute k largest eigenvalues and corresponding eigenvectors with numpy

【发布时间】：2012-08-23 10:57:23

【问题描述】：

我有一个很大的 NxN 密集对称矩阵，并且想要对应于 k 个最大特征值的特征向量。找到它们的最佳方法是什么（最好使用 numpy，但如果这是唯一的方法，通常可能使用 blas/atlas/lapack）？一般来说，N 比 k 大得多（比如 N > 5000，k

如果我的起始矩阵是稀疏的，Numpy 似乎只有找到 k 个最大特征值的函数。

【参考方案1】：

在 SciPy 中，您可以使用 linalg.eigh 函数和 eigvals 参数。

eigvals : tuple (lo, hi) 最小和最大的索引（in 升序）特征值和对应的特征向量为 返回：0

在你的情况下应该设置为(N-k,N-1)。

【讨论】：

非稀疏方法对我来说是最快的方法。使用 k=2 的 Giuliano 的基准脚本，我得到 eigh 已用时间：93.704689 eigsh 已用时间：353.433379 eig 已用时间：870.060089 最后一次是 numpy.linalg.eig。这是在我的 macbook pro 上。

【参考方案2】：

实际上，稀疏例程也适用于密集的 numpy 数组，我认为它们使用了一些 一种 Krylov 子空间迭代，因此他们需要计算几个矩阵向量 产品，这意味着如果你的 k

查看文档 http://docs.scipy.org/doc/scipy/reference/tutorial/arpack.html

和下面的代码（去和朋友喝杯好咖啡直到结束）

import numpy as np
from time import clock
from scipy.linalg import eigh as largest_eigh
from scipy.sparse.linalg.eigen.arpack import eigsh as largest_eigsh

np.set_printoptions(suppress=True)
np.random.seed(0)
N=5000
k=10
X = np.random.random((N,N)) - 0.5
X = np.dot(X, X.T) #create a symmetric matrix

# Benchmark the dense routine
start = clock()
evals_large, evecs_large = largest_eigh(X, eigvals=(N-k,N-1))
elapsed = (clock() - start)
print "eigh elapsed time: ", elapsed

# Benchmark the sparse routine
start = clock()
evals_large_sparse, evecs_large_sparse = largest_eigsh(X, k, which='LM')
elapsed = (clock() - start)
print "eigsh elapsed time: ", elapsed

【讨论】：

通过几个我尝试过的“小”k 示例，我得到 44 秒和 18 秒（eigsh 更快），当k=2 时它们大致相同，当@987654326 @（奇怪）或k 是“大”eigsh 相当慢，在所有情况下eigh 大约需要 44 秒。必须有更有效的算法来做到这一点，你会期望找到最大的特征值对，其时间比许多/全部少几个数量级......

这就是为什么我说'可以'......我真的不知道！ AFAIK 用于确定主要特征值的大多数方法都是旧的幂迭代方法的演变，这意味着它们需要多次执行 A*x，并且在 N=5000 和 A 满的情况下，这可能并不理想。无论如何，OP 都在询问 numpy/scipy 中可用的内容，结果证明他可以在 2 种方法之间进行选择。

我认为答案是：在 OP 的情况下，稀疏例程 is 更快！ :)

没有说这么小的边距，其他考虑可能会起作用，例如机器上可用的特定 BLAS/LINPACK 实现。就我而言，在 k=10 的情况下，我在 cpu 时间之间的比率与您相同（在 2008 MacBookPro 上），但速度慢了四倍。但是如果 OP 有一个超快的 BLAS/LINPACK，他甚至可能不会注意到其中的区别......

在我的笔记本电脑上，结果证明非稀疏方法最快。