0.引言

论文参考

A survey on trajectory clustering analysis	Trajectory Clustering: A Partition-and-Group Framework
Distribution-Based Trajectory Clustering	A principled distributional approach to trajectory similarity measurement
A survey of trajectory distance measures and performance evaluation

博客参考

轨迹相似度计算方法汇总	JUST技术：如何通过轨迹相似性度量方法，发现新冠易感人群
TraClus轨迹聚类算法原理及java版实现	Trajectory Clustering（DBSCAN算法进行轨迹聚类）
论文笔记： Trajectory Clustering: A Partition-and-Group Framework	轨迹聚类之分段聚类方法总结
轨迹路线相似度计算	衡量轨迹相似度算法

轨迹作为一种时空数据，指的是某物体在空间中的移动路径，通常表示为GPS点的序列，例如tr=<p1→p2→…pn>，其中点pi=(lat,lng,t)，表示该物体在t时刻位于地理坐标位置(lat,lng)上，lat和lng分别表示纬度和经度。

轨迹数据的分析处理非常具有挑战性，主要包含三个方面：1）轨迹数据量大；2）轨迹数据噪音多；3）轨迹数据获取途径多样。其中，轨迹相似性作为一项基础算法服务，衡量两条轨迹之间的距离大小，可为其上层应用提供支持，也是目前研究的热点之一。

相对于点与点或点与轨迹之间的距离度量，轨迹之间的距离度量更加的复杂，需要考虑的因素也更多，例如轨迹的采样率、考虑轨迹的时间信息和轨迹自身的噪音等。常见的轨迹相似性度量方法大致分类如下图所示。

$D_{owd}(T_1,T_2)=\frac{1}{|T_1|}({\sum }_{p\in T_1}{D}_{point}(p,T_2))$
判断两条轨迹的相似性方法:

• 基于点方法： EDR，LCSS，DTW等
• 基于形状的方法： Frechet， Hausdorff
• 基于分段的方法：One Way Distance, LIP distance
• 基于特定任务的方法：TRACLUS， Road Network，grid等

我们定义如下两条轨迹，长度分别为n和m，则：

1.欧式距离（Euclidean Distance）

欧式距离要求两条轨迹的长度相同一一对应，其数学定义为：

欧氏距离的定义简单明了，就是两条轨迹对应点的空间距离的平均值，但是缺点也很明显，就是不能度量不同长度的轨迹相似性，且对噪音点敏感。

.欧式距离示意图（黑点与红点一一对应）

2.动态时间归整（Dynamic Time Warping，DTW）

如上所述，欧式距离的一个明显的限制是要求两条轨迹长度相等，这在实际情况中是不太可能的，更理想的情况应该在轨迹长度上具有一定的灵活性。
动态时间归整DTW的思想是自动扭曲两个序列，并在时间轴上进行局部的缩放对齐，以使其形态尽可能一致，从而得到最大可能的相似性。DTW 将两条轨迹的点进行多对多的映射，从而较为高效地解决了数据不齐的问题，其动态规划算法如下：

其中Head(tr)=<p1>表示该轨迹的第一个点；Rest(tr)=<p2…pn>表示除第一个点之外的所有点组成的子序列。

DTW示意图（黑点与红点进行多对多映射）

动态时间归整算法灵活，对轨迹长度无限制，且效果较好，但是其并未对噪音点进行处理，离群点也会对结果造成较大影响。

• bilibili
  

3.最长公共字串（Longest Common Sub-Sequence，LCSS）

有一个经典的算法问题：求解两序列的最长公共子序列，不要求公共子序列中的两个连续相连，例如BDCABA和ABCBDAB的最大公共子序列为BCBA。在此基础上，很自然提出了基于最长公共子序列的轨迹相似性度量方法，即LCSS，其值代表最多可被视为同一点的点数，也就是两条轨迹中满足最小距离阈值限制的轨迹点的对数。其基于动态规划的算法如下：

其中，参数是最小距离阈值，两点之间距离小于该值时将被认为是同一点，此外，该算法对轨迹长度没有限制。

LCSS示意图，有两对点的距离小于 e，则被视为同一点，LCSS值为2

最长公共子串距离对噪音点进行了处理，即因噪音点的偏离没有与其相近的轨迹点故不会被计算在最终结果内，这一步骤有效对抗噪音。但与此同时，该算法的最小距离阈值不好定义，还有可能返回并不相似的轨迹。

1.确定dp数组（dp table）以及下标的含义
dp[i][j]：长度为[0, i - 1]的字符串text1与长度为[0, j - 1]的字符串text2的最长公共子序列为dp[i][j]

2.确定递推公式

• 如果text1[i - 1] 与 text2[j - 1]相同，那么找到了一个公共元素，所以dp[i][j] = dp[i - 1][j - 1] + 1;

• 如果text1[i - 1] 与 text2[j - 1]不相同，那就看看text1[0, i - 2]与text2[0, j - 1]的最长公共子序列 和 text1[0, i - 1]与text2[0, j - 2]的最长公共子序列，取最大的。即：dp[i][j] = max(dp[i - 1][j], dp[i][j - 1]);

  if (text1[i - 1] == text2[j - 1]) {
      dp[i][j] = dp[i - 1][j - 1] + 1;
  } else {
      dp[i][j] = max(dp[i - 1][j], dp[i][j - 1]);
  }
  

4.编辑距离（Edit Distance on Real sequence，EDR）

给定两个长度分别为n和m的轨迹tr1和tr2，最小距离的匹配阈值，则两条轨迹之间的EDR距离就是需要对tr1及逆行插入、删除或替换使其变为tr2的操作次数，其基于动态规划的算法如下。

其中，

EDR示意图，在p1处“插入”一点、将p2’“替换”为p3和在p5处“插入”一点，共计3个操作使两条轨迹相等（即对应点距离均小于阈值），故其EDR值为3， 值越小越相似。

轨迹的编辑距离为轨迹相似新度量提供了一种新的思路，其缺陷也很明显，就是对噪音点敏感。

1.确定dp数组（dp table）以及下标的含义
dp[i][j] 表示以下标i-1为结尾的字符串word1，和以下标j-1为结尾的字符串word2，最近编辑距离为dp[i][j]。
2.确定递推公式
在确定递推公式的时候，首先要考虑清楚编辑的几种操作，整理如下：

if (word1[i - 1] == word2[j - 1])
    不操作
if (word1[i - 1] != word2[j - 1])
    增
    删
    换

5.豪斯多夫距离（Hausdorff Distance）

简单来说，豪斯多夫距离就是两条轨迹最近点距离的最大值。度量了两个点集间的最大不匹配程度

其中，h(tr1,tr2)称为tr1到tr2的单向豪斯多夫距离，其定义如下：

豪斯多夫距离示意图，每个红点都有一个最近的黑点，豪斯多夫距离就是其中的最大值。

• Hausdorff distance 参考

缺点：聚类阈值不太好选取。

6.弗雷歇距离（Fréchet Distance）

直观的理解，弗雷歇距离就是狗绳距离，即主人走路径A，狗走路径B，各自走完两条路径过程中所需要的最短狗绳长度。实际的是计算两个序列之间的最大差值，是距离的最大值。

狗绳距离示意图，虚线表示主人和狗在同一时刻所处位置的对应，弗雷歇距离即为长度最长的虚线

弗雷歇距离基于动态规划思想的算法如下：

其中，d(p,q)是两个GPS点的欧式距离，tr(n-1)=<p1→p2→…pn-1> 是轨迹tr的长度为n-1的子轨迹。
弗雷歇距离为我们提供了一种简单直观的度量相似性的方式，也能达到较好的效果；但可惜的是其并没有对噪音点进行处理，例如若狗的某个轨迹点因为噪声偏离得很远，那么弗雷歇距离也随之增大，这显然是不合理的。

7.单向距离（One Way Distance，OWD）

单向距离OWD[7]的定义如下：

其中，|tr1|表示轨迹tr1的长度，d(p,tr2)表示GPS点p到tr2的距离。为了对称，简单修改上述公式：

单向距离示意图，该距离即为各多边形的面积之和与折线长度的比值

OWD距离的基本思想基于两条轨迹围成的面积，当面积大，说明轨迹之间距离较远，相似度就低；相反，若围成的面积为0，则说明两条轨迹重合，相似度最高。

8.多线位置距离（Locality In-between Polylines，LIP）

多线位置距离LIP定义如下：

其中，$I_i$表示两条轨迹的第i个交点，权重$w_i$定义如下：

LIP方法示意图，该距离为每个区域面积与其权重乘积之和

LIP方法很好理解，当某区域面积的周长占总长比重大时权重也自然就大；

• 当Area均为0时，说明两条轨迹重合没有缝隙，LIP距离为0；
• 当Area加权和大时，则说明两条轨迹之间缝隙较大，LIP距离也就大。
• 此外，权重由区域周长占总长比重大决定，也一定程度对抗了噪音点的干扰。

9.巴氏距离（Bhattacharyya Distance）

在统计中, Bhattacharyya距离测量两个离散或连续概率分布的相似性。它与衡量两个统计样品或种群之间的重叠量的Bhattacharyya系数密切相关。 Bhattacharyya距离和Bhattacharyya系数以20世纪30年代曾在印度统计研究所工作的一个统计学家A. Bhattacharya命名。同时, Bhattacharyya系数可以被用来确定两个样本被认为相对接近的, 它是用来测量中的类分类的可分离性。

对于离散概率分布 $p$ 和 $q$ 在同一域 $X$, 它被定义为:
$$D_B(p,q)=-\ln (BC(p,q))$$

其中：
$$BC(p,q)=\sum _{x\in X}\sqrt{p(x)q(x)}$$

是Bhattacharyya系数。
对于连续概率分布, Bhattacharyya系数被定义为:
$$BC(p,q)=\int \sqrt{p(x)q(x)}dx$$

Bhattacharyya系数是两个统计样本之间的重叠量的近似测量, 可以被用于确定被考虑的两个样本的相对接近。
计算Bhattacharyya系数涉及集成的基本形式的两个样本的重叠的时间间隔的值的两个样本被分裂成一个选定的分区数, 并且在每个分区中的每个样品的成员的数量, 在下面的公式中使用
$$\text{ Bhattacharyya }=\sum _{i=1}^n\sqrt{\left(\sum a_i\cdot \sum b_i\right)}$$

考虑样品 $a$ 和 $b,n$ 是的分区数, $\sum a_i$ 是指样品 $a$ 中落在分区 $i$ 内的个数, $\sum b_i$ 有类似的定义。

10.Levenshtein Distance（编辑距离）

Levenshtein Distance，一般称为编辑距离（Edit Distance，Levenshtein Distance只是编辑距离的其中一种）或者莱文斯坦距离，算法概念是俄罗斯科学家弗拉基米尔·莱文斯坦（Levenshtein · Vladimir I）在1965年提出。此算法的概念很简单：Levenshtein Distance指两个字串之间，由一个转换成另一个所需的最少编辑操作次数，允许的编辑操作包括：

• 将其中一个字符替换成另一个字符（Substitutions）。
• 插入一个字符（Insertions）。
• 删除一个字符（Deletions）。

• 参考

11. 细胞相似性(CSIM)

CSIM是一种仅考虑两条轨迹所经过区域的相似性度量方案。CSIM使用一个正方形（ $3\times 3$ ）网格结构进行膨胀，然后利用类似Jaccard 系数的方式测量网格单元总数中有多少个共同的网格单元。该算法的主要优点是具有线性时间复杂度, 且结果不受点偏移量的影响。
$$\mathrm{Grid}\left(C_A,C_B\right)=\frac{\left|C_A\cap C_B\right|+\left|C_A\cap C_B^d\right|+\left|C_B\cap C_A^d\right|}{\left|C_A\right|+\left|C_B\right|-\left|C_A\cap C_B\right|}$$

CSIM通过将轨迹点映射到规则的⽹格单元上来提取代表性轨迹。不同的单元格⼤⼩对轨迹提取结果的准确性有显著影响。如果单元尺⼨太⼤，轨迹会被过度压缩，关键特征会 变得模糊。相反，过⼩的⽹格单元尺⼨会增加代表性轨迹中的单元数量并增加计算负荷。需要根据用户历史Log数据来分析寻找最佳的细胞⼤⼩。

12.距离归并算法（Merge distance）

这个算法的核心是把两个轨迹合并成一个中间轨迹，然后计算这个中间轨迹和两条原轨迹的差异作为两条轨迹的距离。

对于左边的图，下面的那个路径就是对上面的路径进行归并后的最短路径；对于右边的图，绿色的实线路径就是对蓝色和红色的虚线路径归并后的最短路径。

Merge Distance算法就是为了求出这个最短的路径长度。

$$\{\begin{array}{ll}\mathrm{Length}\left(T_1(1,i),T_2(1,j)\right)=\min (\mathrm{Length}\left(T_1(1,i-1),T_2(1,j)\right)+d\left(T_{1,i-1},T_{1,i}\right),\mathrm{Length}\left(T_2(1,j),T_1(1,i-1)\right)+d\left(T_{2,j},T_{1,i}\right)),& \text{ if }\le i\le m\text{ and }1\le j\le n\\ \mathrm{Length}\left(T_2(1,j),T_1(1,i)\right)=\min (\mathrm{Length}\left(T_1(1,i),T_2(1,j-1)\right)+d\left(T_{1,i},T_{2,j}\right),\mathrm{Length}\left(T_2(1,j-1),T_1(1,i)\right)+d\left(T_{2,j-1},T_{2,j}\right)),& \text{ if }\le i\le m\text{ and }2\le j\le n\\ \mathrm{Length}\left(T_1(1,i),T_2(1,j)\right)=\left({\sum }_{k=1}^{j-1}d\left(T_{2,k},T_{2,k+1}\right)\right)+d\left(T_{2,j},T_{1,1}\right),& \text{ if }i=1\\ \mathrm{Length}\left(T_1(1,i),T_2(1,j)\right)=\left({\sum }_{k=1}^{i-1}d\left(T_{1,k},T_{1,k+1}\right)\right)+d\left(T_{1,i},T_{2,1}\right),& \text{ if }j=1\\ \mathrm{Length}\left(T_{\text{super }}\right)=\min \left(\mathrm{Length}\left(T_1(1,m),T_2(1,n)\right),\mathrm{Length}\left(T_2(1,n),T_1(1,m)\right)\right)& \end{array}$$

$d_{MD}(T_1,T_2)=\frac{Length(T_{super})}{Length(T_1)+Length(T_2)}-1$

• 离散连续均可
• 所有采样点匹配
• 不具有可度量性
• 仅考虑位置信息
• 不需要手动设定参数
• 时间复杂度为$O(mn)$
• 对于误差有一定容忍性