【国庆特辑文章】时间序列~动态时间规整(Dynamic Time Wraping)

时间序列-动态时间规整(Dynamic Time Wraping)

一、解决的问题

测量两端时间序列的相似性。

在语音识别中，特别是单音节的识别中，每个人说话时间长短不同，导致时间序列长度不同。DTW算法就是将某些数据点的时间Wrap到另一个时间序列的某些数据点，以辅助计算相似性。

二、算法设计

1.规则

• 两端对齐
• 一个点可以对应另一序列的多个点（允许重合对应)
• 每一对数据点对齐不可交叉（只能向前对应）

2.成本矩阵

• 有两不同长度时间序列
  $$X=x_1,x_2,......,x_n$$
  $$Y=y_1,y_2,......,y_m$$

• 构建距离矩阵
  $$D_N{\times }_M$$
  矩阵元素
  $$d_{i}j=dist(x_i,y_i)$$
  dist通常采用欧氏距离

• 在矩阵D中搜索从$$d_{1}1$$到$$d_{(}n{m}_{)}$$之间的最短路径（采用动态规划搜索算法，向右、上、右上三个方向搜索)；
• 这个最短路径的和作为X和Y之间的相似度。

三、代码实现

'''
@Author: Classmate.Liu loved Technology
@Date: 2022-10-01 08:33:00
@LastEditTime: 2022-10-01 09:40:22
@FilePath: D:\A\Project_1\main_5.py
'''
# import package
from importlib.resources import path
import numpy as np
import matplotlib.pyplot as plt

np.random.seed(42)

# define operation functions (including lambda expressions)
def DTW(X, Y, distance = lambda a, b: np.linalg.norm(a -b)):
    ''' dynamic time warping '''
    N, M = len(X), len(Y)
    cost = np.ones([N, M]) * np.inf
    for i in range(N):
        for j in range(M):
            dist_ij = distance(X[i], Y[j])
            dist_pr = min(cost[i - 1, j] if i-1>=0 else np.inf, cost[i, j-1] if j-1>=0 else np.inf, cost[i, j-1] if j-1>=0 and i+1>=0 else np.inf)
            cost[i, j] = dist_ij + (dist_pr if dist_pr < np.inf else 0)
        
        # traced back cost matrix to get minimum distance path
        i, j = N - 1, M - 1
        path = [(i, j)]
        while i > 0 or j > 0:
            condidate = []
            if i-1>=0:
                condidate.append((cost[i-1, j], (i-1, j)))
            if j-1>=0:
                condidate.append((cost[i, j-1], (i, j-1)))
            if i-1>=0 or j-1>=0:
                condidate.append((cost[i-1, j-1], (i-1, j-1)))
            i, j = min(condidate)[1]
            path.append((i, j))
        
        return cost[N-1, M-1], path

g_X = np.random.uniform(size=18)
g_Y = np.random.uniform(size=16) + 3.0

dist, path = DTW(g_X, g_Y)
print('Dist(X, Y) = ', dist)

plt.plot(g_X)
plt.plot(g_Y)
for ij in path:
    plt.plot(ij[0], ij[1]),
[g_X[ij[0]], g_Y[ij[1]], 'gray']
plt.show()

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53

四、运行结果

图一：窗体运行结果

图二：终端运行结果

五、性能优化

• 绝大多数最短距离不会偏离对角线太远；

• 先在低精度下确定轮廓，随后在轮廓内计算更高精度路径。