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

Question

我有一个大的NxN密集对称矩阵,并希望特征向量对应于k个最大特征值.找到它们的最佳方法是什么(最好是使用numpy,但也许通常使用blas/atlas/lapack,如果这是唯一的方法)？通常,N比k大得多(比如N> 5000,k <10).

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

Answer 1

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

eigvals:元组(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