如何无循环计算两个矩阵数据集的Euclidean distance matrix？

当然可行，而且这种基于广播的实现效率远高于Python循环，尤其适合处理大规模数据。核心是通过欧式距离的代数变形结合numpy的广播机制来完成计算。

原理推导 #

对于两个d维向量$\boldsymbol{x}$（来自matrix1，维度为$(d,)$）和$\boldsymbol{y}$（来自matrix2，维度为$(d,)$），欧式距离的平方可展开为：
$$
|\boldsymbol{x} - \boldsymbol{y}|^2 = |\boldsymbol{x}|^2 + |\boldsymbol{y}|^2 - 2\boldsymbol{x} \cdot \boldsymbol{y}
$$
将这个公式扩展到矩阵层面：

1. 计算matrix1中每个向量的模长平方，得到形状为$(n_1, 1)$的数组
2. 计算matrix2中每个向量的模长平方，转置后得到形状为$(1, n_2)$的数组
3. 计算两个矩阵的点积，得到形状为$(n_1, n_2)$的数组
4. 利用广播将前两个结果相加，再减去2倍的点积，最后开根号得到最终的欧式距离矩阵

代码实现示例 #

基于numpy的具体实现如下：

import numpy as np

# 示例数据（替换为你的实际矩阵）
n1, n2, d = 50, 30, 10
matrix1 = np.random.rand(n1, d)
matrix2 = np.random.rand(n2, d)

# 计算各向量的模长平方
norm1 = np.sum(matrix1 ** 2, axis=1, keepdims=True)  # shape (n1, 1)
norm2 = np.sum(matrix2 ** 2, axis=1, keepdims=True).T  # shape (1, n2)

# 计算距离平方矩阵，修正数值精度问题后开根号
dist_sq = norm1 + norm2 - 2 * np.dot(matrix1, matrix2.T)
dist_sq = np.maximum(dist_sq, 0.0)  # 避免浮点误差导致的极小负数
dist_matrix = np.sqrt(dist_sq)  # 最终shape为(n1, n2)

关键说明 #

• 广播机制的作用：当norm1（$(n1,1)$）和norm2（$(1,n2)$）相加时，numpy会自动将二者扩展为$(n1,n2)$的形状，实现逐元素匹配相加。
• 数值精度修正：浮点运算可能导致原本为0的距离平方出现极小负数，开根号会产生NaN，因此用np.maximum将其修正为0。
• 效率优势：这种方法的时间复杂度和循环一致（$O(n_1 n_2 d)$），但numpy底层基于C实现，执行速度比Python循环快几个数量级。

内容的提问来源于stack exchange，提问作者guavcxn