时间序列中的动态时间扭曲算法解析

1,190 阅读10分钟

Dynamic Time Warping Algorithm in Time Series, Explained
使用运动捕捉系统记录了两次重复的行走序列。虽然不同的重复在行走速度上存在差异,但四肢的空间路径仍然高度相似。

简介

乍一看,"动态时间扭曲 "这个短语可能会让人联想到《*回到未来》*系列中马蒂-麦克弗里以88英里/小时的速度驾驶他的迪罗瑞安。然而,动态时间扭曲并不涉及时间旅行;相反,当比较数据点之间的时间指数不完全同步时,它是一种用于动态比较时间序列数据的技术。

Introduction

正如我们将在下面探讨的那样,动态时间扭曲最突出的用途之一是在语音识别中--确定一个短语是否与另一个短语相匹配,即使该短语的语速比其对比的快或慢。你可以想象,这在识别用于激活你的谷歌Home或亚马逊Alexa设备的 "唤醒词 "时很方便--即使你的语音很慢,因为你还没有喝完你每天的一杯(几杯)咖啡。

动态时间扭曲(DTW)

时间序列分析中,动态时间扭曲(DTW)是测量两个时间序列序列之间相似性的算法之一,这两个序列的速度可能不同。DTW的主要思想是计算时间序列之间相似元素的匹配距离。它使用动态编程技术来寻找两个时间序列元素之间的最佳时间匹配。

例如,即使一个人比另一个人走得快,或者在观察期间有 加速和减速的情况,也可以用DTW检测出走路的相似性。DTW已经被应用于视频、音频和图形数据的时间序列--事实上,任何可以转化为线性序列的数据都可以用DTW进行分析。一个著名的应用是自动 语音识别,以应对不同的说话速度。其他应用包括 说话人识别和在线 签名识别。它还可以用于部分 形状匹配应用。

digital_geometry

时间序列比较方法的目的是在两个输入时间序列之间产生一个距离度量。两个时间序列的相似性或不相似性通常是通过将数据转换为矢量并计算矢量空间中这些点之间的欧氏距离来计算的。

一般来说,DTW是一种计算两个给定序列(如时间序列)之间最佳匹配的方法,有一定的限制和规则。

  1. 第一个序列中的每个指数必须与另一个序列中的一个或多个指数相匹配,反之亦然

  2. 第一个序列的第一个索引必须与另一个序列的第一个索引相匹配(但它不一定是唯一的匹配)。

  3. 第一个序列的最后一个索引必须与另一个序列的最后一个索引匹配(但它不一定是唯一的匹配)。

  4. 第一序列的指数与另一序列的指数的映射必须是单调递增的,反之亦然,即如果j>i是第一序列的指数,那么在另一序列中不能有两个指数l>k,从而使指数i与指数l匹配,指数j与指数k匹配,反之亦然。

Time normalized distance between Series A and B.
A系列和B系列之间的时间归一化距离。

从数学上理解DTW

让我们假设我们有两个像下面这样的序列。

Understand DTW Mathematically

序列X和Y可以排列成一个n乘m的网格,其中每个点(i,j)是xi和yj之间的排列。

一个翘曲路径W映射X和Y的元素,使它们之间的距离最小。W是一个网格点(i,j)的序列。我们将在后面看到一个翘曲路径的例子。

翘曲路径和DTW距离

通往(i_k, j_k)的最优路径可以通过以下方式计算。

Optimal path

其中d是欧几里得距离。然后,整体路径成本可以计算为:。

overall path cost can be calculated as

对翘曲函数的限制

翘曲路径是用动态编程的方法找到的,以对齐两个序列。遍历所有可能的路径是 "组合爆炸性的"。因此,为了提高效率,限制可能的翘曲路径的数量是很重要的,因此概述了以下限制条件。

  • 边界条件。这个约束条件确保翘曲路径从两个信号的起点开始,并以它们的终点结束。

    Boundary Condition

  • 单调性条件。这个约束条件保留了各点的时间顺序(不回溯)。

    Monotonicity condition

  • 连续性(步骤大小)条件。这个约束条件将路径转换限制在相邻的时间点上(而不是在时间上跳转)。

    Continuity (step size) condition

除了上述三个约束条件外,还有其他一些不太常见的可允许的翘曲路径的条件。

  • 翘曲窗口条件。允许的点可以被限制在一个给定的宽度为Equation (正整数)的翘曲窗口内。

    Warping window condition

  • 坡度条件。可以通过限制坡度来约束翘曲路径,从而避免在一个方向的极端运动。

一个可接受的翘曲路径的棋王动作组合是:。

  • 水平移动:(i,j)->(i,j+1)。
  • 垂直移动:(i, j) -> (i+1, j)
  • 对角线移动:(i, j) -> (i+1, j+1)

为什么我们需要DTW?

任何两个时间序列都可以用欧氏距离或其他类似的距离在时间轴上一对一地进行比较。第一个时间序列在时间T的振幅将与第二个时间序列在时间T的振幅进行比较。即使两个时间序列在形状上非常相似,但在时间上不相一致,这也会导致一个非常差的比较和相似性分数。

DTW compares the amplitude of the first signal at time

DTW将第一个信号在时间T的振幅与第二个信号在时间T+1和T-1或T+2和T-2的振幅进行比较。这确保了它不会对具有相似形状和不同相位的信号给出一个低的相似性分数。

DTW算法是如何工作的?

让我们取两个时间序列信号A和B。

two-time series signals A and B

signals A and B

步骤1:创建空成本矩阵

创建一个空的成本矩阵M,用x和y标签作为两个序列的振幅来进行比较。

Empty Cost Matrix

第2步:成本计算

使用下面提到的公式从左角和底角开始填充成本矩阵。

M(i, j) = |A(i) - B(j)| + min ( M(i-1,j-1), M(i, j-1), M(i-1,j) )

其中,
M是矩阵

i是A系列的迭代器

j是B系列的迭代器

所以在下表中,我们计算了一些数值的成本。

computed the cost of some values

让我们从上表中抽取几个例子(6、9和12)来说明下表中所强调的计算过程

illustrate the calculation as highlighted

对于6。

|8 -5| + min( 6, 6, 3 ) = 6

同样的,对于9来说。

|4 -1| + min( 6 ) = 9

12也是如此。

|1 -5| + min( 8 ) = 12

完整的表格将看起来像这样。

full table

第3步:确定翘曲路径

识别从矩阵的右上角开始,遍历到左下角的翘曲路径。遍历路径的确定是基于具有最小值的邻居。

在我们的例子中,它从17开始,在其邻居18、14和12中寻找最小值,即12。

Warping Path Identification

Warping Path Identification

这个过程一直持续到我们到达表格左轴的底部。

warping traversal path

最后的路径会是这样的。

final path will look like

让我们把这个翘曲路径系列称为d。

d = [17, 12, 9, 9, 9, 7, 7, 6, 3, 2, 1, 0]

第四步:最终距离计算

时间归一化的距离,D

Time normalized distance

其中k是序列d的长度。

在我们的例子中,k=12。

D = ( 17 + 12 + 9 + 9 + 9 +7 + 7 + 6 + 3 + 2 +1 + 0)/ 12

D = 6.8333

Python中的DTW实现

为了实现,我们将使用 fastdtwpython库。

FastDTW一种近似的动态时间扭曲(DTW)算法,它以O(N)的时间和内存复杂度提供最优或接近最优的排列,与标准DTW算法的O(N^2)要求形成对比。

安装库。

pip install fastdtw

导入所有需要的库

import pandas as pd

import numpy as np
# Plotting Packages

import matplotlib.pyplot as plt

import seaborn as sbn
import matplotlib as mpl

mpl.rcParams['figure.dpi'] = 150

savefig_options = dict(format="png", dpi=150, bbox_inches="tight")
# Computation packages

from scipy.spatial.distance import euclidean

from fastdtw import fastdtw

让我们定义一个方法来计算翘曲路径的累积成本矩阵D。成本矩阵使用欧氏距离来计算每两个点之间的距离。计算欧氏距离矩阵和累积成本矩阵的方法定义如下。

def compute_euclidean_distance_matrix(x, y) -> np.array:

    """Calculate distance matrix

    This method calculates the pairwise Euclidean distance between two sequences.

    The sequences can have different lengths.

    """

    dist = np.zeros((len(y), len(x)))

    for i in range(len(y)):

        for j in range(len(x)):

            dist[i,j] = (x[j]-y[i])**2

    return dist


def compute_accumulated_cost_matrix(x, y) -> np.array:

    """Compute accumulated cost matrix for warp path using Euclidean distance

    """

    distances = compute_euclidean_distance_matrix(x, y)
    # Initialization

    cost = np.zeros((len(y), len(x)))

    cost[0,0] = distances[0,0]
    for i in range(1, len(y)):

        cost[i, 0] = distances[i, 0] + cost[i-1, 0]  
    for j in range(1, len(x)):

        cost[0, j] = distances[0, j] + cost[0, j-1]  
    # Accumulated warp path cost

    for i in range(1, len(y)):

        for j in range(1, len(x)):

            cost[i, j] = min(

                cost[i-1, j],    # insertion

                cost[i, j-1],    # deletion

                cost[i-1, j-1]   # match

            ) + distances[i, j] 
    return cost

现在,创建两个序列。

x = [7, 1, 2, 5, 9]

y = [1, 8, 0, 4, 4, 2, 0]

我们不能计算x和y之间的距离,因为它们的长度不相等。

fig, ax = plt.subplots(figsize=(6, 4))
# Remove the border and axes ticks

fig.patch.set_visible(False)

ax.axis('off')

xx = [(i, x[i]) for i in np.arange(0, len(x))]

yy = [(j, y[j]) for j in np.arange(0, len(y))]

for i, j in zip(xx, yy[:-2]):

  ax.plot([i[0], j[0]], [i[1], j[1]], '--k', linewidth=1)

ax.plot(x, '-ro', label='x', linewidth=1, markersize=20, markerfacecolor='lightcoral', markeredgecolor='lightcoral')

ax.plot(y, '-bo', label='y', linewidth=1, markersize=20, markerfacecolor='skyblue', markeredgecolor='skyblue')

ax.set_title("Euclidean Distance", fontsize=10, fontweight="bold")

Euclidean distance between x and y
x和y之间的欧几里得距离

计算DTW距离和翘曲路径

许多Python软件包只需提供序列和距离类型(通常默认为欧几里得)就可以计算DTW。在这里,我们使用一个流行的DTW的Python实现,即 FastDTW,这是一个近似的DTW算法,具有较低的时间和内存复杂度。

dtw_distance, warp_path = fastdtw(x, y, dist=euclidean)

请注意,我们使用的是 SciPy的距离函数Euclidean,我们之前已经导入了这个函数。为了更好地理解翘曲路径,让我们首先计算累积成本矩阵,然后在网格上可视化路径。下面的代码将绘制出累积成本矩阵的热图。

cost_matrix = compute_accumulated_cost_matrix(x, y)

fig, ax = plt.subplots(figsize=(6, 4))

ax = sbn.heatmap(cost_matrix, annot=True, square=True, linewidths=0.1, cmap="YlGnBu", ax=ax)

ax.invert_yaxis()

# Get the warp path in x and y directions

path_x = [p[0] for p in warp_path]

path_y = [p[1] for p in warp_path]

# Align the path from the center of each cell

path_xx = [x+0.5 for x in path_x]

path_yy = [y+0.5 for y in path_y]

ax.plot(path_xx, path_yy, color='blue', linewidth=1, alpha=0.2)

path plot

色条显示了网格中每个点的成本。可以看出,翘曲路径(蓝线)正穿过网格中成本最低的地方。让我们通过打印这两个变量来看看DTW距离和翘曲路径。

print("DTW distance: ", dtw_distance)

print("Warp path: ", warp_path)
>>> DTW distance:  23.0

>>> Warp path:  [(0, 0), (0, 1), (1, 2), (2, 3), (3, 4), (3, 5), (4, 6)]

翘曲路径从点(0,0)开始,通过6次移动在(4,6)结束。让我们也用我们前面定义的函数计算一下累计成本最,并将其值与热图进行比较。

cost_matrix = compute_accumulated_cost_matrix(x, y)

print(np.flipud(cost_matrix))

cost_matrix

上面打印的成本矩阵与热图的数值相似。

现在让我们绘制这两个序列,并连接映射点。下面给出了绘制xy之间DTW距离的代码。

fig, ax = plt.subplots(figsize=(6, 4))

# Remove the border and axes ticks

fig.patch.set_visible(False)

ax.axis('off')

for [map_x, map_y] in warp_path:

    ax.plot([map_x, map_y], [x[map_x], y[map_y]], '--k', linewidth=1)

ax.plot(x, '-ro', label='x', linewidth=1, markersize=20, markerfacecolor='lightcoral', markeredgecolor='lightcoral')

ax.plot(y, '-bo', label='y', linewidth=1, markersize=20, markerfacecolor='skyblue', markeredgecolor='skyblue')

ax.set_title("DTW Distance", fontsize=10, fontweight="bold")

DTW distance between x and y
x和y之间的DTW距离

你可以从 Github下载该代码。

DTW的应用

  1. 检测行走的相似性。如果一个人比另一个人走得快,或者在观察期间有加速和减速的情况。

  2. 口语识别应用。它们用于将样本语音命令与其他人的命令相匹配,即使这个人说话比预先录制的样本语音更快或更慢。

  3. 相关的功率分析

参考文献