为什么使用scipy的稀疏csr_matrix进行向量点积比使用numpy的密集数组慢？

Question

我有一个问题，需要从稀疏矩阵中提取一行并将其与密集矩阵的一行进行点积。使用scipy的csr_matrix，这似乎比使用numpy的密集数组乘法慢得多。这让我感到惊讶，因为我原本以为稀疏点积将涉及更少的操作。以下是一个例子：

import timeit as ti

sparse_setup = 'import numpy as np; import scipy.sparse as si;' + \
               'u = si.eye(10000).tocsr()[10];' + \
               'v = np.random.randint(100, size=10000)'

dense_setup  = 'import numpy as np; u = np.eye(10000)[10];' + \
               'v = np.random.randint(100, size=10000)'

ti.timeit('u.dot(v)', setup=sparse_setup, number=100000)
2.788649031019304

ti.timeit('u.dot(v)', setup=dense_setup, number=100000)
2.179030169005273

对于矩阵-向量乘法，稀疏表示方式胜出，但在这种情况下并非如此。我尝试使用csc_matrix，但性能甚至更差：

>>> sparse_setup = 'import numpy as np; import scipy.sparse as si;' + \
...                'u = si.eye(10000).tocsc()[10];' + \
...                'v = np.random.randint(100, size=10000)'
>>> ti.timeit('u.dot(v)', setup=sparse_setup, number=100000)
7.0045155879925005

为什么在这种情况下 numpy 能胜过 scipy.sparse？是否有一种矩阵格式对此类计算更快？

- castle-bravo

---

稀疏矩阵在计算过程中更复杂，但可以节省内存。 - MaxNoe

你是如何计算“操作”的？只有乘法和加法吗？还是考虑了索引、迭代等因素？在现代处理器上，两个数字相乘并不是一项昂贵的操作。此外，“点积”会被委托给快速数值库。稀疏矩阵乘法也会被编译，但使用的优化库不同。 - hpaulj

只是为了明确，您正在测试一个非常稀疏的向量（10000个元素中只有1个非零元素），与相同大小的密集向量进行比较。它最终使用了sparse._sparsetools.csr_matvec，这是一个已编译的函数。我需要查看scipy github以进一步了解。 - hpaulj

@hpaulj，我认为稀疏矩阵是一组行-列-值元组的集合。k个稀疏行向量与n个密集向量的点积应该只需要O(k)次操作（比密集向量具有更高的常数因子），而密集乘法应该需要O(n)次操作。我会采纳您的好建议并阅读源代码。 - castle-bravo

当我改变u中非零值的数量（以及密集等效项）时，时间并没有发生变化。(u.data * v[u.indices]).sum()更接近于您想象中的情况。它略微更快，但其时间仍不是O(k)。 - hpaulj

是的，这正是我所想象的发生的事情。受到您的评论的启发，我尝试使用DOK格式而不是CSR。这导致了3.5倍的加速，并且比密集的点积更快。 - castle-bravo

Answer 1

CSR/CSC向量乘积调用每次需要几微秒的开销，这是由于执行一小段Python代码和在编译代码（scipy.sparse._sparsetools.csr_matvec）中处理参数导致的。在现代处理器上，计算向量点积非常快，因此在这种情况下，开销占据了计算时间的主要部分。矩阵-向量乘积本身更昂贵，在这里类似的开销不可见。为什么Numpy的开销较小？主要是由于代码的优化更好; csr_matrix的性能可能会在这里得到改善。

- pv.