使用scipy求解大型稀疏矩阵的逆矩阵

6
我需要反转一个大的稀疏矩阵。我不能避免矩阵反转,唯一的捷径就是只关注主对角线元素,忽略非对角线元素(虽然这不是我的首选方案,但作为解决方法也可以接受)。
我需要反转的矩阵通常很大(40000 * 40000),只有一小部分对角线上的元素为非零值。我目前的做法是建立所有的稀疏矩阵,然后使用以下代码进行反转: posterior_covar = np.linalg.inv( hessian.todense() ) 显然,这需要很长时间和大量内存。
请问有什么提示,或者只需要耐心等待或减小问题规模即可吗?

scipy文档中提到您不需要将矩阵稠密化,所以我有点困惑。http://docs.scipy.org/doc/scipy-0.8.x/reference/sparse.html - BenDundee
3
SciPy的0.12版本(目前正在进行测试)具有函数scipy.sparse.linalg.inv。 - Warren Weckesser
1个回答

8

我不认为稀疏模块有明确的逆方法,但它确实有稀疏求解器。类似这个玩具例子的东西可以工作:

>>> a = np.random.rand(3, 3)
>>> a
array([[ 0.31837307,  0.11282832,  0.70878689],
       [ 0.32481098,  0.94713997,  0.5034967 ],
       [ 0.391264  ,  0.58149983,  0.34353628]])
>>> np.linalg.inv(a)
array([[-0.29964242, -3.43275347,  5.64936743],
       [-0.78524966,  1.54400931, -0.64281108],
       [ 1.67045482,  1.29614174, -2.43525829]])

>>> a_sps = scipy.sparse.csc_matrix(a)
>>> lu_obj = scipy.sparse.linalg.splu(a_sps)
>>> lu_obj.solve(np.eye(3))
array([[-0.29964242, -0.78524966,  1.67045482],
       [-3.43275347,  1.54400931,  1.29614174],
       [ 5.64936743, -0.64281108, -2.43525829]])

请注意,结果已经转置了!
如果您希望反函数也是稀疏的,并且从最后一个求解返回的密集矩阵无法适合内存,那么您也可以逐行(列)生成它,提取非零值,并从中构建稀疏逆矩阵。
>>> for k in xrange(3) :
...     b = np.zeros((3,))
...     b[k] = 1
...     print lu_obj.solve(b)
... 
[-0.29964242 -0.78524966  1.67045482]
[-3.43275347  1.54400931  1.29614174]
[ 5.64936743 -0.64281108 -2.43525829]

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接