仅使用存储元素的 Scipy 稀疏矩阵分配答案

【问题标题】：Scipy sparse matrix assignment using only stored elements仅使用存储元素的 Scipy 稀疏矩阵分配
【发布时间】：2016-04-15 15:31:47
【问题描述】：

我有一个大的稀疏矩阵 globalGrid (lil_matrix) 和一个较小的矩阵 localGrid (coo_matrix)。 localGrid 代表 globalGrid 的一个子集，我想用 localGrid 更新 globalGrid。为此，我使用以下代码（在 Python Scipy 中）：

globalGrid[xLocalgrid:xLocalgrid + localGrid.shape[0], yLocalgrid: yLocalgrid + localGrid.shape[1]] = localGrid

其中 xLocalGrid 和 yLocalGrid 是 localGrid 原点相对于 globalGrid 的偏移量。

问题是localGrid是稀疏的，而且零元素也分配给globalGrid。有没有办法只能分配存储的元素而不是 0 元素？

我在 numpy 中发现了掩码数组，但这似乎不适用于稀疏 scipy 矩阵。

编辑：针对下面的cmets，这里有一个例子来说明我的意思：

首先设置矩阵：

M=sparse.lil_matrix(2*np.ones([5,5]))
m = sparse.eye(3)

M.todense()
matrix([[ 2.,  2.,  2.,  2.,  2.],
    [ 2.,  2.,  2.,  2.,  2.],
    [ 2.,  2.,  2.,  2.,  2.],
    [ 2.,  2.,  2.,  2.,  2.],
    [ 2.,  2.,  2.,  2.,  2.]])

m.todense()
matrix([[ 1.,  0.,  0.],
    [ 0.,  1.,  0.],
    [ 0.,  0.,  1.]])

然后赋值：

M[1:4, 1:4] = m

现在结果是：

M.todense()
matrix([[ 2.,  2.,  2.,  2.,  2.],
    [ 2.,  1.,  0.,  0.,  2.],
    [ 2.,  0.,  1.,  0.,  2.],
    [ 2.,  0.,  0.,  1.,  2.],
    [ 2.,  2.,  2.,  2.,  2.]])

而我需要的结果是：

matrix([[ 2.,  2.,  2.,  2.,  2.],
    [ 2.,  1.,  2.,  2.,  2.],
    [ 2.,  2.,  1.,  2.,  2.],
    [ 2.,  2.,  2.,  1.,  2.],
    [ 2.,  2.,  2.,  2.,  2.]])

【问题讨论】：

此分配是否有效，globalGrid 中的值正确？您是否担心结果中有不必要的 0 项？你看到那些我的.data 属性了吗？您可能需要显示示例矩阵。
您能给我们一个minimal reproducible example 和确切的预期输出吗？
我在问题中添加了一个我想要实现的示例。

标签： python numpy matrix scipy sparse-matrix

【解决方案1】：

应该这条线

问题是localGrid是稀疏的，而且非零元素也被分配给globalGrid。有没有办法只能分配存储的元素而不是 0 元素？

改成？

问题在于localGrid 是稀疏的，而且零元素也被分配给了globalGrid。有没有办法只能分配存储的元素而不是 0 元素？

您的问题不太清楚，但我猜是因为 globalGrid[a:b, c:d] 索引跨越两个数组中应该为 0 的值，所以您担心 0 被复制。

让我们用真实的矩阵试试这个。

In [13]: M=sparse.lil_matrix((10,10))
In [14]: m=sparse.eye(3)
In [15]: M[4:7,5:8]=m
In [16]: m
Out[16]: 
<3x3 sparse matrix of type '<class 'numpy.float64'>'
    with 3 stored elements (1 diagonals) in DIAgonal format>
In [17]: M
Out[17]: 
<10x10 sparse matrix of type '<class 'numpy.float64'>'
    with 3 stored elements in LInked List format>
In [18]: M.data
Out[18]: array([[], [], [], [], [1.0], [1.0], [1.0], [], [], []], dtype=object)
In [19]: M.rows
Out[19]: array([[], [], [], [], [5], [6], [7], [], [], []], dtype=object)

M 没有任何不必要的 0。

如果稀疏矩阵中有不必要的 0，则应该使用csr 格式来处理它们

M.tocsr().tolil()

csr 格式也有一个就地的.eliminate_zeros() 方法。

所以您担心的是过度写入目标数组的非零值。

对于密集数组，使用nonzero（或where）可以解决这个问题：

In [87]: X=np.ones((10,10),int)*2
In [88]: y=np.eye(3)
In [89]: I,J=np.nonzero(y)
In [90]: X[I+3,J+2]=y[I,J]
In [91]: X
Out[91]: 
array([[2, 2, 2, 2, 2, 2, 2, 2, 2, 2],
       [2, 2, 2, 2, 2, 2, 2, 2, 2, 2],
       [2, 2, 2, 2, 2, 2, 2, 2, 2, 2],
       [2, 2, 1, 2, 2, 2, 2, 2, 2, 2],
       [2, 2, 2, 1, 2, 2, 2, 2, 2, 2],
       [2, 2, 2, 2, 1, 2, 2, 2, 2, 2],
       [2, 2, 2, 2, 2, 2, 2, 2, 2, 2],
       [2, 2, 2, 2, 2, 2, 2, 2, 2, 2],
       [2, 2, 2, 2, 2, 2, 2, 2, 2, 2],
       [2, 2, 2, 2, 2, 2, 2, 2, 2, 2]])

尝试稀疏等效：

In [92]: M=sparse.lil_matrix(X)
In [93]: M
Out[93]: 
<10x10 sparse matrix of type '<class 'numpy.int32'>'
    with 100 stored elements in LInked List format>
In [94]: m=sparse.coo_matrix(y)
In [95]: m
Out[95]: 
<3x3 sparse matrix of type '<class 'numpy.float64'>'
    with 3 stored elements in COOrdinate format>
In [96]: I,J=np.nonzero(m)
In [97]: I
Out[97]: array([0, 1, 2], dtype=int32)
In [98]: J
Out[98]: array([0, 1, 2], dtype=int32)
In [99]: M[I+3,J+2]=m[I,J]
...
TypeError: 'coo_matrix' object is not subscriptable

我本可以使用自己的稀疏矩阵nonzero。

In [106]: I,J=m.nonzero()

对于coo 格式，与此相同

In [109]: I,J=m.row, m.col

在这种情况下，我也可以使用data 属性：

In [100]: M[I+3,J+2]=m.data
In [101]: M.A
Out[101]: 
array([[2, 2, 2, 2, 2, 2, 2, 2, 2, 2],
       [2, 2, 2, 2, 2, 2, 2, 2, 2, 2],
       [2, 2, 2, 2, 2, 2, 2, 2, 2, 2],
       [2, 2, 1, 2, 2, 2, 2, 2, 2, 2],
       [2, 2, 2, 1, 2, 2, 2, 2, 2, 2],
       [2, 2, 2, 2, 1, 2, 2, 2, 2, 2],
       [2, 2, 2, 2, 2, 2, 2, 2, 2, 2],
       [2, 2, 2, 2, 2, 2, 2, 2, 2, 2],
       [2, 2, 2, 2, 2, 2, 2, 2, 2, 2],
       [2, 2, 2, 2, 2, 2, 2, 2, 2, 2]], dtype=int32)

m.nonzero 的代码可能具有指导意义

    A = self.tocoo()
    nz_mask = A.data != 0
    return (A.row[nz_mask],A.col[nz_mask])

所以你需要小心确保稀疏矩阵的索引和数据属性匹配。

还要注意哪些稀疏格式允许索引。 lil 适合更改值。 csr 允许逐个元素索引，但如果您尝试将零值更改为非零（或 v.v.），则会引发效率警告。 coo 有很好的索引和数据配对，但不允许索引。

另一个微妙之处：在构造coo 时，您可能会重复坐标。当转换为csr 格式时，这些值会相加。但是我建议的分配只会使用最后一个值，而不是总和。因此，请确保您了解您的 local 矩阵是如何构建的，并知道它是否是数据的“干净”表示。

【讨论】：

您建议的更改是正确的。我已经相应地编辑了这个问题。你的答案不适用于我想要的。我已经编辑了问题并添加了一个示例，以使我的愿望更加清晰。
所以您担心global 中的零会覆盖非零。然后你需要使用一些版本的nonzero 索引。用dense 试试，看看sparse 有什么用。
查看我的编辑以使用nonzero 来索引插入点。
我接受了您的回答，因为它有效且简单。但是，我在下面发布的方法在我的（小）测试用例上快了大约 30%。

【解决方案2】：

我找到了一个可行的解决方案。我没有使用赋值，而是遍历稀疏矩阵中的数据（使用 M.data 和 M.rows）并逐个替换元素。

for idx, row in enumerate(m.rows):
    for idy, col in enumerate(row):
        M[yOffset+col, xOffset+idx] = m.data[idx][idy]

我仍然很好奇，如果没有更简单/更快的方法来实现这个结果。

我对上述答案的实现：

I,J = m.nonzero()
M[I+yOffset,J+xOffset] = m[I,J]
return M'

然而，这稍微慢了一点。

【讨论】：