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

Question

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

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

Answer 1

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

eigvals : tuple (lo, hi) 要返回的最小和最大（按升序）特征值和相应特征向量的索引：0 <= lo < hi <= M-1。如果省略，则返回所有特征值和特征向量。

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

Answer 2

实际上，稀疏例程也适用于密集的 numpy 数组，我认为它们使用某种 Krylov 子空间迭代，因此他们需要计算几个矩阵向量乘积，这意味着如果您的 k << N，稀疏例程可能是（勉强？ ） 快点。

查看文档 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