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

Scipy sparse matrix assignment using only stored elements

我有一个大的稀疏矩阵 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 矩阵。

编辑:针对下面的评论,这里有一个例子来说明我的意思:

首先设置矩阵:

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.]])

这条线要不要

The problem is that the localGrid is sparse, but also the non-zero elements are assigned to the globalGrid. Is there a way I can only assign the stored elements and not the 0-elements?

改成?

The problem is that the localGrid is sparse, but also the zero elements are assigned to the globalGrid. Is there a way I can only assign the stored elements and not the 0-elements?

你的问题不是很清楚,但我猜是因为 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 格式也有一个 inplace .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 矩阵的构建方式,并了解它是否是数据的 'clean' 表示。

我找到了可行的解决方案。我没有使用赋值,而是遍历稀疏矩阵中的数据(使用 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]

我还是很好奇有没有simpler/faster的方法可以达到这个效果

我对上面答案的实现:

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

但是速度稍慢。