如何通过广播密集1d数组元素乘法scipy.sparse矩阵？

假设我有一个2d稀疏数组。在我的实际用例的行和列两者的数量要大得多（说20000和50000），因此它不适合在内存中，当使用密表示：如何通过广播密集1d数组元素乘法scipy.sparse矩阵？

>>> import numpy as np 
>>> import scipy.sparse as ssp 

>>> a = ssp.lil_matrix((5, 3)) 
>>> a[1, 2] = -1 
>>> a[4, 1] = 2 
>>> a.todense() 
matrix([[ 0., 0., 0.], 
     [ 0., 0., -1.], 
     [ 0., 0., 0.], 
     [ 0., 0., 0.], 
     [ 0., 2., 0.]])

现在假设我有所有的密集维数组与（在我的现实生活中的情况下，或50000）大小3个非零组件：

>>> d = np.ones(3) * 3 
>>> d 
array([ 3., 3., 3.])

我想计算的的元素单元的乘法和使用numpy的通常的广播语义天。然而，在SciPy的稀疏矩阵是np.matrix的：在“*”操作符被重载把它像一个矩阵乘法，而不是按元素乘法：

>>> a * d 
array([ 0., -3., 0., 0., 6.])

一个解决办法是做“一个”切换到阵列的语义为‘*’操作，这将产生预期的结果：

>>> a.toarray() * d 
array([[ 0., 0., 0.], 
     [ 0., 0., -3.], 
     [ 0., 0., 0.], 
     [ 0., 0., 0.], 
     [ 0., 6., 0.]])

但我不能这样做，因为调用ToArray的（）将兑现的密集版本的‘一个’，这不适合内存（并且结果也会很密集）：

>>> ssp.issparse(a.toarray()) 
False

任何想法如何建立这个，而只保留稀疏的数据结构，而不必在'a'的列上做一个不够高效的python循环？

来源

2010-07-14 ogrisel

如果'D'是相同尺寸的稀疏矩阵为'A'可以使用'a.multiply（d）'。也许你可以做一个长度为N行的“d”，并且每次循环N行“a”？ – mtrw 2010-07-14 19:43:44

但是d密集并且不能在存储器中明确地被广播以满足多重形状要求。循环一个批处理是一种选择，但我觉得这有点hackish。我原以为有没有python循环这样做的香草vectorized/scipy的方式来做到这一点。 – ogrisel 2010-07-14 22:15:41

我想问题是你想要一个（稀疏）矩阵的表示，但是一个数组的多重操作。我认为你不得不不幸地推出自己的产品。 – mtrw 2010-07-15 19:33:00

我在scipy.org上也回复了，但我想我应该在这里添加一个答案，以防其他人在搜索时找到这个页面。

您可以将矢量转换为稀疏对角矩阵，然后使用矩阵乘法（使用*）与广播做同样的事情，但效率很高。

>>> d = ssp.lil_matrix((3,3)) 
>>> d.setdiag(np.ones(3)*3) 
>>> a*d 
<5x3 sparse matrix of type '<type 'numpy.float64'>' 
with 2 stored elements in Compressed Sparse Row format> 
>>> (a*d).todense() 
matrix([[ 0., 0., 0.], 
     [ 0., 0., -3.], 
     [ 0., 0., 0.], 
     [ 0., 0., 0.], 
     [ 0., 6., 0.]])

希望有帮助！

来源

2010-07-28 21:27:58 mitmatt

感谢它看起来会解决我的问题。 – ogrisel 2010-07-30 09:05:43

关于这一点的好处在于，当'X'是一个'ndarray'或一个密集矩阵时它也可以工作。 +1。 – 2012-09-08 16:40:41

这可以使用['scipy.sparse.diags（d，0）']进一步简化（http://docs.scipy.org/doc/scipy-0.16.1/reference/generated/scipy.sparse.diags。 html）而不是'lil_matrix' – 2016-01-21 23:05:31

那么，这里有一个简单的代码，将做你想做的。我不知道这是否是有效率，你想，所以要么接受，要么离开它：

import scipy.sparse as ssp 
def pointmult(a,b): 
    x = a.copy() 
    for i in xrange(a.shape[0]): 
     if x.data[i]: 
      for j in xrange(len(x.data[i])): 
       x.data[i] *= b[x.rows[i]] 
    return x

它只能与律矩阵的工作，所以你必须做出一些改变，如果你想要的工作与其他格式。

来源

2010-07-18 22:34:09

谢谢我希望避免python中的循环。但是，对于这个用例来说，当前的scipy.sparse类可能没有出路。 – ogrisel 2010-07-19 17:49:07

我认为A.multiply（B）应该在scipy稀疏工作。该方法乘法执行“逐点”乘法，而不是矩阵乘法。

HTH

来源

2010-12-21 20:01:07 salma

结果是一个密集的矩阵。不好。 – 2016-07-15 11:45:58

@ K3 --- rnc只有B密集时结果才是密集的。如果您将B转换为任何稀疏格式，它将会执行此操作。例如。 A.multiply（csc_matrix（B）） – markhor 2017-01-23 18:27:52

如何通过广播密集1d数组元素乘法scipy.sparse矩阵？

回答

相关问题