稀疏矩阵外积求和

Question

我正在尝试使用scipy的稀疏包来实现以下方程：

W = x[:,1] * y[:,1].T + x[:,2] * y[:,2].T + ...

这里的x和y是一个nxm的csc_matrix。简单来说，我想把x的每一列和y的每一列相乘，然后把得到的nxn矩阵加在一起。最后，我想把所有非零的元素都变成1。

这是我现在的实现方式：

    c = sparse.csc_matrix((n, n))
    for i in xrange(0,m):
        tmp = bam.id2sym_thal[:,i] * bam.id2sym_cort[:,i].T
        minimum(tmp.data,ones_like(tmp.data),tmp.data)
        maximum(tmp.data,ones_like(tmp.data),tmp.data)

        c = c + tmp

这个实现有以下几个问题：

1. 内存使用量似乎暴涨。根据我的理解，内存应该只会在c变得不那么稀疏时增加，但我发现循环开始消耗超过20GB的内存，n=10,000，m=100,000（x和y的每一行大约只有60个非零元素）。

2. 我使用的是Python循环，这样效率不高。

我的问题是：有没有更好的方法来做到这一点？控制内存使用是我最关心的，但如果能让速度更快就更好了！

谢谢！

Answer 1

从内存和性能的角度来看，这可能是使用 Cython 的一个很好的选择。

下面的论文中有一部分专门讲述了它在稀疏 scipy 矩阵中的应用：

http://folk.uio.no/dagss/cython_cise.pdf

Answer 2

请注意，你提到的外积求和其实和把两个矩阵相乘是一样的。换句话说，

sum_i X[:,i]*Y[:,i].T == X*Y.T

所以只需要把这两个矩阵相乘就可以了。

Z = X*Y.T

对于n=10000和m=100000的情况，并且在X和Y的每一列中都有一个非零元素，我的笔记本电脑几乎瞬间就能计算出来。