了解NumPy的einsum

Question

我很难理解einsum是如何工作的。我看过文档和一些例子，但它似乎并不坚持。

这是我们在课堂上看到的一个例子：

C = np.einsum("ij,jk->ki", A, B)

对于两个阵列A和B

我认为这需要A^T * B，但我不确定（它正在将其中一个的转置正确吗？）。任何人都可以告诉我这里发生了什么（一般情况下使用einsum）？

Answer 1

（注意：这个答案基于我之前写的关于blog post的简短einsum。）

`einsum`做什么？

想象一下，我们有两个多维数组，A和B。现在让我们假设我们想......

将A与B以特定方式相乘以创建新的产品阵列;然后也许吧
将这个新阵列与特定轴相加;然后也许吧
以特定顺序转置新阵列的轴。

einsum很有可能帮助我们更快，更有效地记忆，如multiply，sum和transpose等NumPy函数的组合将允许。

`einsum`如何运作？

这是一个简单（但不是完全无关紧要）的例子。采用以下两个数组：

A = np.array([0, 1, 2])

B = np.array([[ 0,  1,  2,  3],
              [ 4,  5,  6,  7],
              [ 8,  9, 10, 11]])

我们将元素乘以A和B，然后沿新数组的行求和。在“正常”的NumPy中，我们写道：

>>> (A[:, np.newaxis] * B).sum(axis=1)
array([ 0, 22, 76])

所以在这里，A上的索引操作将两个数组的第一个轴对齐，以便可以广播乘法。然后将产品数组的行相加以返回答案。

现在，如果我们想要使用einsum，我们可以写：

>>> np.einsum('i,ij->i', A, B)
array([ 0, 22, 76])

签名字符串'i,ij->i'是关键，需要一点解释。你可以把它分成两半。在左侧（->的左侧），我们标记了两个输入数组。在->的右边，我们已经标记了我们想要最终得到的数组。

接下来会发生什么：

A有一个轴;我们把它标记为i。 B有两个轴;我们将轴0标记为i，将轴1标记为j。
通过在两个输入数组中重复标签i，我们告诉einsum这两个轴应该相乘。换句话说，我们将数组A与B数组的每一列相乘，就像A[:, np.newaxis] * B一样。
请注意，j不会在我们想要的输出中显示为标签;我们刚刚使用了i（我们希望最终得到一维数组）。通过省略标签，我们告诉einsum沿着这个轴求和。换句话说，我们正在总结产品的行，就像.sum(axis=1)一样。

这基本上是你需要知道使用einsum。它有助于玩一点;如果我们在输出中留下两个标签'i,ij->ij'，我们会得到一个2D数组的产品（与A[:, np.newaxis] * B相同）。如果我们说没有输出标签，'i,ij->，我们得到一个数字（与做(A[:, np.newaxis] * B).sum()相同）。

然而，关于einsum的伟大之处在于，它不是首先构建一系列临时产品;它只是对产品进行总结。这可以大大节省内存使用量。

一个稍大的例子

为了解释点积，这里有两个新的数组：

A = array([[1, 1, 1],
           [2, 2, 2],
           [5, 5, 5]])

B = array([[0, 1, 0],
           [1, 1, 0],
           [1, 1, 1]])

我们将使用np.einsum('ij,jk->ik', A, B)计算点积。这是一张图片，显示了A和B的标签以及我们从函数中获得的输出数组：

你可以看到重复标签j - 这意味着我们将A的行与B的列相乘。此外，标签j不包括在输出中 - 我们正在总结这些产品。保留标签i和k用于输出，因此我们返回2D阵列。

将此结果与标签j未求和的数组进行比较可能更清楚。下面，在左侧，您可以看到编写np.einsum('ij,jk->ijk', A, B)的3D数组（即我们保留了标签j）：

求和轴j给出了预期的点积，如右图所示。

一些练习

为了更好地了解einsum，使用下标符号实现熟悉的NumPy数组操作会很有用。任何涉及乘法和求和轴组合的东西都可以使用einsum编写。

设A和B为两个长度相同的1D阵列。例如，A = np.arange(10)和B = np.arange(5, 15)。

A的总和可以写成： np.einsum('i->', A)
元素乘法，A * B，可以写成： np.einsum('i,i->i', A, B)
内部产品或点产品np.inner(A, B)或np.dot(A, B)可以写成： np.einsum('i,i->', A, B) # or just use 'i,i'
外部产品np.outer(A, B)可以写成： np.einsum('i,j->ij', A, B)

对于2D数组，C和D，假设轴是兼容的长度（两个长度相同或长度为1），这里有几个例子：

C（主对角线的总和），np.trace(C)的痕迹可以写成： np.einsum('ii', C)
C的元素乘法和qazxswpoi，D的转置可以写成： C * D.T
将qazxsw poi的每个元素乘以数组qazxsw poi（制作4D数组），np.einsum('ij,ji->ij', C, D)，可以写成： C

Answer 2

如果你直观地理解D的想法就很容易理解。作为一个例子，让我们从涉及矩阵乘法的简单描述开始。

要使用C[:, :, None, None] * D，您所要做的就是将所谓的下标字符串作为参数传递，然后输入您的输入数组。

假设你有两个二维数组，np.einsum('ij,kl->ijkl', C, D)和numpy.einsum()，你想做矩阵乘法。所以你也是：

numpy.einsum()

下标字符串A对应于数组B，而下标字符串np.einsum("ij, jk -> ik", A, B)对应于数组ij。此外，最重要的是要注意每个下标字符串中的字符数必须与数组的尺寸相匹配。（即2D阵列的两个字符，3D数组的三个字符，等等。）如果你重复下标字符串之间的字符（在我们的例子中是A），那么这意味着你希望jksum沿着那些维度发生。因此，他们将减少总和。（即那个维度将会消失）

这个B之后的下标字符串将是我们的结果数组。如果将其保留为空，则将对所有内容求和，并返回标量值作为结果。否则，结果数组将根据下标字符串具有维度。在我们的例子中，它将是j。这是直观的，因为我们知道对于矩阵乘法，数组ein中的列数必须与数组->中的行数相匹配，这是在这里发生的事情（即我们通过在下标字符串中重复char ik来编码这些知识）

下面是一些例子，说明A在实现一些常见的张量或nd阵列操作时的用法/能力，简洁明了。

输入

1）矩阵乘法（类似于j）

np.einsum()

2）沿主对角线提取元素（类似于# a vector In [197]: vec Out[197]: array([0, 1, 2, 3]) # an array In [198]: A Out[198]: array([[11, 12, 13, 14], [21, 22, 23, 24], [31, 32, 33, 34], [41, 42, 43, 44]]) # another array In [199]: B Out[199]: array([[1, 1, 1, 1], [2, 2, 2, 2], [3, 3, 3, 3], [4, 4, 4, 4]])）

np.matmul(arr1, arr2)

3）Hadamard产品（即两个阵列的元素产品）（类似于In [200]: np.einsum("ij, jk -> ik", A, B) Out[200]: array([[130, 130, 130, 130], [230, 230, 230, 230], [330, 330, 330, 330], [430, 430, 430, 430]])）

np.diag(arr)

4）元素方形（类似于In [202]: np.einsum("ii -> i", A) Out[202]: array([11, 22, 33, 44])或arr1 * arr2）

In [203]: np.einsum("ij, ij -> ij", A, B)
Out[203]: 
array([[ 11,  12,  13,  14],
       [ 42,  44,  46,  48],
       [ 93,  96,  99, 102],
       [164, 168, 172, 176]])

5）跟踪（即主对角元素的总和）（类似于np.square(arr)）

arr ** 2

6）矩阵转置（类似于In [210]: np.einsum("ij, ij -> ij", B, B) Out[210]: array([[ 1, 1, 1, 1], [ 4, 4, 4, 4], [ 9, 9, 9, 9], [16, 16, 16, 16]])）

np.trace(arr)

7）外部产品（矢量）（类似于In [217]: np.einsum("ii -> ", A) Out[217]: 110）

np.transpose(arr)

8）内积（向量）（类似于In [221]: np.einsum("ij -> ji", A) Out[221]: array([[11, 21, 31, 41], [12, 22, 32, 42], [13, 23, 33, 43], [14, 24, 34, 44]])）

np.outer(vec1, vec2)

9）沿轴0的总和（类似于In [255]: np.einsum("i, j -> ij", vec, vec) Out[255]: array([[0, 0, 0, 0], [0, 1, 2, 3], [0, 2, 4, 6], [0, 3, 6, 9]])）

np.inner(vec1, vec2)

10）沿轴1的总和（类似于In [256]: np.einsum("i, i -> ", vec, vec) Out[256]: 14）

np.sum(arr, axis=0)

11）批量矩阵乘法

In [260]: np.einsum("ij -> j", B)
Out[260]: array([10, 10, 10, 10])

12）沿轴2的总和（类似于np.sum(arr, axis=1)）

In [261]: np.einsum("ij -> i", B)
Out[261]: array([ 4,  8, 12, 16])

13）对数组中的所有元素求和（类似于In [287]: BM = np.stack((A, B), axis=0) In [288]: BM Out[288]: array([[[11, 12, 13, 14], [21, 22, 23, 24], [31, 32, 33, 34], [41, 42, 43, 44]], [[ 1, 1, 1, 1], [ 2, 2, 2, 2], [ 3, 3, 3, 3], [ 4, 4, 4, 4]]]) In [289]: BM.shape Out[289]: (2, 4, 4) # batch matrix multiply using einsum In [292]: BMM = np.einsum("bij, bjk -> bik", BM, BM) In [293]: BMM Out[293]: array([[[1350, 1400, 1450, 1500], [2390, 2480, 2570, 2660], [3430, 3560, 3690, 3820], [4470, 4640, 4810, 4980]], [[ 10, 10, 10, 10], [ 20, 20, 20, 20], [ 30, 30, 30, 30], [ 40, 40, 40, 40]]]) In [294]: BMM.shape Out[294]: (2, 4, 4)）

np.sum(arr, axis=2)

14）多轴的和（即边缘化）（类似于In [330]: np.einsum("ijk -> ij", BM) Out[330]: array([[ 50, 90, 130, 170], [ 4, 8, 12, 16]])）

np.sum(arr)

15）In [335]: np.einsum("ijk -> ", BM) Out[335]: 480（类似于np.sum（hadamard-product）cf. 3）

np.sum(arr, axis=(axis0, axis1, axis2, axis3, axis4, axis6, axis7))

16）2D和3D阵列乘法

当求解要验证结果的线性方程组（Ax = b）时，这种乘法可能非常有用。

# 8D array
In [354]: R = np.random.standard_normal((3,5,4,6,8,2,7,9))

# marginalize out axis 5 (i.e. "n" here)
In [363]: esum = np.einsum("ijklmnop -> n", R)

# marginalize out axis 5 (i.e. sum over rest of the axes)
In [364]: nsum = np.sum(R, axis=(0,1,2,3,4,6,7))

In [365]: np.allclose(esum, nsum)
Out[365]: True

相反，如果必须使用Double Dot Products进行此验证，我们必须执行几个In [772]: A Out[772]: array([[1, 2, 3], [4, 2, 2], [2, 3, 4]]) In [773]: B Out[773]: array([[1, 4, 7], [2, 5, 8], [3, 6, 9]]) In [774]: np.einsum("ij, ij -> ", A, B) Out[774]: 124操作以实现相同的结果，如：

# inputs
In [115]: A = np.random.rand(3,3)
In [116]: b = np.random.rand(3, 4, 5)

# solve for x
In [117]: x = np.linalg.solve(A, b.reshape(b.shape[0], -1)).reshape(b.shape)

# 2D and 3D array multiplication :)
In [118]: Ax = np.einsum('ij, jkl', A, x)

# indeed the same!
In [119]: np.allclose(Ax, b)
Out[119]: True

额外奖励：在这里阅读更多数学：np.matmul()，绝对在这里：reshape

Answer 3

让我们制作2个阵列，具有不同但兼容的尺寸，以突出它们的相互作用

# reshape 3D array `x` to 2D, perform matmul
# then reshape the resultant array to 3D
In [123]: Ax_matmul = np.matmul(A, x.reshape(x.shape[0], -1)).reshape(x.shape)

# indeed correct!
In [124]: np.allclose(Ax, Ax_matmul)
Out[124]: True

你的计算，得到a（2,3）与a（3,4）的'dot'（乘积之和）来产生（4,2）数组。 Einstein-Summation是Tensor-Notation的第一个暗淡，In [43]: A=np.arange(6).reshape(2,3) Out[43]: array([[0, 1, 2], [3, 4, 5]]) In [44]: B=np.arange(12).reshape(3,4) Out[44]: array([[ 0, 1, 2, 3], [ 4, 5, 6, 7], [ 8, 9, 10, 11]])的最后一个; i是A的最后一个C。 k被总和'消耗'。

这与C相同 - 它是转置的最终输出。

要查看更多j发生的事情，请将In [45]: C=np.einsum('ij,jk->ki',A,B) Out[45]: array([[20, 56], [23, 68], [26, 80], [29, 92]])下标更改为np.dot(A,B).T：

这也可以通过以下方式生成：

也就是说，在ijk末尾添加In [46]: np.einsum('ij,jk->ijk',A,B) Out[46]: array([[[ 0, 0, 0, 0], [ 4, 5, 6, 7], [16, 18, 20, 22]], [[ 0, 3, 6, 9], [16, 20, 24, 28], [40, 45, 50, 55]]])维度，在A[:,:,None]*B[None,:,:]前面添加k，得到（2,3,4）阵列。

A，i等;汇总B并转置以获得更早的结果：

0 + 4 + 16 = 20

Answer 4

我发现9 + 28 + 55 = 92具有指导意义

我们使用 - >来表示输出数组的顺序。因此，将'ij，i-> j'视为左手侧（LHS）和右手侧（RHS）。 LHS上的任何重复标签都会计算产品元素，然后求和。通过改变RHS（输出）侧的标签，我们可以定义我们想要相对于输入数组进行的轴，即沿轴0,1等的总和。

注意有三个轴，i，j，k，并且重复j（在左侧）。 np.sum(A[:,:,None] * B[None,:,:], axis=1).T # C[k,i] = sum(j) A[i,j (,k) ] * B[(i,) j,k]代表NumPy: The tricks of the trade (Part II)的行和列。 import numpy as np >>> a array([[1, 1, 1], [2, 2, 2], [3, 3, 3]]) >>> b array([[0, 1, 2], [3, 4, 5], [6, 7, 8]]) >>> d = np.einsum('ij, jk->ki', a, b)为i,j。

为了计算产品并对齐a轴，我们需要在j,k上添加一个轴。（b将沿（？）第一轴播出）

a不在右侧，因此我们总结了b，它是3x3x3阵列的第二轴

a[i, j, k]
   b[j, k]

>>> c = a[:,:,np.newaxis] * b
>>> c
array([[[ 0,  1,  2],
        [ 3,  4,  5],
        [ 6,  7,  8]],

       [[ 0,  2,  4],
        [ 6,  8, 10],
        [12, 14, 16]],

       [[ 0,  3,  6],
        [ 9, 12, 15],
        [18, 21, 24]]])

最后，指数在右侧（按字母顺序）反转，因此我们进行转置。

了解NumPy的einsum

问题描述投票：133回答：4

4个回答

`einsum`做什么？

`einsum`如何运作？

一个稍大的例子

一些练习

最新问题

了解NumPy的einsum

问题描述 投票：133回答：4

4个回答

einsum做什么？

einsum如何运作？

一个稍大的例子

一些练习

最新问题

问题描述投票：133回答：4

`einsum`做什么？

`einsum`如何运作？