沿轴 0 重复一个 scipy csr 稀疏矩阵

Question

我想重复 scipy csr 稀疏矩阵的行，但是当我尝试调用 numpy 的 repeat 方法时，它只是将稀疏矩阵视为对象，并且只会将其作为 ndarray 中的对象重复。我查看了文档，但找不到任何实用程序来重复 scipy csr 稀疏矩阵的行。

我写了以下对内部数据进行操作的代码，似乎可以工作

def csr_repeat(csr, repeats):
    if isinstance(repeats, int):
        repeats = np.repeat(repeats, csr.shape[0])
    repeats = np.asarray(repeats)
    rnnz = np.diff(csr.indptr)
    ndata = rnnz.dot(repeats)
    if ndata == 0:
        return sparse.csr_matrix((np.sum(repeats), csr.shape[1]),
                                 dtype=csr.dtype)
    indmap = np.ones(ndata, dtype=np.int)
    indmap[0] = 0
    rnnz_ = np.repeat(rnnz, repeats)
    indptr_ = rnnz_.cumsum()
    mask = indptr_ < ndata
    indmap -= np.int_(np.bincount(indptr_[mask],
                                  weights=rnnz_[mask],
                                  minlength=ndata))
    jumps = (rnnz * repeats).cumsum()
    mask = jumps < ndata
    indmap += np.int_(np.bincount(jumps[mask],
                                  weights=rnnz[mask],
                                  minlength=ndata))
    indmap = indmap.cumsum()
    return sparse.csr_matrix((csr.data[indmap],
                              csr.indices[indmap],
                              np.r_[0, indptr_]),
                             shape=(np.sum(repeats), csr.shape[1]))

并且效率相当高，但我宁愿不猴子修补课程。有没有更好的方法来做到这一点？

编辑

当我重新审视这个问题时，我想知道为什么我首先发布它。几乎所有我能想到的用重复矩阵做的事情都会更容易用原始矩阵做，然后再应用重复。我的假设是，帖子重复总是比任何潜在答案更好的方法来解决这个问题。

Answer 1

重复稀疏矩阵的最有效方法之一是 OP 建议的方式。我修改了indptr，使其不输出0行。

## original sparse matrix
indptr = np.array([0, 2, 3, 6])
indices = np.array([0, 2, 2, 0, 1, 2])
data = np.array([1, 2, 3, 4, 5, 6])
x = scipy.sparse.csr_matrix((data, indices, indptr), shape=(3, 3))
x.toarray()

array([[1, 0, 2],
       [0, 0, 3],
       [4, 5, 6]])

要重复此操作，您需要重复数据和索引，并且需要修复 indptr。这不是最优雅的方式，但确实有效。

## repeated sparse matrix
repeat = 5 
new_indptr = indptr
for r in range(1,repeat):
    new_indptr = np.concatenate((new_indptr, new_indptr[-1]+indptr[1:]))
x = scipy.sparse.csr_matrix((np.tile(data,repeat), np.tile(indices,repeat), new_indptr))
x.toarray()

array([[1, 0, 2],
       [0, 0, 3],
       [4, 5, 6],
       [1, 0, 2],
       [0, 0, 3],
       [4, 5, 6],
       [1, 0, 2],
       [0, 0, 3],
       [4, 5, 6],
       [1, 0, 2],
       [0, 0, 3],
       [4, 5, 6],
       [1, 0, 2],
       [0, 0, 3],
       [4, 5, 6]])

Answer 2

在有人发布了关于如何最好地做到这一点的非常聪明的回应后，我重新审视了我最初的问题，看看是否有更好的方法。我想出了另一种有利有弊的方法。我们可以指示 scipy 重用重复行的数据，而不是重复所有数据（就像接受的答案一样），创建类似于原始稀疏数组视图的东西（就像您可能对 @987654321 所做的那样@)。这可以通过简单地平铺indptr 字段来完成。

repeated = sparse.csr_matrix((orig.data, orig.indices, np.tile(orig.indptr, repeat_num)))

此技术重复向量repeat_num 次，同时仅修改indptr。不利的一面是，由于 csr 矩阵对数据进行编码的方式，它不是创建一个维度为 repeat_num x n 的矩阵，而是创建一个 (2 * repeat_num - 1) x n 的矩阵，其中每个奇数行都是 0。这应该没什么大不了的因为任何操作都会很快，因为每一行都是 0，并且之后它们应该很容易被切掉（使用类似 [::2] 的东西），但这并不理想。

我认为标记的答案可能仍然是做到这一点的“最佳”方式。

Answer 3

from scipy.sparse import csr_matrix
repeated_row_matrix = csr_matrix(np.ones([repeat_number,1])) * sparse_row

Answer 4

np.repeat 不起作用并不奇怪。它将操作委托给硬编码的a.repeat 方法，如果失败，首先将a 转换为一个数组（如果需要则为对象）。

在开发稀疏代码的线性代数世界中，大部分组装工作在创建稀疏矩阵之前在row、col、data 数组上完成。重点是高效的数学运算，而不是添加/删除/索引行和元素。

我还没有完成你的代码，但我并不感到惊讶 csr 格式矩阵需要这么多工作。

我为lil 格式制定了一个类似的函数（从lil.copy 工作）：

def lil_repeat(S, repeat):
    # row repeat for lil sparse matrix
    # test for lil type and/or convert
    shape=list(S.shape)
    if isinstance(repeat, int):
        shape[0]=shape[0]*repeat
    else:
        shape[0]=sum(repeat)
    shape = tuple(shape)
    new = sparse.lil_matrix(shape, dtype=S.dtype)
    new.data = S.data.repeat(repeat) # flat repeat
    new.rows = S.rows.repeat(repeat)
    return new

但也可以使用索引重复。 lil 和 csr 都支持接近常规 numpy 数组的索引（至少在足够新的版本中）。因此：

S = sparse.lil_matrix([[0,1,2],[0,0,0],[1,0,0]])
print S.A.repeat([1,2,3], axis=0)
print S.A[(0,1,1,2,2,2),:]
print lil_repeat(S,[1,2,3]).A
print S[(0,1,1,2,2,2),:].A

给出相同的结果

最重要的是？

print S[np.arange(3).repeat([1,2,3]),:].A