前言
时间序列相似性度量是时间序列相似性检索、时间序列无监督聚类、时间序列分类以及其他时间序列分析的基础。给定时间序列的模式表示之后,需要给出一个有效度量来衡量两个时间序列的相似性。时间序列的相似性可以分为如下三种:
1、 时序相似性
- 时序相似性是指时间序列点的增减变化模式相同,即在同一时间点增加或者减少,两个时间序列呈现一定程度的相互平行。这个一般使用闵可夫斯基距离即可进行相似性度量。
2、 形状相似性
- 形状相似性是指时间序列中具有共同的形状,它通常包含在不同时间点发生的共同的趋势形状或者数据中独立于时间点相同的子模式。两个时间序列整体上使用闵可夫斯基距离刻画可能不相似,但是他们具有共同相似的模式子序列,相似的模式子序列可能出现在不同的时间点。这个一般使用DTW动态时间规整距离来进行相似性刻画。
3、变化相似性
- 变化相似性指的是时间序列从一个时间点到下一个时间点的变化规律相同,两个时间序列在形状上可能并不一致,但是可能来自于同一个模型。这个一般使用ARMA或者HMM等模型匹配方法进行评估。
时间序列相似性度量可能会受到如下因素影响:
- 时间序列作为真实世界的系统输出或者测量结果,一般会夹杂着不同程度的 噪声扰动 ;
- 时间序列一般会呈现各种变形,如
- 振幅平移
- 振幅压缩
- 时间轴伸缩
- 线性漂移
- 不连续点等
- 时间序列之间可能存在不同程度的关联;
以上因素在衡量时间序列相似性度量的时候要根据具体情况进行具体分析。
闵可夫斯基距离
给定两条时间序列:
P
=
(
x
1
,
x
2
,
.
.
.
x
n
)
,
Q
(
y
1
,
y
2
,
.
.
.
y
n
)
P=(x_1,x_2,...x_n),\ \ Q(y_1,y_2,...y_n)
P=(x1,x2,...xn), Q(y1,y2,...yn)
闵可夫斯基距离的定义如下:
d
i
s
t
(
P
,
Q
)
=
(
∑
i
=
1
n
∣
x
i
−
y
i
∣
p
)
1
p
dist(P,Q) = \left(\sum\limits_{i=1}^n|x_i-y_i|^p\right)^{\frac{1}{p}}
dist(P,Q)=(i=1∑n∣xi−yi∣p)p1
注:
-
当 p = 1 p=1 p=1时,闵可夫斯基距离又称为 曼哈顿距离 :
d i s t ( P , Q ) = ∑ i = 1 n ∣ x i − y i ∣ dist(P,Q)=\sum\limits_{i=1}^n |x_i-y_i| dist(P,Q)=i=1∑n∣xi−yi∣ -
当 p = 2 p=2 p=2时,闵可夫斯基距离又称为 欧氏距离 :
d i s t ( P , Q ) = ( ∑ i = 1 n ∣ x i − y i ∣ 2 ) 1 2 dist(P,Q) = \left(\sum\limits_{i=1}^n|x_i-y_i|^2\right)^{\frac{1}{2}} dist(P,Q)=(i=1∑n∣xi−yi∣2)21 -
当 p → ∞ p\rightarrow\infty p→∞时,闵可夫斯基距离又称为 切比雪夫距离 :
lim p → ∞ ( ∑ i = 1 n ∣ x i − y i ∣ p ) 1 p = max i ∣ x i − y i ∣ \lim\limits_{p\rightarrow\infty}\left(\sum\limits_{i=1}^n|x_i-y_i|^p\right)^{\frac{1}{p}} = \max\limits_{i}|x_i-y_i| p→∞lim(i=1∑n∣xi−yi∣p)p1=imax∣xi−yi∣ -
闵可夫斯基距离模型简单,运算速度快。
-
闵可夫斯基距离比较直观,但是它与数据的分布无关,具有一定的局限性
-
在时间轴伸缩和弯曲情况下无法给出确切的度量
-
只能针对长度相等的时间序列
-
如果时间序列 P P P 的幅值远远大于 Q Q Q 的幅值,那么得到的计算结果就会被时间序列 P P P 所过度作用。
-
计算之前需要对数据进行尺度变化,例如减去均值并除以标准差,或者极大极小归一化
-
-
这种方法在假设数据各个维度不相关的情况下,利用数据分布的特性计算出不同的距离,如果数据维度之间数据相关,这时该类距离就不合适了!
余弦相似度
给定两个时间序列:
P
=
(
x
1
,
x
2
,
.
.
.
x
n
)
,
Q
(
y
1
,
y
2
,
.
.
.
y
n
)
P=(x_1,x_2,...x_n),\ \ Q(y_1,y_2,...y_n)
P=(x1,x2,...xn), Q(y1,y2,...yn)
余弦相似度如下:
D
c
o
s
(
P
,
Q
)
=
∑
i
=
1
n
x
i
y
i
∑
i
=
1
n
x
i
2
∑
i
=
1
n
y
i
2
D_{cos}(P,Q) = \frac{\sum\limits_{i=1}^n x_iy_i}{\sqrt{\sum\limits_{i=1}^nx_i^2}\sqrt{\sum\limits_{i=1}^ny_i^2}}
Dcos(P,Q)=i=1∑nxi2i=1∑nyi2i=1∑nxiyi
动态时间规整距离
动态时间规整距离是一种动态规划方法,通过寻找最小的对齐匹配路径使两序列之间距离最小,因而DTW可计算不同长度序列的距离,对时间序列偏移不敏感。
给定两个时间序列(两个时间序列并一定等长,即m和n不一定相等):
P
=
(
x
1
,
x
2
,
.
.
.
x
m
)
,
Q
(
y
1
,
y
2
,
.
.
.
y
n
)
P=(x_1,x_2,...x_m),\ \ Q(y_1,y_2,...y_n)
P=(x1,x2,...xm), Q(y1,y2,...yn)
令
A
m
×
n
=
(
a
i
j
)
m
×
n
A_{m\times n}=(a_{ij})_{m\times n}
Am×n=(aij)m×n 表示
P
P
P 和
Q
Q
Q 的距离矩阵,其中
a
i
j
a_{ij}
aij 基于基距离方法进行计算,一般采用欧氏距离的平方:
a
i
j
=
(
x
i
−
y
j
)
2
a_{ij} = (x_i-y_j)^2
aij=(xi−yj)2
然后在距离矩阵中,令
W
=
w
1
,
w
2
,
.
.
.
,
w
k
W=w_1,w_2,...,w_k
W=w1,w2,...,wk为时间序列
P
P
P 和
Q
Q
Q 的动态时间弯曲路径,其中
w
k
=
(
a
i
j
)
k
w_k=(a_{ij})_k
wk=(aij)k 是路径的第k个元素,且路径
W
W
W 需要满足以下条件:
1、 max { m , n } ≤ K ≤ m + n − 1 \max\{m,n\}\leq K\leq m+n-1 max{m,n}≤K≤m+n−1
2、 w 1 = a 11 , w k = a m n w_1=a_{11}, w_k=a_{mn} w1=a11,wk=amn
3、单调性:若
w
k
=
a
i
j
,
w
k
+
1
=
a
k
l
w_k=a_{ij}, w_{k+1}=a_{kl}
wk=aij,wk+1=akl,则需保证
0
≤
k
−
i
≤
1
,
0
≤
l
−
j
≤
1
0\leq k-i\leq 1, \ \ 0\leq l-j\leq 1
0≤k−i≤1, 0≤l−j≤1
根据上述三个条件可知,动态时间弯曲路径并不唯一,在所有的动态时间弯曲路径中使得
∑
i
=
1
K
w
i
\sqrt{\sum\limits_{i=1}^K w_i}
i=1∑Kwi 的值最小的路径称为最佳动态时间弯曲路径,所对应的距离即为动态时间弯曲距离。记
D
T
W
(
P
,
Q
)
DTW (P,Q)
DTW(P,Q) 为P和Q之间的DTW距离,则有:
D
T
W
(
P
,
Q
)
=
min
(
∑
i
=
1
K
w
i
)
DTW(P,Q) = \min\left(\sqrt{\sum\limits_{i=1}^{K}w_i}\right)
DTW(P,Q)=min⎝⎛i=1∑Kwi⎠⎞
为什么称DTW为动态时间规整距离呢?
因为DTW距离的计算过程是一个 动态规划过程, 若令
P
P
P 和
Q
Q
Q 在
i
i
i 和
j
j
j 处累积距离为
L
(
i
,
j
)
L(i,j)
L(i,j),则由条件3可知,累积距离存在如下递归关系:
L
(
i
,
j
)
=
min
[
L
(
i
−
1
,
j
−
1
)
,
L
(
i
−
1
,
j
)
,
L
(
i
,
j
−
1
)
]
+
a
i
j
L(i,j) = \min\left[L(i-1,j-1), L(i-1,j), L(i,j-1)\right]+a_{ij}
L(i,j)=min[L(i−1,j−1),L(i−1,j),L(i,j−1)]+aij
易知:
L
(
1
,
1
)
=
a
11
L(1,1)=a_{11}
L(1,1)=a11, 且
D
T
W
(
P
,
Q
)
=
L
(
m
,
n
)
DTW(P,Q)=L(m,n)
DTW(P,Q)=L(m,n),因此DTW的计算过程就是:
- 先使用基距离,如欧式距离初始化距离矩阵 A m × n = ( a i j ) m × n A_{m\times n}=(a_{ij})_{m\times n} Am×n=(aij)m×n
- 然后使用累积距离的递归公式求解得到 L ( m , n ) L(m,n) L(m,n) 即为DTW距离
注:
- DTW距离通过允许序列在时间轴上偏移,克服了闵可夫斯基距离不能度量不等长序列、对于噪声和振幅伸缩无法刻画的缺陷
- 可以度量任意两个长度时间序列间的距离,对时间轴上的扭曲不敏感,适用性、可靠性更强
- 计算复杂度高: O ( m n ) O(mn) O(mn) , 时间代价昂贵
- DTW距离不满足三角不等式
- 相应改进方案:下界函数、加权动态时间规整距离、微分动态时间规整距离、FastDTW
基于相关性的相似度度量方法
给定两个时间序列:
P
=
(
x
1
,
x
2
,
.
.
.
x
n
)
,
Q
(
y
1
,
y
2
,
.
.
.
y
n
)
P=(x_1,x_2,...x_n),\ \ Q(y_1,y_2,...y_n)
P=(x1,x2,...xn), Q(y1,y2,...yn)
相关系数定义如下:
C
o
r
(
P
,
Q
)
=
C
o
v
(
P
,
Q
)
V
a
r
(
P
)
V
a
r
(
Q
)
Cor(P,Q) = \frac{Cov(P,Q)}{\sqrt{Var(P)\sqrt{Var(Q)}}}
Cor(P,Q)=Var(P)Var(Q)Cov(P,Q)
相关系数为1,表示完全一致;相关系数为-1,表示完全负相关;相关系数为0,表示没有相关性
基于互信息的相似性度量方法
给定两个时间序列:
P
=
(
x
1
,
x
2
,
.
.
.
x
n
)
,
Q
(
y
1
,
y
2
,
.
.
.
y
n
)
P=(x_1,x_2,...x_n),\ \ Q(y_1,y_2,...y_n)
P=(x1,x2,...xn), Q(y1,y2,...yn)
互信息定义如下:
I
(
P
,
Q
)
=
H
(
P
)
+
H
(
Q
)
−
H
(
P
,
Q
)
I(P,Q)=H(P)+H(Q)-H(P,Q)
I(P,Q)=H(P)+H(Q)−H(P,Q)
互信息取值范围为
[
0
,
∞
]
[0,\infty]
[0,∞],互信息为0表示两者之间没有公共信息,互信息为
∞
\infty
∞表示相关性特别大,类似于相关系数的绝对值为1. 除此之外,可利用如下公式对互信息进行归一化:
Λ
(
P
,
Q
)
=
1
−
e
2
I
(
P
,
Q
)
\Lambda(P,Q) = \sqrt{1-e^{2I(P,Q)}}
Λ(P,Q)=1−e2I(P,Q)
基于互信息的相似性度量方法(可能拓展)
上面基于互信息的相似性度量并没有考虑到不同时间序列之间的时间延迟关系,如果他们存在一定程度的时间延迟,那么可以考虑如下方式的拓展:
- 计算延迟互信息,然后对不同的时间延迟所得互信息求最大值或者直接求和
- 计算转移熵,然后求最大值或者直接求和
马氏距离
上面提到了基于相关性和基于互信息的距离,仅仅考虑到了变量之间的相关性,而马氏距离则是同时考虑到变量不同维度之间存在相关性和尺度变换的关系,它表示的是变量之间的协方差距离。
假设有两个n维空间中的点:
P
=
(
x
1
,
x
2
,
.
.
.
x
n
)
,
Q
(
y
1
,
y
2
,
.
.
.
y
n
)
P=(x_1,x_2,...x_n),\ \ Q(y_1,y_2,...y_n)
P=(x1,x2,...xn), Q(y1,y2,...yn)
那么马氏距离的定义为:
D
M
(
P
,
Q
)
=
(
P
−
Q
)
T
Σ
−
1
(
P
−
Q
)
D_M(P,Q) = \sqrt{(P-Q)^T\Sigma^{-1}(P-Q)}
DM(P,Q)=(P−Q)TΣ−1(P−Q)
其中
Σ
\Sigma
Σ 是个协方差矩阵,但是貌似跟
P
P
P 和
Q
Q
Q 没啥关系,它应该是个
n
×
n
n\times n
n×n 的矩阵, 刻画了n个属性之间的关系,计算方式如下:
- 对n维变量进行采样,采样得到样本集合: { P 1 , P 2 , . . . , P m } \{P_1,P_2,...,P_m\} {P1,P2,...,Pm},这个集合的规模可能比 n n n 大很多, 然后每个属性有一个时间序列,计算各个属性之间的协方差矩阵,以及各个属性的均值组成的均值向量 μ = [ μ 1 , μ 2 , . . . , μ n ] \mu=[\mu_1,\mu_2,...,\mu_n] μ=[μ1,μ2,...,μn]
如果利用矩阵分解,如LU分解对协方差矩阵
Σ
\Sigma
Σ 进行分解的话:
Σ
=
L
L
T
\Sigma=LL^T
Σ=LLT,然后对
P
−
Q
P-Q
P−Q做如下处理:
Z
=
L
−
1
(
P
−
Q
)
Z = L^{-1}(P-Q)
Z=L−1(P−Q)
则马氏距离可以表示为:
D
M
(
P
,
Q
)
=
Z
T
Z
D_M(P,Q) = Z^TZ
DM(P,Q)=ZTZ
特别地,马氏距离可以用来衡量一个点 P 和样本集合
{
P
1
,
P
2
,
.
.
.
,
P
m
}
\{P_1,P_2,...,P_m\}
{P1,P2,...,Pm} 之间的距离,这是通过比较 P 和均值向量
μ
\mu
μ 得到的:
D
M
(
P
)
=
(
P
−
μ
)
Σ
−
1
(
P
−
μ
)
D_M(P) = \sqrt{(P-\mu)\Sigma^{-1}(P-\mu)}
DM(P)=(P−μ)Σ−1(P−μ)
注:
- 马氏距离计算的不是时间序列之间的距离,而是对欧氏距离的推广,计算的是两个点之间的距离
- 马氏距离需要一个样本集合来估算样本均值和协方差。
- 马氏距离相当于将变量按照主成分进行旋转,让维度之间相互独立,然后进行标准化,让维度同分布
- 主成分就是特征向量方向,每个方向的方差就是对应的特征值,所以其实就是按照特征向量的方向进行旋转,然后缩放特征值倍
优点:
- 马氏距离不受量纲的影响,且由标准化数据和中心化数据计算出的两点之间的马氏距离相同
- 马氏距离可以排除变量之间相关性的干扰
缺点:
- 需要计算协方差矩阵,协方差矩阵还需要可逆
- 涉及到矩阵求逆
- 不能处理非线性流形上的问题
参考文献
时间序列聚类相关知识的总结与梳理: https://mp.weixin.qq.com/s/iusTb9UwKybwJBqd2kSvUA
机器学习基础 | 互相关系数和互信息异同探讨:https://www.cnblogs.com/fangsf/p/15000465.html
马氏距离:https://blog.csdn.net/qq_28087491/article/details/114370480?ops_request_misc=%257B%2522request%255Fid%2522%253A%2522162639931516780274116169%2522%252C%2522scm%2522%253A%252220140713.130102334…%2522%257D&request_id=162639931516780274116169&biz_id=0&utm_medium=distribute.pc_search_result.none-task-blog-2alltop_click~default-2-114370480.first_rank_v2_pc_rank_v29&utm_term=%E9%A9%AC%E6%B0%8F%E8%B7%9D%E7%A6%BB&spm=1018.2226.3001.4187
荣馨兵. 基于相似搜索的时间序列预测方法及其应用[D]. 大连海事大学, 2018.
张勇. 时间序列模式匹配技术研究[D]. 华中科技大学. 2012