Numpy Broadcast 执行欧式距离矢量化
Numpy Broadcast to perform euclidean distance vectorized
我有 2 x 4 和 3 x 4 的矩阵。我想找到跨行的欧氏距离，并在最后得到一个 2 x 3 的矩阵。下面是带有一个 for 循环的代码，它针对所有 b 行向量计算 a 中每个行向量的欧氏距离。如何在不使用 for 循环的情况下做同样的事情？
import numpy as np
a = np.array([[1,1,1,1],[2,2,2,2]])
b = np.array([[1,2,3,4],[1,1,1,1],[1,2,1,9]])
dists = np.zeros((2, 3))
for i in range(2):
dists[i] = np.sqrt(np.sum(np.square(a[i] - b), axis=1))
只需在正确的地方使用np.newaxis：
np.sqrt((np.square(a[:,np.newaxis]-b).sum(axis=2)))
此功能已包含在 scipy's spatial module 中，我建议使用它，因为它将在后台进行矢量化和高度优化。但是，正如另一个答案所表明的那样，您可以通过多种方式自己做到这一点。
import numpy as np
a = np.array([[1,1,1,1],[2,2,2,2]])
b = np.array([[1,2,3,4],[1,1,1,1],[1,2,1,9]])
np.sqrt((np.square(a[:,np.newaxis]-b).sum(axis=2)))
# array([[ 3.74165739, 0. , 8.06225775],
# [ 2.44948974, 2. , 7.14142843]])
from scipy.spatial.distance import cdist
cdist(a,b)
# array([[ 3.74165739, 0. , 8.06225775],
# [ 2.44948974, 2. , 7.14142843]])
我最近在使用深度学习（stanford cs231n，Assignment1）时遇到了同样的问题，但是当我使用
np.sqrt((np.square(a[:,np.newaxis]-b).sum(axis=2)))
出现错误
MemoryError
这意味着我 运行 内存不足（事实上，在 middle.It 中产生了一个 500*5000*1024 的数组，这么大！）
为防止该错误，我们可以使用一个公式来简化：
%5E2&space;=&space;a%5E2&space;-&space;2ab+b%5E2 "(a-b)^2 = a^2 - 2ab+b^2")
代码：
import numpy as np
aSumSquare = np.sum(np.square(a),axis=1);
bSumSquare = np.sum(np.square(b),axis=1);
mul = np.dot(a,b.T);
dists = np.sqrt(aSumSquare[:,np.newaxis]+bSumSquare-2*mul)
这里是原始输入变量：
A = np.array([[1,1,1,1],[2,2,2,2]])
B = np.array([[1,2,3,4],[1,1,1,1],[1,2,1,9]])
A
# array([[1, 1, 1, 1],
# [2, 2, 2, 2]])
B
# array([[1, 2, 3, 4],
# [1, 1, 1, 1],
# [1, 2, 1, 9]])
A 是一个 2x4 数组。
B 是一个 3x4 数组。
我们想在一个完全向量化的运算中计算欧氏距离矩阵运算，其中 dist[i,j] 包含 A 中第 i 个实例和 B 中第 j 个实例之间的距离。因此 dist 是 2x3这个例子。
距离
表面上可以用 numpy 编写为
dist = np.sqrt(np.sum(np.square(A-B))) # DOES NOT WORK
# Traceback (most recent call last):
# File "<stdin>", line 1, in <module>
# ValueError: operands could not be broadcast together with shapes (2,4) (3,4)
但是，如上所示，问题在于元素减法运算A-B涉及不兼容的数组大小，特别是第一维中的2和3。
A has dimensions 2 x 4
B has dimensions 3 x 4
为了进行逐元素减法，我们必须填充 A 或 B 以满足 numpy 的广播规则。我将选择用额外的维度填充 A，使其变为 2 x 1 x 4，这允许数组的维度对齐以进行广播。有关 numpy 广播的更多信息，请参阅 tutorial in the scipy manual and the final example in this tutorial.
您可以使用 np.newaxis 值或 np.reshape 命令执行填充。我在下面显示：
# First approach is to add the extra dimension to A with np.newaxis
A[:,np.newaxis,:] has dimensions 2 x 1 x 4
B has dimensions 3 x 4
# Second approach is to reshape A with np.reshape
np.reshape(A, (2,1,4)) has dimensions 2 x 1 x 4
B has dimensions 3 x 4
如您所见，使用任何一种方法都可以使尺寸对齐。我将使用 np.newaxis 的第一种方法。所以现在，这将用于创建 A-B，这是一个 2x3x4 数组：
diff = A[:,np.newaxis,:] - B
# Alternative approach:
# diff = np.reshape(A, (2,1,4)) - B
diff.shape
# (2, 3, 4)
现在我们可以将差分表达式代入dist等式语句得到最终结果：
dist = np.sqrt(np.sum(np.square(A[:,np.newaxis,:] - B), axis=2))
dist
# array([[ 3.74165739, 0. , 8.06225775],
# [ 2.44948974, 2. , 7.14142843]])
注意 sum 超过了 axis=2，这意味着对 2x3x4 数组的第三个轴（其中轴 id 从 0 开始）求和。
如果你的数组很小，那么上面的命令就可以正常工作。然而，如果你有大数组，那么你可能 运行 进入内存问题。请注意，在上面的示例中，numpy 在内部创建了一个 2x3x4 数组来执行广播。如果我们将 A 的维度概括为 a x z，将 B 的维度概括为 b x z，那么 numpy 将在内部创建一个 a x b x z 数组用于广播。
我们可以通过一些数学运算来避免创建这个中间数组。因为您将欧几里德距离计算为平方差和，所以我们可以利用平方差和可以重写的数学事实。
请注意，中间项涉及 逐元素 乘法的总和。这个乘法求和被称为点积。因为A和B各自都是一个矩阵，那么这个操作其实就是一个矩阵乘法。因此，我们可以将上面的内容重写为：
然后我们可以编写如下的numpy代码：
threeSums = np.sum(np.square(A)[:,np.newaxis,:], axis=2) - 2 * A.dot(B.T) + np.sum(np.square(B), axis=1)
dist = np.sqrt(threeSums)
dist
# array([[ 3.74165739, 0. , 8.06225775],
# [ 2.44948974, 2. , 7.14142843]])
请注意，上面的答案与之前的实现完全相同。同样，这里的优点是我们不需要为广播创建中间 2x3x4 数组。
为了完整性，让我们再次检查 threeSums 中每个被加数的维度是否允许广播。
np.sum(np.square(A)[:,np.newaxis,:], axis=2) has dimensions 2 x 1
2 * A.dot(B.T) has dimensions 2 x 3
np.sum(np.square(B), axis=1) has dimensions 1 x 3
因此，正如预期的那样，最终的 dist 数组的尺寸为 2x3。
中也讨论了使用点积代替逐元素乘法之和的方法
使用 numpy.linalg.norm 也适用于广播。为 axis 指定整数值将使用向量范数，默认为欧几里德范数。
import numpy as np
a = np.array([[1,1,1,1],[2,2,2,2]])
b = np.array([[1,2,3,4],[1,1,1,1],[1,2,1,9]])
np.linalg.norm(a[:, np.newaxis] - b, axis = 2)
# array([[ 3.74165739, 0. , 8.06225775],
# [ 2.44948974, 2. , 7.14142843]])
我有 2 x 4 和 3 x 4 的矩阵。我想找到跨行的欧氏距离，并在最后得到一个 2 x 3 的矩阵。下面是带有一个 for 循环的代码，它针对所有 b 行向量计算 a 中每个行向量的欧氏距离。如何在不使用 for 循环的情况下做同样的事情？
import numpy as np
a = np.array([[1,1,1,1],[2,2,2,2]])
b = np.array([[1,2,3,4],[1,1,1,1],[1,2,1,9]])
dists = np.zeros((2, 3))
for i in range(2):
dists[i] = np.sqrt(np.sum(np.square(a[i] - b), axis=1))
只需在正确的地方使用np.newaxis：
np.sqrt((np.square(a[:,np.newaxis]-b).sum(axis=2)))
此功能已包含在 scipy's spatial module 中，我建议使用它，因为它将在后台进行矢量化和高度优化。但是，正如另一个答案所表明的那样，您可以通过多种方式自己做到这一点。
import numpy as np
a = np.array([[1,1,1,1],[2,2,2,2]])
b = np.array([[1,2,3,4],[1,1,1,1],[1,2,1,9]])
np.sqrt((np.square(a[:,np.newaxis]-b).sum(axis=2)))
# array([[ 3.74165739, 0. , 8.06225775],
# [ 2.44948974, 2. , 7.14142843]])
from scipy.spatial.distance import cdist
cdist(a,b)
# array([[ 3.74165739, 0. , 8.06225775],
# [ 2.44948974, 2. , 7.14142843]])
我最近在使用深度学习（stanford cs231n，Assignment1）时遇到了同样的问题，但是当我使用
np.sqrt((np.square(a[:,np.newaxis]-b).sum(axis=2)))
出现错误
MemoryError
这意味着我 运行 内存不足（事实上，在 middle.It 中产生了一个 500*5000*1024 的数组，这么大！）
为防止该错误，我们可以使用一个公式来简化：
代码：
import numpy as np
aSumSquare = np.sum(np.square(a),axis=1);
bSumSquare = np.sum(np.square(b),axis=1);
mul = np.dot(a,b.T);
dists = np.sqrt(aSumSquare[:,np.newaxis]+bSumSquare-2*mul)
这里是原始输入变量：
A = np.array([[1,1,1,1],[2,2,2,2]])
B = np.array([[1,2,3,4],[1,1,1,1],[1,2,1,9]])
A
# array([[1, 1, 1, 1],
# [2, 2, 2, 2]])
B
# array([[1, 2, 3, 4],
# [1, 1, 1, 1],
# [1, 2, 1, 9]])
A 是一个 2x4 数组。 B 是一个 3x4 数组。
我们想在一个完全向量化的运算中计算欧氏距离矩阵运算，其中 dist[i,j] 包含 A 中第 i 个实例和 B 中第 j 个实例之间的距离。因此 dist 是 2x3这个例子。
距离
表面上可以用 numpy 编写为
dist = np.sqrt(np.sum(np.square(A-B))) # DOES NOT WORK
# Traceback (most recent call last):
# File "<stdin>", line 1, in <module>
# ValueError: operands could not be broadcast together with shapes (2,4) (3,4)
但是，如上所示，问题在于元素减法运算A-B涉及不兼容的数组大小，特别是第一维中的2和3。
A has dimensions 2 x 4
B has dimensions 3 x 4
为了进行逐元素减法，我们必须填充 A 或 B 以满足 numpy 的广播规则。我将选择用额外的维度填充 A，使其变为 2 x 1 x 4，这允许数组的维度对齐以进行广播。有关 numpy 广播的更多信息，请参阅 tutorial in the scipy manual and the final example in this tutorial.
您可以使用 np.newaxis 值或 np.reshape 命令执行填充。我在下面显示：
# First approach is to add the extra dimension to A with np.newaxis
A[:,np.newaxis,:] has dimensions 2 x 1 x 4
B has dimensions 3 x 4
# Second approach is to reshape A with np.reshape
np.reshape(A, (2,1,4)) has dimensions 2 x 1 x 4
B has dimensions 3 x 4
如您所见，使用任何一种方法都可以使尺寸对齐。我将使用 np.newaxis 的第一种方法。所以现在，这将用于创建 A-B，这是一个 2x3x4 数组：
diff = A[:,np.newaxis,:] - B
# Alternative approach:
# diff = np.reshape(A, (2,1,4)) - B
diff.shape
# (2, 3, 4)
现在我们可以将差分表达式代入dist等式语句得到最终结果：
dist = np.sqrt(np.sum(np.square(A[:,np.newaxis,:] - B), axis=2))
dist
# array([[ 3.74165739, 0. , 8.06225775],
# [ 2.44948974, 2. , 7.14142843]])
注意 sum 超过了 axis=2，这意味着对 2x3x4 数组的第三个轴（其中轴 id 从 0 开始）求和。
如果你的数组很小，那么上面的命令就可以正常工作。然而，如果你有大数组，那么你可能 运行 进入内存问题。请注意，在上面的示例中，numpy 在内部创建了一个 2x3x4 数组来执行广播。如果我们将 A 的维度概括为 a x z，将 B 的维度概括为 b x z，那么 numpy 将在内部创建一个 a x b x z 数组用于广播。
我们可以通过一些数学运算来避免创建这个中间数组。因为您将欧几里德距离计算为平方差和，所以我们可以利用平方差和可以重写的数学事实。
请注意，中间项涉及 逐元素 乘法的总和。这个乘法求和被称为点积。因为A和B各自都是一个矩阵，那么这个操作其实就是一个矩阵乘法。因此，我们可以将上面的内容重写为：
然后我们可以编写如下的numpy代码：
threeSums = np.sum(np.square(A)[:,np.newaxis,:], axis=2) - 2 * A.dot(B.T) + np.sum(np.square(B), axis=1)
dist = np.sqrt(threeSums)
dist
# array([[ 3.74165739, 0. , 8.06225775],
# [ 2.44948974, 2. , 7.14142843]])
请注意，上面的答案与之前的实现完全相同。同样，这里的优点是我们不需要为广播创建中间 2x3x4 数组。
为了完整性，让我们再次检查 threeSums 中每个被加数的维度是否允许广播。
np.sum(np.square(A)[:,np.newaxis,:], axis=2) has dimensions 2 x 1
2 * A.dot(B.T) has dimensions 2 x 3
np.sum(np.square(B), axis=1) has dimensions 1 x 3
因此，正如预期的那样，最终的 dist 数组的尺寸为 2x3。
使用 numpy.linalg.norm 也适用于广播。为 axis 指定整数值将使用向量范数，默认为欧几里德范数。
import numpy as np
a = np.array([[1,1,1,1],[2,2,2,2]])
b = np.array([[1,2,3,4],[1,1,1,1],[1,2,1,9]])
np.linalg.norm(a[:, np.newaxis] - b, axis = 2)
# array([[ 3.74165739, 0. , 8.06225775],
# [ 2.44948974, 2. , 7.14142843]])