Scipy稀疏矩阵和Cython

Question

我需要在一个Cython方法中对scipy稀疏矩阵执行一组操作。

为了高效地应用这些操作，我需要访问lil_matrix表示。 Python中的lil（链接列表稀疏矩阵）数据表示使用不同长度的列表的列表。

如何有效地将不同长度的列表传递给cython（而不进行复制）？ 是否有其他方法可以在cython中访问lil矩阵？

- SKV

---

也许这个答案会给你一些启示：https://dev59.com/ZYLba4cB1Zd3GeqPe3Yc#25295789 - Saullo G. P. Castro

小矩阵的问题在于它们的数据内容不像CSR或COO那样是Numpy数组。 - SKV

我知道了...你应该考虑仅在切片和转换为其他稀疏矩阵类型时使用lil_matrix...然后你可以使用那个答案中提出的解决方案，其中数据类型是已知的... - Saullo G. P. Castro

问题在于子矩阵在循环中使用，我需要保留lil-matrix并避免每次支付转换成本。这归结为有效地将列表传递给cython。 - SKV

Answer 1

下面的示例遍历一个lil_matrix并计算每行的总和。

请注意，我没有进行任何声明，尽管它非常快因为Cython已经针对内置类型（如列表）进行了优化。时间也显示在下面...

import time

import numpy as np
cimport numpy as np

from scipy.sparse import lil_matrix

cdef iter_over_lil_matrix(m):
    cdef list sums, data_row
    sums = []
    for data_row in m.data:
        s = 0
        for value in data_row:
            s += value
        sums.append(s)
    return sums

def main():
    a = np.random.random((1e4*1e4))
    a[a>0.1] = 0
    a = a.reshape(1e4,1e4)
    m = lil_matrix(a)

    t0 = time.clock()
    sums = iter_over_lil_matrix(m)
    t1 = time.clock()
    print 'Cython lil_matrix Time', t1-t0

    t0 = time.clock()
    array_sums = a.sum(axis=1)
    t1 = time.clock()
    print 'Numpy ndarray Time', t1-t0

    t0 = time.clock()
    lil_sums = m.sum(axis=1)
    t1 = time.clock()
    print 'lil_matrix Time', t1-t0

    mcsr = m.tocsr()
    t0 = time.clock()
    csr_sums = mcsr.sum(axis=1)
    t1 = time.clock()
    print 'csr_matrix Time', t1-t0

    assert np.allclose(array_sums, sums)
    assert np.allclose(array_sums, np.asarray(lil_sums).flatten())
    assert np.allclose(array_sums, np.asarray(csr_sums).flatten())

秒数 - 仅比超级优化的NumPy慢大约2倍:D，比lil_matrix.sum()方法快得多，因为它在转换为csr_matrix()之前，如@hpaulj所澄清的并由下面的结果确认。请注意，csr_matrix.sum()在列上的总和几乎比密集矩阵的总和快一个数量级。

Cython lil_matrix Time 0.183935034665
Numpy ndarray Time 0.106583238273
lil_matrix Time 2.47158218631
csr_matrix Time 0.0140050888745

会减慢代码的东西：

• 使用for i in range(len(m.data)):与data_row = m.data[i]
• 声明缓冲区，如np.ndarray[object, ndim=1] data，并使用data=m.data

不会影响的东西：

• boundscheck或wraparound

- Saullo G. P. Castro

---

非常感谢，但这里的比较并不公平。您应该将lil-matrix求和与此cython函数进行比较：即m.sum(axis = 1)而不是a.sum(axis=1)。您能否也发布一下那个结果（最好是在一个方阵上）？此外，这里的矩阵a并不是非常稀疏。 - SKV

@Siamak 我已经更新了答案，使用了更稀疏的矩阵，并与 m.sum(axis=1) 进行比较...请注意，所提出的方法越来越接近于 NumPy 在稀疏矩阵上的性能，最终会比 NumPy 更快... - Saullo G. P. Castro

谢谢Saullo，这些数字令人惊讶！我猜scipy性能不佳的一个原因是矩阵不是非常稀疏。 - SKV

“稀疏和”是与适当的全1向量进行点积运算。为了进行这种数学运算，“lil”格式被转换成“csr”格式。对于您提供的示例矩阵，“csr”格式的加法比密集型加法快5倍。“lil”加法速度较慢，因为需要将其转换为“csr”格式。 - hpaulj

@hpaulj 感谢您宝贵的评论。我之前不知道这一点。我已经更新了答案，添加了 csr_matrix 的时间来确认您的评论。 - Saullo G. P. Castro