为什么 scipy csr 矩阵的行索引比 numpy 数组慢

Question

我不确定我做错了什么，但与 numpy 数组相比，索引 scipy csr_matrix 的行似乎慢了大约 2 倍（见下面的代码）。

csr 矩阵的行索引不应该比密集矩阵更快，因为像下面的情况一样只提取了很少的非零元素吗？

是否有技巧可以使 scipy csr 矩阵的行索引更快？

import numpy as np
import timeit
from scipy.sparse import csr_matrix

# Generate random matrix
A = np.random.rand(5000, 1000)

# Make A sparse
A[:, 4:] =0

# Create sparse matrix
A_sparse = csr_matrix(A)

# Create row indexing functions
def index_A_dense():
    A[4]

def index_A_dense_copy():
    a = A[4].copy()

def index_A_sparse():
    A_sparse[4]

timeit.timeit(index_A_sparse, number=10000)
Out: 0.6668063667304978
timeit.timeit(index_A_dense, number=10000)
Out: 0.0029007720424942818
timeit.timeit(index_A_dense_copy, number=10000)
Out: 0.00979283193282754

提前致谢！

Answer 1

我在下面演示的简短回答是构建一个新的稀疏矩阵是昂贵的。有一个显着的开销不依赖于行数或特定行中非零元素的数量。

稀疏矩阵的数据表示与密集数组的数据表示有很大不同。数组将数据存储在一个一个的连续缓冲区中，并有效地使用 shape 和 strides 来迭代选定的值。这些值，加上索引，准确地定义了它将在缓冲区中找到数据的位置。将那些 N 字节从一个位置复制到另一个位置是整个操作中相对较小的部分。

稀疏矩阵将数据存储在多个数组（或其他结构）中，包含索引和数据。然后选择一行需要查找相关索引，并使用选定的索引和数据构建一个新的稀疏矩阵。 sparse 包中有编译后的代码，但没有 numpy 数组那么多的底层代码。

为了说明，我将制作小矩阵，并且不那么密集，所以我们没有很多空行：

In [259]: A = (sparse.rand(5,5,.4,'csr')*20).floor()
In [260]: A
Out[260]: 
<5x5 sparse matrix of type '<class 'numpy.float64'>'
    with 10 stored elements in Compressed Sparse Row format>

密集等价，一行复制：

In [262]: Ad=A.A
In [263]: Ad
Out[263]: 
array([[  0.,   0.,   0.,   0.,  10.],
       [  0.,   0.,   0.,   0.,   0.],
       [ 17.,  16.,  14.,  19.,   6.],
       [  0.,   0.,   1.,   0.,   0.],
       [ 14.,   0.,   9.,   0.,   0.]])
In [264]: Ad[4,:]
Out[264]: array([ 14.,   0.,   9.,   0.,   0.])
In [265]: timeit Ad[4,:].copy()
100000 loops, best of 3: 4.58 µs per loop

一个矩阵行：

In [266]: A[4,:]
Out[266]: 
<1x5 sparse matrix of type '<class 'numpy.float64'>'
    with 2 stored elements in Compressed Sparse Row format>

查看此 csr 矩阵的数据表示，3 个一维数组：

In [267]: A.data
Out[267]: array([  0.,  10.,  17.,  16.,  14.,  19.,   6.,   1.,  14.,   9.])  
In [268]: A.indices
Out[268]: array([3, 4, 0, 1, 2, 3, 4, 2, 0, 2], dtype=int32)
In [269]: A.indptr
Out[269]: array([ 0,  2,  2,  7,  8, 10], dtype=int32)

以下是选择行的方式（但在编译代码中）：

In [270]: A.indices[A.indptr[4]:A.indptr[5]]
Out[270]: array([0, 2], dtype=int32)
In [271]: A.data[A.indptr[4]:A.indptr[5]]
Out[271]: array([ 14.,   9.])

'row' 是另一个稀疏矩阵，具有相同类型的数据数组：

In [272]: A[4,:].indptr
Out[272]: array([0, 2])
In [273]: A[4,:].indices
Out[273]: array([0, 2])
In [274]: timeit A[4,:]

是的，稀疏矩阵的计时很慢。不知道有多少时间花在了实际选择数据上，又有多少时间花在了构建新矩阵上。

10000 loops, best of 3: 145 µs per loop
In [275]: timeit Ad[4,:].copy()
100000 loops, best of 3: 4.56 µs per loop

lil 格式可能更容易理解，因为数据和索引存储在子列表中，每行一个。

In [276]: Al=A.tolil()
In [277]: Al.data
Out[277]: array([[0.0, 10.0], [], [17.0, 16.0, 14.0, 19.0, 6.0], [1.0], [14.0, 9.0]], dtype=object)
In [278]: Al.rows
Out[278]: array([[3, 4], [], [0, 1, 2, 3, 4], [2], [0, 2]], dtype=object)
In [279]: Al[4,:].data
Out[279]: array([[14.0, 9.0]], dtype=object)
In [280]: Al[4,:].rows
Out[280]: array([[0, 2]], dtype=object)

在处理紧凑的编译代码时，这样的速度比较是有意义的，其中字节从内存的一部分移动到另一部分是重要的时间消耗者。在 numpy 和 scipy 中混合使用 Python 和编译代码，您不能只计算 O(n) 操作。

=============================

这是从 A 中选择一行所需的时间估计值，以及 return 一个新的稀疏矩阵所需的时间：

刚拿到数据：

In [292]: %%timeit
d1=A.data[A.indptr[4]:A.indptr[5]]
i1=A.indices[A.indptr[4]:A.indptr[5]]
   .....: 
100000 loops, best of 3: 4.92 µs per loop

加上制作矩阵所需的时间：

In [293]: %%timeit
d1=A.data[A.indptr[4]:A.indptr[5]]
i1=A.indices[A.indptr[4]:A.indptr[5]]
sparse.csr_matrix((d1,([0,0],i1)),shape=(1,5))
   .....: 
1000 loops, best of 3: 445 µs per loop

尝试更简单的 coo 矩阵

In [294]: %%timeit
d1=A.data[A.indptr[4]:A.indptr[5]]
i1=A.indices[A.indptr[4]:A.indptr[5]]
sparse.coo_matrix((d1,([0,0],i1)),shape=(1,5))
   .....: 
10000 loops, best of 3: 135 µs per loop

为什么 scipy csr 矩阵的行索引比 numpy 数组慢

why is row indexing of scipy csr matrices slower compared to numpy arrays

python

time

numpy

scipy

sparse-matrix