使用广播与稀疏的scipy矩阵。

Question

我有一个形状为(k,N)的numpy数组Z和第二个形状为(N,n)的数组X。 利用numpy广播，我可以轻松地获得一个新的数组H，形状为(n,k,N)，其切片是通过将Z数组的行乘以X的列而得到的。

H = Z.reshape((1, k, N)) * X.T.reshape((n, 1, N))

这个运行得很好，而且非常快。 现在，X 非常稀疏，我想使用稀疏矩阵操作进一步加速此操作。

但是如果我执行以下操作:

import scipy.sparse as sprs
spX = sprs.csr_matrix(X)
H = (Z.reshape((1,k,N))*spX.T.reshape((n,1,N))).dot(Z.T)

I get the following error:

Traceback (most recent call last):
  File "<input>", line 1, in <module>
  File "C:\Python27\lib\site-packages\scipy\sparse\base.py", line 126, in reshape
    self.__class__.__name__)
NotImplementedError: Reshaping not implemented for csc_matrix.

有没有一种方法可以在使用稀疏的scipy矩阵时进行广播？

- Ulderique Demoitre

---

除了二维限制外，稀疏数学只能与其他稀疏矩阵一起使用。在使用密集数组时，稀疏矩阵会被转换为密集矩阵。尝试使用“spX”进行各种数学运算，看看会发生什么。 - hpaulj

Answer 1

Scipy稀疏矩阵仅限于2D形状。但是您可以以“稀疏”的方式使用Numpy：

H = np.zeros((n,k,N), np.result_type(Z, X))
I, J = np.nonzero(X)
Z_ = np.broadcast_to(Z, H.shape)
H[J,:,I] = Z_[J,:,I] * X[I,J,None]

不幸的是，结果H仍然是一个密集数组。

N.b. 使用None进行索引是在所需轴上添加单位长度维度的方便方法。当将高级索引与切片组合时，结果的顺序在文档中有解释。

- user7138814