距离矩阵 | 会飞的贼

万壑松风知客来，摇扇抚琴待留声

1. 文起

好久没有写文章了，这段时间其实累积了很多想分享和记录的事，很想每一件都写下来，不过由于一篇文章需要仔细筹备相关的知识所以时间抽不出来，不过我会慢慢补上。

先从哪篇文章开始呢，我决定从最近刚遇到的一个问题开始做分享——距离矩阵，所以这是一篇讲通过矩阵来优化计算机运算的知识分享。

距离矩阵能够帮助你减少 for 循环的使用从而达到计算两个矩阵所有元素之间的距离。尤其是 Python 中的一个科学计算包 Numpy，这种矩阵运算简直就是为它准备的。

2. 问题描述

机器学习的 KNN 算法实现原理：计算测试集每一个向量元素与训练集的每一个向量元素的距离，从中选择 K 个与之相近的样本，投票选择这些样本的目标值，从而得到预测结果。

这个距离一般采用欧氏距离（勾股定理计算方式）。如果不采用一定优化措施，其实 KNN 是一个很耗费时间的算法，通过 for 循环测试集与训练集的每个样本之间都需要计算一次。当数据集达到一定数量后运算就会十分的缓慢，由此我们需要一个方法来取代循环，使用更高效的方式达到这个距离的运算。

当两个过大的数据列表计算时，我们很容易联想到使用矩阵运算来代替。矩阵运算能够十分高效的处理这种由多个向量组成的数据。所以现在已经有了尝试解决的方案，下面来以具体数据来展开思考。

3. 距离矩阵

3.1. 场景

（起点随机）两个人在一个棋盘上随意跳动，他们每一步的行走位置由横纵坐标表示。两个人跳一段时间后可以得到两个列表（长度可能不同），列表中的每个元素是一个长度为 2 的元组或列表，由此每个元素表示一个二维坐标信息，两个列表分别记录了两个人跳的每一步，长度则是代表了分别跳的步数。

3.2. 问题

现在需要计算这两组步数两两之间的欧几里得距离，从中找到这两个人相距最近的位置、最远的位置，或者其它按需求寻找的位置信息。

3.3. 解决问题

以上是一个比较简单的场景（这只是一个二维空间），一个位置信息仅需要两个数字就能表示。不过我们再简化一下这个问题，两个人跳的步长并不远都是三步。以列表形式表现为：A: [[1,2],[3,4],[5,6]]、B:[[2,2],[3,3],[4,4]]。

根据上文的思路我们决定以矩阵的方式来解决该问题，所以首先就得将两个列表转换为矩阵的形式。每行样本表示一个位置坐标，由于是二维空间，所以每行只有两列（如果在 N 维空间，则每行有 N 列）

import numpy as np

A = np.array([[1,2],[3,4],[5,6]])
B = np.array([[2,2],[3,3],[4,4]])
A = A.reshape((-1,2))
B = B.reshape((-1,2))
# A:
# [[1 2]
#  [3 4]
#  [5 6]]
# B:
# [[2 2]
#  [3 3]
#  [4 4]]

每行表示一个坐标信息，两两之间求距离，则共有 9 种组合情况，即可以表示成一个 3X3 的矩阵。矩阵的 0、1、2 分别代表了 A、B 矩阵的行索引，i 表示 A 矩阵的任意一行，j 表示 B 矩阵的任意一行。

从单元格入手：

那么首先从最简单的一个单元格入手寻找规律。以该黄色单元格为例表示 A 矩阵的第 0 行（[1,2]）与 B 矩阵的第 1 行（[3,3]）求欧式距离：

$dist(A_{0}\, B_{1}) = \sqrt{(A_{00}-B_{10})^{2} + (A_{01}-B_{11})^{2}}$ $=$ $\sqrt{(A_{00}^{2} + A_{01}^{2}) + (B_{10}^{2} + B_{11}^{2}) - 2\times (A_{00}B_{10} + A_{01}B_{11})}$

$A_{00}^{2} + A_{01}^{2}$ 其实表示的是 A 矩阵第 0 行的所有元素求平方和（这里是二维，所以只有两个数[1,2]）。这种方式可以通过数学符号 $||$ (模) 来表示 (hexo 渲染公式有一点点问题)。$A_{00}B_{10} + A_{01}B_{11}$ 表示两个二维向量的对应位置相乘，可以通过矩阵相乘的方式表达。所以距离矩阵的第 0 行第 1 列化简可以表示为：

$dist(A_{0}\, B_{1}) = \sqrt{\left || A_{0} \right ||^{2} + \left || B_{1} \right ||^{2} - 2\times A_{0}B_{1}^T} \Rightarrow \sqrt{\left || A_{i} \right ||^{2} + \left || B_{j} \right ||^{2} - 2\times A_{i}B_{j}^T}$ （i、j 分别代表 A、B 的任意一行）

扩展到某一行：

上面已经求得了距离矩阵 [0,1] 位置处的距离表达公式，该公式已经可以扩展到距离矩阵的任何一个单元格，只需要根据所求行、列的位置，改变 A、B 的下标即可。

这是一个规律，但还不够，我们做一点改变将这个公式从一个单元格表达式扩展到一行表达式，计算 A 矩阵的第 0 个向量与 B 矩阵的所有向量的距离，也就是距离矩阵的第 0 行：

简单解释下面公式扩展的一个理论知识，否则你可能会有疑惑。

1: 距离矩阵大小等于 A.shape[0] x B.shape[0]，所以一行有 B.shape[0] 个元素，这里 B.shape[0] = 3，扩展时需要保证这一行有 3 个元素。

2: 矩阵相加时，是对应位置一一相加。这两点是后面统一公式的条件。

简单解释后继续我们的分析。下面分别列出了这三个位置对应的计算公式。

$\sqrt{\left | A_{0} \right |^{2} + \left | B_{0} \right |^{2} - 2\times A_{0}B_{0}^T}$ $\sqrt{\left || A_{0} \right ||^{2} + \left || B_{1} \right ||^{2} - 2\times A_{0}B_{1}^T}$ $\sqrt{\left | A_{0} \right |^{2} + \left | B_{2} \right |^{2} - 2\times A_{0}B_{2}^T}$

接着统一这一行的元素公式（也就是把单独计算的每个元素放到一个公式中，还记得上面所说的矩阵对应位置相加的原理吗）：

$dist(A_{0}\, B) = \sqrt{(\left || A_{0} \right ||^{2} \: \left || A_{0} \right ||^{2} \: \left || A_{0} \right ||^{2}) + (\left || B_{0} \right ||^{2} \: \left || B_{1} \right ||^{2} \: \left || B_{2} \right ||^{2}) - 2\times A_{0}\left ( B_{0}^{T} \: B_{1}^{T} \: B_{2}^{T} \right )}$

到此我们得到了计算距离矩阵一行的公式，如果要表示任意一行的话，稍微改变上面公式的元素下标即可（j=2 ，如果多维则 j 表示最高维数）：

$dist(A_{i}\, B) = \sqrt{(\left || A_{i} \right ||^{2} \: \left || A_{i} \right ||^{2} \: \cdots \: \left || A_{i} \right ||^{2}) + (\left || B_{0} \right ||^{2} \: \left || B_{1} \right ||^{2} \: \cdots \: \left || B_{j} \right ||^{2}) - 2\times A_{i}B^{T}}$

扩展到整个距离矩阵：

上面得到了某一行的距离矩阵公式，接下来要做的就是将行公式按列扩展，从而得到整个距离矩阵的表达式（将 i 从 0 表达到 A.shape[0]，这里也就是 3）。

同样以 A、B 矩阵为例，求这个 (3,3) 的距离矩阵的表达式如下：

$$
dist(A \, B) = \sqrt{\begin{pmatrix}
\left || A_{0} \right ||^{2} & \left || A_{0} \right ||^{2} & \left || A_{0} \right ||^{2}\\
\left || A_{1} \right ||^{2} & \left || A_{1} \right ||^{2} & \left || A_{1} \right ||^{2}\\
\left || A_{2} \right ||^{2} & \left || A_{2} \right ||^{2} & \left || A_{2} \right ||^{2}
\end{pmatrix}
+
\begin{pmatrix}
\left || B_{0} \right ||^{2} & \left || B_{1} \right ||^{2} & \left || B_{2} \right ||^{2}\\
\left || B_{0} \right ||^{2} & \left || B_{1} \right ||^{2} & \left || B_{2} \right ||^{2}\\
\left || B_{0} \right ||^{2} & \left || B_{1} \right ||^{2} & \left || B_{2} \right ||^{2}
\end{pmatrix}
-
2\times AB^{T}
}
$$
其中 $AB^{T}$ 得到的也是像前面两个 (3,3) 大小的矩阵。照应前面的说明，3 个矩阵对应位置相运算，最后开根得到距离矩阵。

上面距离矩阵的 i、j 仅表示三个维度，同样我们可以以 i、j 表示任何维度大小，从而得到任意大小两个矩阵的距离矩阵表达式。

$$
dist(A \, B) = \sqrt{\begin{pmatrix}
\left || A_{0} \right ||^{2} & \left || A_{0} \right ||^{2} & \cdots & \left || A_{0} \right ||^{2}\\
\left || A_{1} \right ||^{2} & \left || A_{1} \right ||^{2} & \cdots & \left || A_{1} \right ||^{2}\\
\vdots & \vdots & \ddots & \vdots \\
\left || A_{i} \right ||^{2} & \left || A_{i} \right ||^{2} & \cdots & \left || A_{i} \right ||^{2}
\end{pmatrix}
+
\begin{pmatrix}
\left || B_{0} \right ||^{2} & \left || B_{1} \right ||^{2} & \cdots & \left || B_{j} \right ||^{2}\\
\left || B_{0} \right ||^{2} & \left || B_{1} \right ||^{2} & \cdots & \left || B_{j} \right ||^{2}\\
\vdots & \vdots & \ddots & \vdots \\
\left || B_{0} \right ||^{2} & \left || B_{1} \right ||^{2} & \cdots & \left || B_{j} \right ||^{2}
\end{pmatrix}
-2 \times AB^{T}
}
$$

到此终于完成了任务，我们已经得到了关于两个矩阵（跳出这里 (3,2) 矩阵的思维，两个矩阵可以表示任意长度，任意维度）的距离矩阵表达式。

4. 实现公式

4.1. 手动实现公式

通过上面一步一步来，我们得到了最终的距离矩阵公式，不过使用代码如何实现呢？这里使用 Numpy 按照公式来对 A、B 矩阵进行运算：

A_matrix = np.transpose([np.sum(np.square(A), axis=1)])
B_matrix = np.sum(np.square(B), axis=1)
C_matrix = -2*np.dot(A,B.T)

D_matrix = np.sqrt(A_matrix + B_matrix + C_matrix)
print(D_matrix)
# [[1.         2.23606798 3.60555128]
#  [2.23606798 1.         1.        ]
#  [5.         3.60555128 2.23606798]]

答疑说明：A_matrix、B_matrix、C_matrix 与上面的公式一一对应。A_matrix、B_matrix 如何理解呢？这是简化了矩阵的表示方式可以看到第一个矩阵表达中每一行的元素都相等，第二个矩阵中每一列元素都相等。则可以只用一列、一行来代替，而在用它们和矩阵做运算时，会展开成对应大小的矩阵再对应位置运算。
1
2
3
4
> [[ 5]		[[5 , 5 , 5 ]
>  [25]	->	 [25, 25, 25]
>  [61]]	 [61, 61, 61]]
>

> [ 8 18 32] ->	[[8 , 18, 32]
> 		 [8 , 18, 32]
> 		 [8 , 18, 32]
>

4.2. 使用 scipy 公式

scipy 是 python 的一个数学计算包，包含了大量线性代数以及其它知识。后来我才发现原来 scipy 中早已包含了计算距离矩阵的函数，真是让人省心，我们使用它来验证上面的结果。

from scipy.spatial.distance import cdist

result = cdist(A, B, metric='euclidean')
print(result)
# [[1.         2.23606798 3.60555128]
#  [2.23606798 1.         1.        ]
#  [5.         3.60555128 2.23606798]]

验证结果完全没有问题，说明我们的公式完全正确。虽然有了现成的函数，但如果能掌握推算原理对我们理解问题将会有更进一步的帮助。

4.3. 使用距离矩阵

这里我们简单来使用一下距离矩阵，解决正如前面 3.2. 中提到的问题。

# 返回距离矩阵中的最大值、最小值
max_distance = D_matrix.max()
min_distance = D_matrix.min()

# 通过位置查询找到最大值、最小值在距离矩阵中对应的位置坐标
max_locate = np.where(D_matrix == max_distance)
min_locate = np.where(D_matrix == min_distance)
# (array([2], dtype=int64), array([0], dtype=int64))
# (array([0, 1, 1], dtype=int64), array([0, 1, 2], dtype=int64))

# 第一个元素对应A矩阵的行，第二个元素对应B矩阵的行
# 比如最大距离5，由A的2行和B的0行组成->[5 6]、[2 2]

可以利用距离矩阵做的事，远多于这个简单的例子。

5. 文末

关于向量的运算，很多时候都可以通过矩阵的方式来进行优化。上面的距离矩阵运算帮助我提高了程序 30 倍以上的运算效率，所以我觉得该知识点确实值得去学习。

放假抽了一个晚上记录完了这篇文章，还是很欣慰的。