数组序列的 numpy 外积

numpy outerproduct of sequence of arrays

我有一个矩阵 A (nXm) 。我的最终目标是获得维度的 Z (nXmXm) 目前我正在使用它来完成它但是可以在不使用 for 循环的情况下使用一些 matrix.tensordot 或 matrix.multiply.outer

 for i in range(0,A.shape[0]):
       Z[i,:,:] = np.outer(A[i,:],A[i,:])

您可以像这样使用 numpy's Einstein summation

np.einsum('ij, ik -> ijk', a, a)

为了完整起见,与来自 unutbu 的同样出色的答案 (+1) 进行时间比较:

In [39]: A = np.random.random((1000,50))

In [40]: %timeit using_einsum(A)
100 loops, best of 3: 11.6 ms per loop

In [41]: %timeit using_broadcasting(A)
100 loops, best of 3: 10.2 ms per loop

In [42]: %timeit orig(A)
10 loops, best of 3: 27.8 ms per loop

教我的

  1. unutbu 的机器比我的快
  2. 广播会比np.einsum
  3. 稍快
for i in range(0,A.shape[0]):
    Z[i,:,:] = np.outer(A[i,:],A[i,:])

表示

Z_ijk = A_ij * A_ik

可以使用NumPy broadcasting计算:

Z = A[:, :, np.newaxis] * A[:, np.newaxis, :]

A[:, :, np.newaxis] has shape (n, m, 1)A[:, np.newaxis, :] 有形状 (n, 1, m)。将两者相乘导致两个阵列被广播到 形状 (n, m, m).

NumPy 乘法总是按元素执行。沿线的价值观 广播轴在任何地方都相同,因此元素乘法结果 在 Z_ijk = A_ij * A_ik.


import numpy as np

def orig(A):
    Z = np.empty(A.shape+(A.shape[-1],), dtype=A.dtype)
    for i in range(0,A.shape[0]):
        Z[i,:,:] = np.outer(A[i,:],A[i,:])
    return Z

def using_broadcasting(A):
    return A[:, :, np.newaxis] * A[:, np.newaxis, :]

这是一个完整性检查,显示这会产生正确的结果:

A = np.random.random((1000,50))
assert np.allclose(using_broadcasting(A), orig(A))

通过选择 A.shape[0] 变大,我们得到一个展示 在 Python:

中广播优于循环的优势
In [107]: %timeit using_broadcasting(A)
10 loops, best of 3: 6.12 ms per loop

In [108]: %timeit orig(A)
100 loops, best of 3: 16.9 ms per loop