Scipy Sparse:意外的身份行为

Question

下面是我创建一个以c为对角线的稀疏矩阵的尝试。我知道还有其他方法可以解决这个问题，但是我很好奇为什么下面的代码不能像预期的那样工作：

import numpy as np
import scipy.sparse as sparse

c = np.arange(0,5)
>>> np.identity(5)*c
array([[ 0.,  0.,  0.,  0.,  0.],
       [ 0.,  1.,  0.,  0.,  0.],
       [ 0.,  0.,  2.,  0.,  0.],
       [ 0.,  0.,  0.,  3.,  0.],
       [ 0.,  0.,  0.,  0.,  4.]])
>>> sparse.identity(5)*c
array([ 0.,  1.,  2.,  3.,  4.])
#expected output:
<5x5 sparse matrix of type '<type 'numpy.float64'>'
with 5 stored elements (1 diagonals) in DIAgonal format>
# and (sparse.identity(5)*c).todense() == np.identity(5)*c

Answer 1

在表达式sparse.identity(5)*c中，使用了稀疏矩阵的乘法运算符，它是代数矩阵乘法(即，矩阵乘以向量给出一个向量)。

可以使用scipy.sparse.diags创建具有给定对角线的稀疏对角线矩阵

In [18]: from scipy import sparse

In [19]: c = np.arange(5)

In [20]: d = sparse.diags(c, 0)

In [21]: d
Out[21]: 
<5x5 sparse matrix of type '<type 'numpy.float64'>'
    with 5 stored elements (1 diagonals) in DIAgonal format>

In [22]: d.A
Out[22]: 
array([[ 0.,  0.,  0.,  0.,  0.],
       [ 0.,  1.,  0.,  0.,  0.],
       [ 0.,  0.,  2.,  0.,  0.],
       [ 0.,  0.,  0.,  3.,  0.],
       [ 0.,  0.,  0.,  0.,  4.]])

Answer 2

有两种常见类型的数组乘法，逐个元素和矩阵。

在MATLAB中，*是矩阵形式，.*是逐个元素。

在*中，numpy是一个元素接一个元素(带有广播)，'np.dot‘是矩阵乘法的基本形式。Python开发人员已经批准@作为一个运算符，可以(最终)用于矩阵乘法。

对于matrix子类，*是矩阵乘法，np.multiply用于逐个元素。(np.multiply也适用于ndarray。)

scipy遵循np.matrix约定。*是矩阵乘法。sparse.identity(5).multiply(c)通过元素乘法进行元素乘法(尽管它返回一个np.matrix，而不是稀疏的)。

至于原因，归根结底是开发人员习惯的约定。对于线性代数问题，矩阵乘法很常见，因此在sparse中使用矩阵乘法。np.matrix复制了MATLAB的约定。MATLAB是用来访问FORTRAN矩阵库的。

在物理学中还有另一个惯例，爱因斯坦符号。这是一个广义的矩阵乘法，扩展到更多的维度。np.einsum可以做到这一点。它可以实现一个元素一个元素的乘法，尽管它的核心使用了“乘积和”方法。但是它不是为np.matrix或sparse实现的(实际上也不需要，因为它们总是2D的)。

请注意，您的示例的规范是多么相似

np.einsum('ij,j->ij',np.identity(5),c) # element by element
np.einsum('ij,j->i',np.identity(5),c)  # matrix (sum on j)