Numpy/Scipy 稀疏矩阵与稠密矩阵相乘

Question

似乎scipy稀疏矩阵类型和普通的numpy矩阵类型存在一些差异

import scipy.sparse as sp
A = sp.dia_matrix(tri(3,4))
vec = array([1,2,3,4])

print A * vec                        #array([ 1.,  3.,  6.])

print A * (mat(vec).T)               #matrix([[ 1.],
                                     #        [ 3.],
                                     #        [ 6.]])

print A.todense() * vec              #ValueError: matrices are not aligned

print A.todense() * (mat(vec).T)     #matrix([[ 1.],
                                     #        [ 3.],
                                     #        [ 6.]])

为什么稀疏矩阵可以将数组解释为列向量，而普通矩阵却无法做到呢？

- John_C

---

稀疏矩阵不是numpy矩阵的子类，甚至都不是ndarray。 - hpaulj

Answer 1

在spmatrix类中（可以在scipy/sparse/base.py中查看），__mul__()函数中有一组“ifs”可以回答你的问题：

class spmatrix(object):
    ...
    def __mul__(self, other):
        ...
        M,N = self.shape
        if other.__class__ is np.ndarray:
            # Fast path for the most common case
            if other.shape == (N,):
                return self._mul_vector(other)
            elif other.shape == (N, 1):
                return self._mul_vector(other.ravel()).reshape(M, 1)
            elif other.ndim == 2  and other.shape[0] == N:
                return self._mul_multivector(other)

对于一维数组，它将始终从compressed.py中的_cs_matrix类中的_mul_vector()进行处理。以下是代码：

def _mul_vector(self, other):
    M,N = self.shape

    # output array
    result = np.zeros(M, dtype=upcast_char(self.dtype.char,
                                           other.dtype.char))

    # csr_matvec or csc_matvec
    fn = getattr(sparsetools,self.format + '_matvec')
    fn(M, N, self.indptr, self.indices, self.data, other, result)

    return result

请注意，该函数假定稀疏矩阵的输出行数已知。基本上，它将您的输入一维数组视为适合于稀疏矩阵的列数（没有转置或非转置）。但对于ndim==2的ndarray，它无法做出这样的假设，因此如果您尝试：

vec = np.array([[1,2,3,4],
                [1,2,3,4]])

A * vec.T 是唯一可行的选项。

对于一个1维矩阵，稀疏模块也不假设它适合列数。你可以尝试检查一下:

A * mat(vec)
#ValueError: dimension mismatch

而 A * mat(vec).T 将是你唯一的选择。

- Saullo G. P. Castro

---

因为有“快速路径”的存在，所以A*vec比A*mvec更快，其中mvec=mat(vec).T。 - hpaulj