在NumPy中计算行-行叉积的欧氏距离？

Question

我有两个包含n个坐标的numpy数组（每行包含两个项目）。

coords_a = np.random.random((20, 2))
coords_b = np.random.random((20, 2))

现在，对于每一行的组合，我想要计算一个函数，并将返回值保存为矩阵中的项。因此，结果数组应该具有形状(20, 20)，并且可以像下面展示的那样“延迟”计算。作为示例函数，使用欧几里得距离。

def euclidean_dist(x1: float, y1: float, x2: float, y2: float) -> float:
    """Return the euclidean distance between two the points (x1, y1) and (x2, y2)."""
    return np.sqrt(np.square(x1 - x2) + np.square(y1 - y2))

matrix = []
for a in coords_a:
    row = []
    for b in coords_b:
        row.append(euclidean_dist(*a, *b))
    matrix.append(row)
    
matrix = np.array(matrix)

正如您所想象的那样，这个嵌套的for循环非常耗时，仅使用2000个坐标对就需要超过25秒钟。 是否有建议的向量化方法来处理这种叉积运算？

提前感谢。

- BBQuercus

---

为什么需要计算变得懒惰？您能详细说明一下实际问题吗？您能够进行每列计算吗？您的数组有多大？ - Ehsan

“懒惰”是指下面的示例不是矢量化的。我的问题是，这个计算必须在几千组数据上执行，每组数据包含2000到5000个坐标。 - BBQuercus

一个5000x5000（X4floats）的计算应该很容易使用向量化在几乎任何系统中进行计算。那么你希望向量化的确切函数是什么？这取决于函数。例如，对于欧几里得距离，有内置函数可用。 - Ehsan

很抱歉没有理解。我将使用欧几里得距离（两个坐标之间的数学距离） - 应用于所有行的组合。 - BBQuercus

我认为向量化实现取决于您需要应用的函数。 - FBruzzesi

对于在NumPy中作为ufunc可用的函数，我敢打赌，在ufunc上进行外部乘积运算会快得多。但任何任意的函数都不一定能以这种方式使用。 - Curt F.

Answer 1

我想补充一下我的意见，因为并不是每个函数都已经在numpy或scipy中实现了。一般而言，您可以使用numpy广播来实现向量化解决方案。 对于欧几里得距离的特定情况，以下是如何完成它：

import numpy as np

# Define the arrays of coordinates
coords_a = np.random.random((20, 2))
coords_b = np.random.random((20, 2))

# Expand their dimensions
a = coords_a[:, None]
b = coords_b[None, None]

# Use broadcasting to compute pairwise difference
d = a-b

# Apply formula for euclidean distance
r = np.sqrt(np.sum(d**2, axis=-1)) 

就时间性能而言，scipy.spatial.distance.cdist在这种特定情况下速度更快，但并非所有函数都可用：

import numpy as np
from scipy.spatial.distance import cdist

a = np.random.random((10_000, 2))
b = np.random.random((10_000, 2))

euc_broadcast = lambda a,b: np.sqrt(np.sum(np.square(a[:, None]-b[None, :]), axis=-1))

%timeit euc_broadcast(a, b)
3.39 s ± 149 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)

%timeit cdist(a, b)
603 ms ± 13.2 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)

- FBruzzesi

---

这是一个非常有帮助的答案。我正需要这个，但找了一段时间才找到！ - seeker_after_truth

Answer 2

针对您的具体示例，您可以执行以下操作：

from scipy.spatial.distance import cdist
cdist(coords_b,coords_a)

一般而言，向量化取决于您的函数。

- Ehsan