更有效的方法将3D矩阵的每列与3D矩阵的每个片相乘

我有一个8x8x25000阵列W和一个8 x 25000阵列r。我想在每列（8x1）的每个8x8切片的W中多取一个，并将结果保存在Wres中，最终将成为一个8x25000矩阵。更有效的方法将3D矩阵的每列与3D矩阵的每个片相乘

我完成这个使用for循环这样：

for i in range(0,25000): 
    Wres[:,i] = np.matmul(W[:,:,i],res[:,i])

但这是缓慢的，我希望有做到这一点更快的方法。

任何想法？

来源

2017-07-07 gisuser_14213

任何例如输入/输出，你可以分享？ –

另请参见[Numpy元素明智的点积]（https：// stackoverflow。COM /问题/ 41443444/numpy的元素，明智点产品/ 41443497＃41443497） –

只要2个阵列共享相同的1轴长度，Matmul就可以传播。从文档：

如果任一参数是N-D，N> 2，它将被视为驻留在最后两个索引中的矩阵堆栈并相应地进行广播。

因此，你必须之前matmul执行2个操作：

import numpy as np 
a = np.random.rand(8,8,100) 
b = np.random.rand(8, 100)

转置a和b使得第一轴线是100片
添加一个额外的维度b所以b.shape = (100, 8, 1)

Then：

at = a.transpose(2, 0, 1) # swap to shape 100, 8, 8 
bt = b.T[..., None] # swap to shape 100, 8, 1 
c = np.matmul(at, bt)

c现在100, 8, 1，重塑回8, 100：

c = np.squeeze(c).swapaxes(0, 1)

或

c = np.squeeze(c).T

而在去年，一个班轮只是conveniende：

c = np.squeeze(np.matmul(a.transpose(2, 0, 1), b.T[..., None])).T

来源

2017-07-07 22:45:21

使用np.matmul的替代方法是np.einsum，它可以在1个较短且可以说是更适宜的代码行中完成，而且没有方法链接。

实例阵列：

np.random.seed(123) 
w = np.random.rand(8,8,25000) 
r = np.random.rand(8,25000) 
wres = np.einsum('ijk,jk->ik',w,r) 

# a quick check on result equivalency to your loop 
print(np.allclose(np.matmul(w[:, :, 1], r[:, 1]), wres[:, 1])) 
True

定时相当于@了Imanol的解决方案，从而把你的两个选秀权。两者比循环快30倍。在这里，einsum将因为阵列的大小而具有竞争力。如果数组大于这些数组，它可能会胜出，而对于较小的数组则会失败。有关更多信息，请参阅this讨论。

def solution1(): 
    return np.einsum('ijk,jk->ik',w,r) 

def solution2(): 
    return np.squeeze(np.matmul(w.transpose(2, 0, 1), r.T[..., None])).T 

def solution3(): 
    Wres = np.empty((8, 25000)) 
    for i in range(0,25000): 
     Wres[:,i] = np.matmul(w[:,:,i],r[:,i]) 
    return Wres 

%timeit solution1() 
100 loops, best of 3: 2.51 ms per loop 

%timeit solution2() 
100 loops, best of 3: 2.52 ms per loop 

%timeit solution3() 
10 loops, best of 3: 64.2 ms per loop

Credit到：@Divakar

来源

2017-07-07 23:54:31

更有效的方法将3D矩阵的每列与3D矩阵的每个片相乘

回答

相关问题