假设我有两个密集矩阵U(10000×50)和V(50×10000),以及一个稀疏矩阵A(10000×10000). A中的每个元素都是1或0.我希望找到A *(UV),注意’*’是逐元素乘法.为了解决这个问题,Scipy / numpy将首先计算一个密集的矩阵UV.但紫外线密集而且很大(10000×10000)所以它非常慢.
因为我只需要A指示的一些UV元素,如果只计算必要的元素而不是计算所有元素然后用A过滤,它应该节省很多时间.有没有办法指示scipy这样做?
顺便说一句,我使用Matlab来解决这个问题,Matlab非常聪明,可以找到我正在尝试做的事情,并且有效地工作.
更新:
我发现Matlab完全按照scipy计算了UV.我的scipy安装太慢了……
最佳答案 这是一个测试脚本和可能的加速.基本思想是使用A的非零坐标来选择U和V的行和列,然后使用einsum来执行可能的点积的子集.
import numpy as np
from scipy import sparse
#M,N,d = 10,5,.1
#M,N,d = 1000,50,.1
M,N,d = 5000,50,.01 # about the limit for my memory
A=sparse.rand(M,M,d)
A.data[:] = 1 # a sparse 0,1 array
U=(np.arange(M*N)/(M*N)).reshape(M,N)
V=(np.arange(M*N)/(M*N)).reshape(N,M)
A1=A.multiply(U.dot(V)) # the direct solution
A2=np.einsum('ij,ik,kj->ij',A.A,U,V)
print(np.allclose(A1,A2))
def foo(A,U,V):
# use A to select elements of U and V
A3=A.copy()
U1=U[A.row,:]
V1=V[:,A.col]
A3.data[:]=np.einsum('ij,ji->i',U1,V1)
return A3
A3 = foo(A,U,V)
print(np.allclose(A1,A3.A))
3个解决方案匹配.对于大型阵列,foo比直接解决方案快约2倍.对于小尺寸,纯净的einsum是有竞争力的,但是对于大型阵列来说是困难的.
在foo中使用dot会计算太多的产品,ij,jk-> ik而不是ij,ji-> i.