python - Scipy 稀疏矩阵特殊减法

Question

我正在做一个项目，我在其中进行了大量的矩阵计算。

我正在寻找一种聪明的方法来加速我的代码。在我的项目中，我正在处理一个大小为 100Mx1M 且具有大约 10M 非零值的稀疏矩阵。下面的例子只是为了证明我的观点。

假设我有:

• 大小为 (2) 的向量 v
• 大小为 (3) 的向量 c

• 大小为 (2,3) 的稀疏矩阵 X

  v = np.asarray([10, 20])
  c = np.asarray([ 2,  3,  4])
  data = np.array([1, 1, 1, 1])
  row  = np.array([0, 0, 1, 1])
  col  = np.array([1, 2, 0, 2])
  X = coo_matrix((data,(row,col)), shape=(2,3))
  X.todense()
  # matrix([[0, 1, 1],
  #         [1, 0, 1]])
  

目前我在做:

result = np.zeros_like(v)
d = scipy.sparse.lil_matrix((v.shape[0], v.shape[0]))
d.setdiag(v)
tmp = d * X

print tmp.todense()
#matrix([[  0.,  10.,  10.],
#        [ 20.,   0.,  20.]])
# At this point tmp is csr sparse matrix

for i in range(tmp.shape[0]):
     x_i = tmp.getrow(i)
     result += x_i.data * ( c[x_i.indices] - x_i.data)
     # I only want to do the subtraction on non-zero elements    

print result
# array([-430, -380])

我的问题是 for 循环，尤其是减法。 我想找到一种通过仅减去非零元素来向量化此操作的方法。

在减法上直接得到稀疏矩阵的东西:

matrix([[  0.,  -7.,  -6.],
       [ -18.,   0.,  -16.]])

有没有办法聪明地做到这一点？

Answer 1

您不需要遍历行来执行您已经在执行的操作。您可以使用类似的技巧将行乘以第一个向量:

import scipy.sparse as sps

# number of nonzero entries per row of X
nnz_per_row = np.diff(X.indptr)
# multiply every row by the corresponding entry of v
# You could do this in-place as:
# X.data *= np.repeat(v, nnz_per_row)
Y = sps.csr_matrix((X.data * np.repeat(v, nnz_per_row), X.indices, X.indptr),
                   shape=X.shape)

# subtract from the non-zero entries the corresponding column value in c...
Y.data -= np.take(c, Y.indices)
# ...and multiply by -1 to get the value you are after
Y.data *= -1

要查看它是否有效，请设置一些虚拟数据

rows, cols = 3, 5
v = np.random.rand(rows)
c = np.random.rand(cols)
X = sps.rand(rows, cols, density=0.5, format='csr')

运行上面的代码后:

>>> x = X.toarray()
>>> mask = x == 0
>>> x *= v[:, np.newaxis]
>>> x = c - x
>>> x[mask] = 0
>>> x
array([[ 0.79935123,  0.        ,  0.        , -0.0097763 ,  0.59901243],
       [ 0.7522559 ,  0.        ,  0.67510109,  0.        ,  0.36240006],
       [ 0.        ,  0.        ,  0.72370725,  0.        ,  0.        ]])
>>> Y.toarray()
array([[ 0.79935123,  0.        ,  0.        , -0.0097763 ,  0.59901243],
       [ 0.7522559 ,  0.        ,  0.67510109,  0.        ,  0.36240006],
       [ 0.        ,  0.        ,  0.72370725,  0.        ,  0.        ]])

您累积结果的方式要求每一行中有相同数量的非零条目，这似乎是一件很奇怪的事情。你确定那是你想要的吗？如果这真的是您想要的，您可以通过以下方式获得该值(value):

result = np.sum(Y.data.reshape(Y.shape[0], -1), axis=0)

但我很难相信这真的是您所追求的......