【问题标题】:Numpy equivalent of dot(A,B,3)点(A,B,3)的numpy等价物
【发布时间】:2014-07-26 19:46:52
【问题描述】:

假设我有两个 3 维矩阵,就像这样(取自这个 matlab 示例 http://www.mathworks.com/help/matlab/ref/dot.html):

A = cat(3,[1 1;1 1],[2 3;4 5],[6 7;8 9])
B = cat(3,[2 2;2 2],[10 11;12 13],[14 15; 16 17])

如果我想沿第三维取成对点积,我可以在 matlab 中这样做:

C = dot(A,B,3)

这会给出结果:

C =
  106   140
  178   220

numpy 中的等效操作是什么,最好是矢量化选项,以避免必须编写一个遍历整个数组的双 for 循环。我似乎无法理解 np.tensordotnp.inner 应该做什么,但它们可能是选项。

【问题讨论】:

    标签: python matlab numpy matrix dot-product


    【解决方案1】:
    In [169]:
    
    A = np.dstack([[[1, 1],[1 ,1]],[[2 ,3],[4, 5]],[[6, 7],[8, 9]]])
    B = np.dstack([[[2, 2],[2, 2]],[[10, 11],[12, 13]],[[14, 15], [16, 17]]])
    c=np.tensordot(A, B.T,1)
    np.vstack([np.diag(c[:,i,i]) for i in range(A.shape[0])]).T
    Out[169]:
    array([[106, 140],
           [178, 220]])
    

    但令人惊讶的是它是最慢的:

    In [170]:
    
    %%timeit
    c=np.tensordot(A, B.T,1)
    np.vstack([np.diag(c[:,i,i]) for i in range(A.shape[0])]).T
    10000 loops, best of 3: 95.2 µs per loop
    In [171]:
    
    %timeit np.einsum('i...,i...',a,b)
    100000 loops, best of 3: 6.93 µs per loop
    In [172]:
    
    %timeit inner1d(A,B)
    100000 loops, best of 3: 4.51 µs per loop
    

    【讨论】:

    • 两个n x n 矩阵的一般情况是什么?
    • @flebool,与 joshAdel 的其他答案几乎一样,请参阅编辑。
    • 请注意 A,B 此处的尺寸不正确。在 Matlab 代码中,数组的形状为 2,2,3,而您的数组的形状为 3,2,2,需要更改程序的所有逻辑。
    • @flebool,已更新。令人惊讶的是,tensordot 是最慢的。我实际上避免了 Einsterin 求和,假设它会很慢。
    • 我猜它会因为 for 循环和对diag 的调用而变慢。显然我期待inner1d 更快,但它甚至没有记录,或者至少我找不到它的文档。
    【解决方案2】:

    这里有一个解决方案:

    A = dstack([[[1, 1],[1 ,1]],[[2 ,3],[4, 5]],[[6, 7],[8, 9]]])
    B = dstack([[[2, 2],[2, 2]],[[10, 11],[12, 13]],[[14, 15], [16, 17]]])
    
    C = einsum('...k,...k',A,B)
    

    基本上dstack沿第三个轴串联,(docs),然后你使用numpy提供的强大的爱因斯坦求和工具einsum(docs)

    【讨论】:

    • 假设我已经有了我的 2 个 3-d n x n x k 数组,我需要做的就是 C = einsum('...k,...k',A,B)?这会比双 for 循环更有效吗?
    • 是的,如果你习惯了爱因斯坦求和,它会更高效,更易读。
    【解决方案3】:

    使用 np.einsum:

    In [9]: B = np.array([[[2, 2],[2, 2]],[[10, 11],[12, 13]],[[14, 15],[16, 17]]])
    
    In [10]: A = np.array([[[1, 1],[1, 1]],[[2, 3],[4, 5]],[[6, 7],[8, 9]]]) 
    
    In [11]: np.einsum('i...,i...',A,B)
    Out[11]:
    array([[106, 140],
           [178, 220]])
    

    或者这是另一个有趣的:

    In [37]: from numpy.core.umath_tests import inner1d
    
    In [38]: inner1d(A,B)
    Out[38]:
    array([[106, 140],
           [178, 220]])
    

    编辑响应@flebool 的评论,inner1d 适用于 (2,2,3) 和 (3,2,2) 形状的数组:

    In [41]: A = dstack([[[1, 1],[1 ,1]],[[2 ,3],[4, 5]],[[6, 7],[8, 9]]])
    
    In [42]: B = dstack([[[2, 2],[2, 2]],[[10, 11],[12, 13]],[[14, 15], [16, 17]]])
    
    In [43]: inner1d(A,B)
    Out[43]:
    array([[106, 140],
           [178, 220]])
    

    【讨论】:

    • 请注意A,B 此处的尺寸不正确。在 Matlab 代码中,数组的形状为2,2,3,而您的数组的形状为3,2,2,需要更改程序的所有逻辑
    • 这似乎是最干净的解决方案。与 einsum 相比,它是否有任何性能劣势?
    • @Stankalank 用于 OP 中的测试阵列,einsum = 5.33 µs 每循环,inner1d = 4.3 µs 每循环。对于 (2,2,10000) einsum = 42.6 µs 每个循环和 inner1d = 41.6 µs 每个循环,尽管时间会因硬件和 numpy 版本/编译设置而异
    • inner1d 的文档在哪里?它的__doc__ 很简陋。
    • @hpaulj 我不记得我第一次偶然发现这个方法的地方,但你可以从导入中看到它不是一个明显的方法。您可以在以下位置查看源代码(C 代码):github.com/numpy/numpy/blob/master/numpy/core/src/umath/…
    猜你喜欢
    • 1970-01-01
    • 2010-12-08
    • 2020-09-05
    • 1970-01-01
    • 1970-01-01
    • 2013-08-08
    • 2019-06-29
    • 2016-02-26
    • 1970-01-01
    相关资源
    最近更新 更多