Scipy getrow() 的稀疏矩阵替代方案
Scipy sparse matrix alternative for getrow()
我正在处理大型稀疏二进制矩阵。我使用 Scipy 稀疏矩阵实现对它们进行了压缩。从scipy.spatial.distance计算Jaccard distance不支持对稀疏矩阵直接运算，所以要么：
将整个稀疏矩阵转换为密集矩阵，然后将每一行作为向量进行操作，这内存不足
或
遍历稀疏，使用getrow()抓取每一行并进行操作。
或
编写我们自己的实现来处理稀疏矩阵。
为了正确看待事物，这里是示例代码：
import scipy.spatial.distance as d
import numpy as np
from scipy.sparse import csr_matrix
# benchmark performance
X = np.random.random((3000, 3000))
# binarize
X[X > 0.3] = 0
X[X>0] = 1
mat = csr_matrix(X)
a = np.zeros(3000)
a[4] = a[100] = a[22] =1
a = csr_matrix(a)
def jaccard_fast(v1,v2):
common = v1.dot(v2.T)
dis = (v1 != v2).getnnz()
if common[0,0]:
return 1.0-float(common[0,0])/float(common[0,0]+dis)
else:
return 0.0
def benchmark_jaccard_fast():
for i in range(mat.shape[0]):
jaccard_fast(mat.getrow(i),a)
def benchmark_jaccard_internal_todense():
for v1,v2 in zip(mat.todense(),a.todense()):
d.jaccard(v1,v2)
def benchmark_jaccard_internal_getrow():
for i in range(mat.shape[0]):
d.jaccard(mat.getrow(i),a)
print "Jaccard Fast:"
%time benchmark_jaccard_fast()
print "Jaccard Scipy (expanding to dense):"
%time benchmark_jaccard_internal_todense()
print "Jaccard Scipy (using getrow):"
%time benchmark_jaccard_internal_getrow()
其中 jaccard_fast 是我自己的实现。在 scipy 稀疏矩阵上，我的实现似乎比内部实现更快，但是 getrow() 似乎减慢了我的实现速度。当我将 jaccard_fast 与 scipy.spatial.distance.jaccard 进行基准测试时，结果是：
Jaccard Fast:
CPU times: user 1.28 s, sys: 0 ns, total: 1.28 s
Wall time: 1.28 s
Jaccard Scipy (expanding to dense):
CPU times: user 28 ms, sys: 8 ms, total: 36 ms
Wall time: 37.2 ms
Jaccard Scipy (using getrow):
CPU times: user 1.82 s, sys: 0 ns, total: 1.82 s
Wall time: 1.81 s
任何有关如何避免 getrow 瓶颈的帮助将不胜感激。由于内存限制，我无法使用 todense() 扩展我的稀疏矩阵。
稀疏索引以速度较慢而著称，例如
In [33]: timeit for row in mat: x=row # sparse iteration
1 loops, best of 3: 510 ms per loop
In [35]: timeit for row in mat.todense(): x=row # dense iteration
10 loops, best of 3: 175 ms per loop
但我发现你的 d.jacard 在使用稀疏矩阵时速度也较慢
In [36]: ad=a.todense()
In [37]: timeit for row in mat.todense(): d.jaccard(row,ad) # all dense
1 loops, best of 3: 734 ms per loop
In [38]: timeit for row in mat: d.jaccard(row.todense(),ad) # inner dense
1 loops, best of 3: 1.69 s per loop
In [39]: timeit for row in mat: d.jaccard(row,a) # all sparse
1 loops, best of 3: 4.61 s per loop
消除 getrow 因素
In [40]: mrow=mat.getrow(0)
In [41]: mrowd=mrow.todense()
In [42]: timeit d.jaccard(mrow, a) # one sparse row
1000 loops, best of 3: 1.32 ms per loop
In [43]: timeit d.jaccard(mrow.todense(), a.todense()) # with conversion
1000 loops, best of 3: 539 µs per loop
In [44]: timeit d.jaccard(mrowd, ad) # dense
10000 loops, best of 3: 173 µs per loop
======================
我需要重新运行 那些测试，因为 d.jaccard 不适用于稀疏（并且 jaccard_fast 不适用于密集）。因此，将稀疏行迭代问题与 jaccard 计算分开将需要更多工作。
我修改了 jaccard_fast 一点：
def my_jaccard(mat, a):
common = mat*a.T # sparse does the large matrix product well
dis=np.array([(row!=a).getnnz() for row in mat]) # iterative
cA = common.A.ravel()
return 1 - cA/(cA + dis)
它 returns 值与密集行上的 d.jaccard 运行 匹配。 d.jaccard returns 1 对于 common 为 0 的行。我似乎不需要 cA 测试（除非 cA 和 dis 在同一槽位为 0）。
In [141]: r=np.array([d.jaccard(row,ad) for row in mat.todense()])
In [142]: r1=my_jaccard(mat,a)
In [143]: np.allclose(r,r1)
Out[143]: True
而且速度只有一半。如果我可以返工 dis 计算器应该有类似的速度。
In [144]: timeit r=np.array([d.jaccard(row,ad) for row in mat.todense()])
1 loops, best of 3: 783 ms per loop
In [145]: timeit r1=my_jaccard(mat,a)
1 loops, best of 3: 1.29 s per loop
进一步调整计算。我屏蔽掉了 0 的 common 值。这有两个目的 - 它确保我们没有除以 0 的问题，并且它减少了 dis 迭代的次数，从而略微提高了速度.
def my_jaccard(mat, a):
common=mat*a.T
cA = common.A.ravel()
mask = cA!=0
cA = cA[mask]
dis = np.array([(row!=a).getnnz() for row, b in zip(mat,mask) if b])
ret = np.ones(mat.shape[0])
ret[mask] = 1 - cA/(cA+dis)
return ret
有了这个，时间会减少一点。
In [188]: timeit my_jaccard(mat,a)
1 loops, best of 3: 1.04 s per loop
==================
的问题有重叠
在那个问题中，我比较了稀疏矩阵和 1 行矩阵，发现使用 sparse.kron 复制行是复制 numpy 广播的最快方法。
使用 jaccard 中的想法计算 dis 数组
def my_jaccard1(mat, a):
common = mat*a.T
cA = common.A.ravel()
aM = sparse.kron(a,np.ones((mat.shape[0],1),int))
dis = (mat!=aM).sum(1)
ret = 1-cA/(cA+dis.A1)
return ret
有了这个时间显着改善（10 倍）：
In [318]: timeit my_jaccard1(mat,a)
1 loops, best of 3: 97.1 ms per loop
我可以像以前一样应用掩码来防止被零除；但它实际上减慢了计算速度（到 140 毫秒）。
def my_jaccard3(mat, a):
common = mat*a.T
cA = common.A.ravel()
mask = cA!=0
cA = cA[mask]
aM = sparse.kron(a,np.ones((len(cA),1),int))
dis = (mat[mask,:]!=aM).sum(1)
ret = np.ones(mat.shape[0])
ret[mask] = 1 - cA/(cA+dis.A1)
return ret
========================
编辑-疑似病例检测
In [75]: x,y= np.array([1,1,0,0,1,0]), np.array([0,0,1,0,1,0])
In [76]: d.jaccard(x,y)
Out[76]: 0.75
In [78]: jaccard_fast(sparse.csr_matrix(x),sparse.csr_matrix(y))
Out[78]: 0.75
我的版本：
In [79]: my_jaccard(sparse.csr_matrix(x),sparse.csr_matrix(y))
Out[79]: array([ 0.75])
...
In [82]: my_jaccard3(sparse.csr_matrix(x),sparse.csr_matrix(y))
Out[82]: array([ 0.75])
（编辑 - 明确使用 sparse.kron）
试试这个：
def getrow_using_underlying_repr(a):
for i,j in enumerate(a.indptr[:-1]):
a.data[a.indptr[i]:a.indptr[i+1]]
比getrow()快100多倍：
Jaccard Fast:
CPU times: user 1.41 s, sys: 0 ns, total: 1.41 s
Wall time: 1.41 s
Jaccard Scipy (expanding to dense):
CPU times: user 19.9 ms, sys: 29.8 ms, total: 49.6 ms
Wall time: 49.3 ms
Jaccard Scipy (using getrow):
CPU times: user 1.48 s, sys: 1.27 ms, total: 1.48 s
Wall time: 1.48 s
getrow_using_underlying_repr:
CPU times: user 11 µs, sys: 1e+03 ns, total: 12 µs
Wall time: 13.8 µs
看懂代码见文档中上一版实例化中的解释：
https://docs.scipy.org/doc/scipy/reference/generated/scipy.sparse.csr_matrix.html#scipy.sparse.csr_matrix
要更好地理解稀疏 CSR 的基础表示，请参阅：
https://www.geeksforgeeks.org/sparse-matrix-representations-set-3-csr/amp/?ref=rp
通过查看源码中built-ingetrow()的实现，发现慢的原因是没有使用底层表示，而是CSR稀疏矩阵的矩阵乘积。
我正在处理大型稀疏二进制矩阵。我使用 Scipy 稀疏矩阵实现对它们进行了压缩。从scipy.spatial.distance计算Jaccard distance不支持对稀疏矩阵直接运算，所以要么：
将整个稀疏矩阵转换为密集矩阵，然后将每一行作为向量进行操作，这内存不足
或
遍历稀疏，使用
getrow()抓取每一行并进行操作。或
编写我们自己的实现来处理稀疏矩阵。
为了正确看待事物，这里是示例代码：
import scipy.spatial.distance as d
import numpy as np
from scipy.sparse import csr_matrix
# benchmark performance
X = np.random.random((3000, 3000))
# binarize
X[X > 0.3] = 0
X[X>0] = 1
mat = csr_matrix(X)
a = np.zeros(3000)
a[4] = a[100] = a[22] =1
a = csr_matrix(a)
def jaccard_fast(v1,v2):
common = v1.dot(v2.T)
dis = (v1 != v2).getnnz()
if common[0,0]:
return 1.0-float(common[0,0])/float(common[0,0]+dis)
else:
return 0.0
def benchmark_jaccard_fast():
for i in range(mat.shape[0]):
jaccard_fast(mat.getrow(i),a)
def benchmark_jaccard_internal_todense():
for v1,v2 in zip(mat.todense(),a.todense()):
d.jaccard(v1,v2)
def benchmark_jaccard_internal_getrow():
for i in range(mat.shape[0]):
d.jaccard(mat.getrow(i),a)
print "Jaccard Fast:"
%time benchmark_jaccard_fast()
print "Jaccard Scipy (expanding to dense):"
%time benchmark_jaccard_internal_todense()
print "Jaccard Scipy (using getrow):"
%time benchmark_jaccard_internal_getrow()
其中 jaccard_fast 是我自己的实现。在 scipy 稀疏矩阵上，我的实现似乎比内部实现更快，但是 getrow() 似乎减慢了我的实现速度。当我将 jaccard_fast 与 scipy.spatial.distance.jaccard 进行基准测试时，结果是：
Jaccard Fast:
CPU times: user 1.28 s, sys: 0 ns, total: 1.28 s
Wall time: 1.28 s
Jaccard Scipy (expanding to dense):
CPU times: user 28 ms, sys: 8 ms, total: 36 ms
Wall time: 37.2 ms
Jaccard Scipy (using getrow):
CPU times: user 1.82 s, sys: 0 ns, total: 1.82 s
Wall time: 1.81 s
任何有关如何避免 getrow 瓶颈的帮助将不胜感激。由于内存限制，我无法使用 todense() 扩展我的稀疏矩阵。
稀疏索引以速度较慢而著称，例如
In [33]: timeit for row in mat: x=row # sparse iteration
1 loops, best of 3: 510 ms per loop
In [35]: timeit for row in mat.todense(): x=row # dense iteration
10 loops, best of 3: 175 ms per loop
但我发现你的 d.jacard 在使用稀疏矩阵时速度也较慢
In [36]: ad=a.todense()
In [37]: timeit for row in mat.todense(): d.jaccard(row,ad) # all dense
1 loops, best of 3: 734 ms per loop
In [38]: timeit for row in mat: d.jaccard(row.todense(),ad) # inner dense
1 loops, best of 3: 1.69 s per loop
In [39]: timeit for row in mat: d.jaccard(row,a) # all sparse
1 loops, best of 3: 4.61 s per loop
消除 getrow 因素
In [40]: mrow=mat.getrow(0)
In [41]: mrowd=mrow.todense()
In [42]: timeit d.jaccard(mrow, a) # one sparse row
1000 loops, best of 3: 1.32 ms per loop
In [43]: timeit d.jaccard(mrow.todense(), a.todense()) # with conversion
1000 loops, best of 3: 539 µs per loop
In [44]: timeit d.jaccard(mrowd, ad) # dense
10000 loops, best of 3: 173 µs per loop
======================
我需要重新运行 那些测试，因为 d.jaccard 不适用于稀疏（并且 jaccard_fast 不适用于密集）。因此，将稀疏行迭代问题与 jaccard 计算分开将需要更多工作。
我修改了 jaccard_fast 一点：
def my_jaccard(mat, a):
common = mat*a.T # sparse does the large matrix product well
dis=np.array([(row!=a).getnnz() for row in mat]) # iterative
cA = common.A.ravel()
return 1 - cA/(cA + dis)
它 returns 值与密集行上的 d.jaccard 运行 匹配。 d.jaccard returns 1 对于 common 为 0 的行。我似乎不需要 cA 测试（除非 cA 和 dis 在同一槽位为 0）。
In [141]: r=np.array([d.jaccard(row,ad) for row in mat.todense()])
In [142]: r1=my_jaccard(mat,a)
In [143]: np.allclose(r,r1)
Out[143]: True
而且速度只有一半。如果我可以返工 dis 计算器应该有类似的速度。
In [144]: timeit r=np.array([d.jaccard(row,ad) for row in mat.todense()])
1 loops, best of 3: 783 ms per loop
In [145]: timeit r1=my_jaccard(mat,a)
1 loops, best of 3: 1.29 s per loop
进一步调整计算。我屏蔽掉了 0 的 common 值。这有两个目的 - 它确保我们没有除以 0 的问题，并且它减少了 dis 迭代的次数，从而略微提高了速度.
def my_jaccard(mat, a):
common=mat*a.T
cA = common.A.ravel()
mask = cA!=0
cA = cA[mask]
dis = np.array([(row!=a).getnnz() for row, b in zip(mat,mask) if b])
ret = np.ones(mat.shape[0])
ret[mask] = 1 - cA/(cA+dis)
return ret
有了这个，时间会减少一点。
In [188]: timeit my_jaccard(mat,a)
1 loops, best of 3: 1.04 s per loop
==================
在那个问题中，我比较了稀疏矩阵和 1 行矩阵，发现使用 sparse.kron 复制行是复制 numpy 广播的最快方法。
使用 jaccard 中的想法计算 dis 数组
def my_jaccard1(mat, a):
common = mat*a.T
cA = common.A.ravel()
aM = sparse.kron(a,np.ones((mat.shape[0],1),int))
dis = (mat!=aM).sum(1)
ret = 1-cA/(cA+dis.A1)
return ret
有了这个时间显着改善（10 倍）：
In [318]: timeit my_jaccard1(mat,a)
1 loops, best of 3: 97.1 ms per loop
我可以像以前一样应用掩码来防止被零除；但它实际上减慢了计算速度（到 140 毫秒）。
def my_jaccard3(mat, a):
common = mat*a.T
cA = common.A.ravel()
mask = cA!=0
cA = cA[mask]
aM = sparse.kron(a,np.ones((len(cA),1),int))
dis = (mat[mask,:]!=aM).sum(1)
ret = np.ones(mat.shape[0])
ret[mask] = 1 - cA/(cA+dis.A1)
return ret
========================
编辑-疑似病例检测
In [75]: x,y= np.array([1,1,0,0,1,0]), np.array([0,0,1,0,1,0])
In [76]: d.jaccard(x,y)
Out[76]: 0.75
In [78]: jaccard_fast(sparse.csr_matrix(x),sparse.csr_matrix(y))
Out[78]: 0.75
我的版本：
In [79]: my_jaccard(sparse.csr_matrix(x),sparse.csr_matrix(y))
Out[79]: array([ 0.75])
...
In [82]: my_jaccard3(sparse.csr_matrix(x),sparse.csr_matrix(y))
Out[82]: array([ 0.75])
（编辑 - 明确使用 sparse.kron）
试试这个：
def getrow_using_underlying_repr(a):
for i,j in enumerate(a.indptr[:-1]):
a.data[a.indptr[i]:a.indptr[i+1]]
比getrow()快100多倍：
Jaccard Fast:
CPU times: user 1.41 s, sys: 0 ns, total: 1.41 s
Wall time: 1.41 s
Jaccard Scipy (expanding to dense):
CPU times: user 19.9 ms, sys: 29.8 ms, total: 49.6 ms
Wall time: 49.3 ms
Jaccard Scipy (using getrow):
CPU times: user 1.48 s, sys: 1.27 ms, total: 1.48 s
Wall time: 1.48 s
getrow_using_underlying_repr:
CPU times: user 11 µs, sys: 1e+03 ns, total: 12 µs
Wall time: 13.8 µs
看懂代码见文档中上一版实例化中的解释：
https://docs.scipy.org/doc/scipy/reference/generated/scipy.sparse.csr_matrix.html#scipy.sparse.csr_matrix
要更好地理解稀疏 CSR 的基础表示，请参阅：
https://www.geeksforgeeks.org/sparse-matrix-representations-set-3-csr/amp/?ref=rp
通过查看源码中built-ingetrow()的实现，发现慢的原因是没有使用底层表示，而是CSR稀疏矩阵的矩阵乘积。