python - scipy.sparse 矩阵的索引操作向量化

Question

尽管一切似乎都已矢量化，但以下代码运行速度太慢。

from numpy import *
from scipy.sparse import *

n = 100000;
i = xrange(n); j = xrange(n);
data = ones(n);

A=csr_matrix((data,(i,j)));

x = A[i,j]

问题似乎是索引操作是作为 python 函数实现的，调用 A[i,j] 会导致以下分析输出

         500033 function calls in 8.718 CPU seconds

   Ordered by: internal time

   ncalls  tottime  percall  cumtime  percall filename:lineno(function)
   100000    7.933    0.000    8.156    0.000 csr.py:265(_get_single_element)
        1    0.271    0.271    8.705    8.705 csr.py:177(__getitem__)
(...)

也就是说，python 函数 _get_single_element 被调用了 100000 次，效率非常低。为什么这不是用纯 C 实现的？有人知道绕过这个限制并加速上述代码的方法吗？我应该使用不同的稀疏矩阵类型吗？

Answer 1

您可以使用 A.diagonal() 更快地检索对角线(0.0009 秒对我机器上的 3.8 秒)。但是，如果您想进行任意索引，那么这是一个更复杂的问题，因为您使用的不是切片而是索引列表。 _get_single_element 函数被调用了 100000 次，因为它只是迭代您传递给它的迭代器(i 和 j)。切片将是 A[30:60,10] 或类似的东西。

另外，为了简单起见，我会使用 csr_matrix(eye(n,n)) 来制作您使用迭代器制作的相同矩阵。

更新:

好的，既然你的问题真的是关于能够快速访问大量随机元素，我会尽可能地回答你的问题。

• 为什么不使用纯 C 实现？

答案很简单:还没有人了解它。据我所知，Scipy 的稀疏矩阵模块领域还有很多工作要做。用C语言实现的一个部分是不同稀疏矩阵格式之间的转换。

• 有没有人知道绕过这个限制并加速上述代码的方法？

您可以尝试深入研究稀疏矩阵模块并尝试加快它们的速度。我这样做了，并且在尝试上面的代码使用 csr 矩阵进行随机访问时，能够将时间减少到原来的三分之一以下。我必须直接访问 _get_single_element 并显着削减代码来做到这一点，包括取消绑定(bind)检查。

但是，使用 lil_matrix 仍然更快(尽管初始化矩阵速度较慢)，但我必须使用列表推导式进行访问，因为 lil 矩阵未针对您正在执行的索引类型设置。顺便说一句，对 csr_matrix 使用列表推导仍然使 lil 矩阵方法遥遥领先。最后，lil 矩阵访问随机元素的速度更快，因为它没有被压缩。

使用原始形式的 lil_matrix 的运行时间大约是上面列出的代码的五分之一。如果我去掉一些绑定(bind)检查并直接调用 lil_matrix 的 _get1() 方法，我可以将时间进一步缩短为原始时间的 7% 左右。为清楚起见，这是从 3.4-3.8 秒加速到大约 0.261 秒。

最后，我尝试制作自己的函数来直接访问 lil 矩阵的数据并避免重复函数调用。这个时间大约是 0.136 秒。这没有利用正在排序的数据，这是另一个潜在的优化(特别是如果您正在访问同一行上的许多元素)。

如果您想要比这更快，那么您可能必须编写自己的 C 代码稀疏矩阵实现。

• 我应该使用不同的稀疏矩阵类型吗？

好吧，如果您打算访问大量元素，我建议使用 lil 矩阵，但这完全取决于您需要做什么。例如，您还需要乘以矩阵吗？请记住，矩阵之间的更改至少有时(在某些情况下)非常快，因此不排除更改为不同的矩阵格式来执行不同的操作。

如果您实际上不需要对矩阵执行任何代数运算，那么也许您应该只使用 defaultdict 或类似的东西。 defaultdicts 的危险在于，无论何时请求一个不在字典中的元素，它都会将该项目设置为默认值并存储它，这样可能会出现问题。