scipy 稀疏矩阵除法答案

【问题标题】：scipy sparse matrix divisionscipy 稀疏矩阵除法
【发布时间】：2017-07-02 16:34:40
【问题描述】：

我一直在尝试将 python scipy 稀疏矩阵除以其行的向量和。这是我的代码

sparse_mat = bsr_matrix((l_data, (l_row, l_col)), dtype=float)
sparse_mat = sparse_mat / (sparse_mat.sum(axis = 1)[:,None])

但是，无论我如何尝试，它都会引发错误

sparse_mat = sparse_mat / (sparse_mat.sum(axis = 1)[:,None])
File "/usr/lib/python2.7/dist-packages/scipy/sparse/base.py", line 381, in __div__
return self.__truediv__(other)
File "/usr/lib/python2.7/dist-packages/scipy/sparse/compressed.py", line 427, in __truediv__
raise NotImplementedError
NotImplementedError

有人知道我哪里出错了吗？

【问题讨论】：

划分调用true_division，这是一个元素划分。这似乎没有为多个值实现。因此，(sparse_mat.sum(axis = 1)[:,None] 的结果很可能不是一个数字。
@Dschoni 是的，结果是一个向量，我的目标是将稀疏矩阵每一行中的每个元素除以行元素的总和。所以如果 M=[[2,4],[1,2]]，我想得到 Ans=[[2/6, 4/6],[1/3, 2/3]]。
你试过sparse_mat = sparse_mat*(1 / (sparse_mat.sum(axis = 1)[:,None]))。稀疏矩阵的划分似乎是问题所在。您可能还需要将除数转换为密集数组sparse_mat = sparse_mat*(1 / (sparse_mat.sum(axis = 1).toarray()[:,None]))
@uchman21 请提供一个独立的小例子。此问题可能与您放入矩阵的数据有关。（或者可能是你的 scipy 太旧了——我尝试过的稀疏矩阵除法适用于 Python 3 和 scipy 0.18。）
我正在使用 python 2.7.13 和 scipy 0.18。该矩阵只是一个 232 x 232 的简单稀疏矩阵

标签： python numpy scipy sparse-matrix

【解决方案1】：

您可以通过从行和的倒数创建一个稀疏对角矩阵然后将其与矩阵相乘来规避该问题。在乘积中，对角矩阵向左，你的矩阵向右。

例子：

>>> a
array([[0, 9, 0, 0, 1, 0],
       [2, 0, 5, 0, 0, 9],
       [0, 2, 0, 0, 0, 0],
       [2, 0, 0, 0, 0, 0],
       [0, 9, 5, 3, 0, 7],
       [1, 0, 0, 8, 9, 0]])
>>> b = sparse.bsr_matrix(a)
>>> 
>>> c = sparse.diags(1/b.sum(axis=1).A.ravel())
>>> # on older scipy versions the offsets parameter (default 0)
... # is a required argument, thus
... # c = sparse.diags(1/b.sum(axis=1).A.ravel(), 0)
...
>>> a/a.sum(axis=1, keepdims=True)
array([[ 0.        ,  0.9       ,  0.        ,  0.        ,  0.1       ,  0.        ],
       [ 0.125     ,  0.        ,  0.3125    ,  0.        ,  0.        ,  0.5625    ],
       [ 0.        ,  1.        ,  0.        ,  0.        ,  0.        ,  0.        ],
       [ 1.        ,  0.        ,  0.        ,  0.        ,  0.        ,  0.        ],
       [ 0.        ,  0.375     ,  0.20833333,  0.125     ,  0.        ,  0.29166667],
       [ 0.05555556,  0.        ,  0.        ,  0.44444444,  0.5       ,  0.        ]])
>>> (c @ b).todense() # on Python < 3.5 replace c @ b with c.dot(b)
matrix([[ 0.        ,  0.9       ,  0.        ,  0.        ,  0.1       ,  0.        ],
        [ 0.125     ,  0.        ,  0.3125    ,  0.        ,  0.        ,  0.5625    ],
        [ 0.        ,  1.        ,  0.        ,  0.        ,  0.        ,  0.        ],
        [ 1.        ,  0.        ,  0.        ,  0.        ,  0.        ,  0.        ],
        [ 0.        ,  0.375     ,  0.20833333,  0.125     ,  0.        ,  0.29166667],
        [ 0.05555556,  0.        ,  0.        ,  0.44444444,  0.5       ,  0.        ]])

【讨论】：

我试过这个解决方案，但它给出了错误 elem_sum = csc_matrix((1/sparse_mat.sum(axis = -1).A.ravel(), numpy.arange(sparse_mat.shape[0] ), numpy.arange(sparse_mat.shape[0]+1))) 文件“/usr/lib/python2.7/dist-packages/scipy/sparse/compressed.py”，第 548 行，总和返回 spmatrix.sum (self,axis) 文件“/usr/lib/python2.7/dist-packages/scipy/sparse/base.py”，第 629 行，总之 raise ValueError("axis out of bounds") ValueError: axis out of bounds
@uchman21 奇怪，然后尝试axis = 1，这似乎在您的代码中有效。
是的，这行得通。但是，对我来说，我需要将其设置为 c = sparse.diags(1/b.sum(axis=1).A.ravel(),0) 以指定主对角线以使其最终起作用。请将其添加到您的答案中。
b.sum(axis=1).A1 应该可以工作。 sum 产生一个 np.matrix,which has a A1`属性。 stackoverflow.com/a/20765358/901925

【解决方案2】：

发生了一些有趣的事情。执行元素分割没有问题。我想知道这是否是 Py2 的问题。我正在使用 Py3。

In [1022]: A=sparse.bsr_matrix([[2,4],[1,2]])
In [1023]: A
Out[1023]: 
<2x2 sparse matrix of type '<class 'numpy.int32'>'
    with 4 stored elements (blocksize = 2x2) in Block Sparse Row format>
In [1024]: A.A
Out[1024]: 
array([[2, 4],
       [1, 2]], dtype=int32)
In [1025]: A.sum(axis=1)
Out[1025]: 
matrix([[6],
        [3]], dtype=int32)
In [1026]: A/A.sum(axis=1)
Out[1026]: 
matrix([[ 0.33333333,  0.66666667],
        [ 0.33333333,  0.66666667]])

或尝试其他示例：

In [1027]: b=sparse.bsr_matrix([[0, 9, 0, 0, 1, 0],
      ...:        [2, 0, 5, 0, 0, 9],
      ...:        [0, 2, 0, 0, 0, 0],
      ...:        [2, 0, 0, 0, 0, 0],
      ...:        [0, 9, 5, 3, 0, 7],
      ...:        [1, 0, 0, 8, 9, 0]])
In [1028]: b
Out[1028]: 
<6x6 sparse matrix of type '<class 'numpy.int32'>'
    with 14 stored elements (blocksize = 1x1) in Block Sparse Row format>
In [1029]: b.sum(axis=1)
Out[1029]: 
matrix([[10],
        [16],
        [ 2],
        [ 2],
        [24],
        [18]], dtype=int32)
In [1030]: b/b.sum(axis=1)
Out[1030]: 
matrix([[ 0.        ,  0.9       ,  0.        ,  0.        ,  0.1       , 0.        ],
        [ 0.125     ,  0.        ,  0.3125    ,  0.        ,  0.        , 0.5625    ],
 ....
        [ 0.05555556,  0.        ,  0.        ,  0.44444444,  0.5       ,     0.        ]])

这种稀疏/密集的结果也是密集的，其中c*b（c 是稀疏对角线）是稀疏的。

In [1039]: c*b
Out[1039]: 
<6x6 sparse matrix of type '<class 'numpy.float64'>'
    with 14 stored elements in Compressed Sparse Row format>

稀疏和是一个密集矩阵。它是二维的，因此无需扩展它的尺寸。事实上，如果我尝试我得到一个错误：

In [1031]: A/(A.sum(axis=1)[:,None])
....
ValueError: shape too large to be a matrix.

【讨论】：

这似乎取决于 scipy 版本。使用过时的版本，这实际上按我的预期工作，其中两个稀疏向量的划分返回一个稀疏向量。毕竟，如果被除数有一个空单元格，那么这个单元格中的结果无论如何都应该是 0。对于较新版本的 scipy，同一行返回一个密集的 numpy 矩阵...

【解决方案3】：

根据this message，为了保持矩阵稀疏，您可以访问 data 值并使用（非零）索引：

sums = np.asarray(A.sum(axis=1)).squeeze()  # this is dense
A.data /= sums[A.nonzero()[0]]

如果除以非零行平均值而不是总和，则可以

nnz = A.getnnz(axis=1)  # this is also dense
means = sums / nnz
A.data /= means[A.nonzero()[0]]

【讨论】：