从DTW到DILATE---时序预测损失函数总结

时间序列预测损失函数

时间序列预测是指通过分析历史数据来预测未来数据的变化趋势。时间序列预测在许多领域都有着广泛的应用，例如金融、气象、交通等。为了能够准确预测未来数据，时间序列预测中损失函数的选择非常关键。损失函数的选择直接影响模型的预测精度，因此，对时间序列预测中常用的损失函数进行综述和评价具有重要的意义。

在时间序列预测中，常用的损失函数包括均方误差（Mean Squared Error, MSE）、平均绝对误差（Mean Absolute Error, MAE）、动态时间规整（Dynamic Time Warping, DTW）等。这些损失函数各有优劣，不同的损失函数适用于不同的预测场景。

本篇综述将对时间序列预测中常用的损失函数进行系统性地介绍和分析，并比较不同损失函数的优缺点和适用场景，以便研究者和实践者在实际应用中选择最适合的损失函数。

一、DTW

DTW是一种常用的时间序列相似性度量方法，主要用于衡量两个时间序列间的相似性。

对于语音序列，由于说话人的发音风格不同，导致同样内容的两段语音即使在序列形状十分相似的情况下，也很难实现在时间点上的一一对齐，出现波动延后，起伏区间不同的情况（图左）。此时通过欧氏距离方法很难实现两序列间在时间轴上的距离最小化。

欧式距离：
$$d=||y_i-{\hat{y}}_i|{|}^2$$
均方误差（MSE）与均方根误差（RMSE）
$$\begin{gathered} L_{MSE}=\frac{1}{n}\sum _{i=1}^n(Y_i-{\hat{Y}}_i{)}^2 \\ {L}_{RMSE}=\sqrt{\frac{1}{n}\sum _{i=1}^n(Y_i-{\hat{Y}}_i{)}^2} \end{gathered}$$
如果将机器学习模型 $Y$ 的输出看作是一个向量，那么 MSE 求和的部分就是在求估计样本 $Y$ 和 $\hat{Y}$ 之间的欧氏距离。即MSE 是机器学习中的一个概念，使用了欧氏距离来衡量模型输出的准确度。RMSE同理，使用欧氏距离来计算误差。

DTW引入 时间轴扭曲 的概念，允许时间轴在不同位置进行弯曲和拉伸，使得序列上的点实现“一对多”、“多对一”的对应（图右），通过这种方法，DTW中定义两序列间的距离为二者间最小的时间轴扭曲代价，该代价通过计算两序列对应时间点之间的距离，并找到一条使两序列间距离最小化的最短路径而得到的。

DTW 算法简介

DTW的本质是一种动态规划算法，对于两个长度分别为 $m,n$ 的时间序列 $A,B$，存在
$$\begin{gathered} A=a_1,a_2,...,a_n \\ B=b_1,b_2,...,b_m \end{gathered}$$
以二者长度构建一个大小为 $m\times n$ 的矩阵，矩阵内位置 $(i,j)$ 负责存储点 $a_i$ 与 $b_j$ 之间的欧氏距离。DTW算法则是负责从矩阵中找出一条从 $(1,1)\to (m,n)$ 代价最小的一条路径 $P$（warping path）

上图中，路径 $P$ 是矩阵中的连续元素，其第 $s$ 个元素定义为 $p_s=(i,j{)}_s$，$P=p_1,p_2,...,p_s,...,p_k$ 。

$A,B$ 间的时间正规化距离（Time-normalized distance） $D_p(A,B)=\left[\frac{{\sum }_{s=1}^k\delta (p_s)\cdot w_s}{{\sum }_{s=1}^k{w}_s}\right]$ ，$i_s,j_s$ 间距离为$\delta (p_s)$ 表示两点间的距离，$w_s>0$ 为加权系数。

以此DTW算法的优化目标是找到矩阵中的最佳路径 $P_0$ ，并满足 $P_0=\underset{P}{\mathrm{argmin}}(D_p(A,B))$，为减小搜索空间大小，需要增添下述限制条件：

• 边界限制：$P_i=(1,1)，P_k=(m,n)$，保证开始与最后时刻对齐，路径从左下角出发，右上角结束
• 连续性：满足 $i_{s-1}-i_s\le 1,j_{s-1}-j_s\le 1$，即路径连续，不存在跳步情况
• 单调性：满足 $i_{s-1}\le i_s,j_{s-1}\le j_s$ ，即路径随时间单调前进，不存在回头的情况
• 窗口规整：$|i_s-j_s|\le window,window\ge 0$ ，确保路径不会在某些特征过度停留，而是尽量沿对角线。
• 斜率限制：$\frac{j_{s_p}-j_{s_0}}{i_{s_p}-i_{s_0}}\le p，\frac{i_{s_q}-i_{s_0}}{j_{s_q}-j_{s_0}}\le q,p,q\ge 0$，其中 $p,q$ 为 $y,x$ 方向步数。

根据上述约束条件，当 $P$ 出发时，仅有上、右上以及右侧三个方向可选，同时由于边界限制，$P$ 仅可在 $window$ 范围内运动，设点 $(a_i,b_j)$ 间距离为 $d(a_i,b_j)=\sqrt{(a_i-b_j{)}^2}={\delta }_{i,j}\cdot w_{i,j}$， 因此对于该动态规划问题 $\tau$ ，$(1,1)\to (i,j)$ 的最小累加距离为：
τi,j=δi,j⋅wi,j+min{τi,j−1,τi−1,j,τi−1,j−1}∴DP0(A,B)=1C⋅τn,m \tau_{i,j}=\delta_{i,j} \cdot w_{i,j}+min\{\tau_{i,j-1}, \tau_{i-1,j}, \tau_{i-1,j-1}\} \\ \therefore D_{P_0}(A,B)=\frac{1}{C} \cdot \tau_{n,m} τi,j​=δi,j​⋅wi,j​+min{τi,j−1​,τi−1,j​,τi−1,j−1​}∴DP0​​(A,B)=C1​⋅τn,m​
在实际应用中，在训练阶段，用户的每段语音被记录下来存入模板库，随后识别阶段中，通过DTW算法将采集到的特征与模板库中进行对比，找到最短距离，实现准确对齐。

由于要在 $m\times n$ 的矩阵中采用暴力动态规划方法求最优解，故DTW 算法时间和空间复杂度为：$O(n^2)$

示例

参考资料：

• https://blog.youkuaiyun.com/qq_39516859/article/details/81705010
• https://zhuanlan.zhihu.com/p/346674141

二、Fast-DTW

文章：FastDTW: Toward Accurate Dynamic Time Warping in Linear Time and Space

贡献：实现了线性时间、空间下DTW算法的精确近似，将DTW算法的时间和空间复杂度由 $O(n^2)$ 成功降低为 $O(n)$

Fast-DTW 算法简介

Fast-DTW采用的多层方法主要受图平分的多层方法启发，通过将图拆分为大小相似的多个子图，在子图中分别找出最优解，随后将其拓展为大图，并对所求最优解进行更新，以求得最终最优解。相较于标准DTW算法，采用多级方法的Fast-DTW无需暴力考虑整个矩阵，而是仅将代价矩阵填充在前一个分辨率所投影的路径内进行考虑即可。

因此，Fast-DTW的加速主要通过三个步骤：

1. 粗粒度化：通过平均序列中相邻点进行数据抽象，对原始时间序列进行压缩，获得不同尺度的序列结果。
2. 投影：在粗粒度化后获取的序列上计算一条初始路径，将其投影到高分辨率时间序列中，作为高分辨率序列中路径的启发
3. 细粒度化：在高分辨率序列中投影路径区域内进行局部微调，以寻找最优解，通常会额外扩张 $K$ 个粒度，$K$ 为半径参数，用于控制邻域大小，一般取值为 $1$ 或 $2$ 。

以下图为例，对于原始矩阵，Fast-DTW算法首先通过粗化创建其四个子分辨率矩阵用于后续计算：

1. 在 $1/8$ 粗粒度空间中运行标准DTW算法，获取临时最优路径
2. 将获取到的临时最优路径投影到 $1/4$ 分辨率空间中作为启发，即阴影部分。将其作为DTW算法的约束空间，以求得临时最优解。同时考虑到最优路径可能不包含在投影路径中，因此通过半径参数设置投影路径边的额外空间参与DTW计算。
3. 在 $1/4$ 分辨率空间中获取到临时最优后，将其再次细化，投影到 $1/2$ 分辨率空间中再次进行细化
4. 将细化后的投影路径映射到全尺寸矩阵 $(1/1)$ 中，进行最后一次细化，以求得最终解。

通过上述步骤，算法的搜索空间由整个矩阵减小为一个随时间序列线性增长的局部空间，从而将时间复杂度由 $O(n^2)$ 降低为 $O(n)$。

但是需要注意的是，Fast-DTW算法并非总能找到最优解，但其找到的路径基本都非常接近最优。同时半径参数的存在也对最优路径的获取存在积极作用，半径参数 $K$ 越大，获取到最优路径的可能性越高，但时间和空间复杂度也会随之上升，当 $K$ 的大小与一条时间序列的长度相同时，时间和空间复杂度增加到 $O(n^2)$ ，同时Fast-DTW被推广到DTW。

三、Soft-DTW

文章：Soft-DTW: a Differentiable Loss Function for Time-Series

贡献：基于动态规划的标准DTW算法属于离散不可微计算，无法用于深度学习中神经网络的损失函数计算，本文采用Soft minimum取代DTW minimum，将DTW由离散不可微计算拓展为连续可微的损失函数，实现了通过梯度下降进行函数结果更新。

Soft-DTW 算法简介

对于两序列 $x=(x_1,x_2,...,x_n)\in {\mathbb{R}}^{p\ast n}$ 与 $y=y_1,y_2,...,y_m\in {\mathbb{R}}^{p\ast m}$ ，与DTW相似，定义其代价矩阵为 $△(x,y):=[\delta (x_i,y_i){]}_{ij}\in {\mathbb{R}}^{n\ast m}$ 。$\delta$ 为可微代价函数，定义 $\delta :{\mathbb{R}}^p\times {\mathbb{R}}^p\to {\mathbb{R}}_{+}$，即时间序列在某时刻的信息均为实数值，且该代价函数表示为两对应点间欧氏距离。

定义集合 $R=[r_{i,j}],R\in {\mathbb{R}}^{n\times m}$ 为路径上的代价和 $r_{i,j}$ 组成的集合，则对于DTW，其动态规划式如下：
ri,jDTW=δi,j⋅wi,j+min{ri,j−1,ri−1,j−1,ri−1,j} r_{i,j}^{DTW}=\delta_{i,j} \cdot w_{i,j} + min\{r_{i,j-1},r_{i-1,j-1},r_{i-1,j}\} ri,jDTW​=δi,j​⋅wi,j​+min{ri,j−1​,ri−1,j−1​,ri−1,j​}
$min$ 为一个离散过程，其存在导致了DTW的离散。

在Soft-min中，采用 $mi{n}^{\gamma }$ 替换了原始公式中的 $min$ 部分，具体公式如下：
minγ{a1,a2,...,an}={mini≤naiγ=0−γlog⁡∑i=1ne−ai/γγ>0 min^{\gamma}\{a_1,a_2,...,a_n\}= \begin{cases} &min_{i \leq n}a_i \qquad \qquad \gamma=0\\ &-\gamma\log\sum^n_{i=1}e^{-a_i/\gamma} \qquad \gamma>0 \end{cases} minγ{a1​,a2​,...,an​}={​mini≤n​ai​γ=0−γlog∑i=1n​e−ai​/γγ>0​
当 $\gamma \to 0$ 时，$mi{n}^{\gamma }\to min$。

已知 $max$ 函数的平滑化处理如下：

根据泰勒展开式，存在 $f(X)\approx f(c)+f^{\prime }(c)(X-c)$

设 $f(X)=\log (X),X=a+x,c=x$

故 $\log (x+a)=f(x+a)\approx f(x)+f^{\prime }(x)(x+a-x)=f(x)+\frac{a}{x}=log(x)+\frac{a}{x}$

即 $log(x+a)\approx log(x)+\frac{a}{x}$

此时，$\log ({\sum }_{i\in [1,n],i\ne j}{e}^{x_i}+e^{x_j})\approx \log (e^{x_j})+\frac{{\sum }_{i\in [1,n],i\ne j}{e}^{x_i}}{e^{x_j}}$

当 $x_j$ 为最大项时，相较于 $log(e^{x_j})$，后半部分可忽略不计

因此存在，log⁡(∑i∈[1,n]exi)=log⁡(∑i∈[1,n],i≠jexi+exj)≈log⁡(exj)=xj=max{x1,x2,...,xn}\log(\sum_{i \in [1,n]}e^{x_i})= \log(\sum_{i \in [1,n],i \not=j}e^{x_i}+e^{x_j}) \approx \log(e^{x_j})=x_j=max\{x_1,x_2,...,x_n\}log(∑i∈[1,n]​exi​)=log(∑i∈[1,n],i=j​exi​+exj​)≈log(exj​)=xj​=max{x1​,x2​,...,xn​}

以此实现 $max$ 函数的平滑化

对于 $min$ 函数的平滑，与 $max$ 函数同理
$$\begin{array}{rl}& min(a_1,a_2,...,a_n)=-max(-a_1,-a_2,...,-a_n)\\ & min(a_1,a_2,...,a_n)=-\log (\sum _{i\in [1,n]}{e}^{-a_i})=-log[\sum _{i\in [1,n]}(e^{\frac{-a_j}{\gamma }}{)}^{\gamma }]\approx -log[\sum _{i\in [1,n]}(e^{\frac{-a_j}{\gamma }}){]}^{\gamma }\\ & \therefore min(a_1,a_2,...,a_n)\approx -\gamma \log \sum _{i=1}^n{e}^{\frac{-a}{r}},\gamma >0\end{array}$$
Soft-DTW的函数计算公式为：
$$\begin{array}{rl}DT{W}_{\gamma }(x,y)& =mi{n}^{\gamma }\{⟨A,△(x,y)⟩,A\in A_{n,m}\}\\ & =-\gamma log(\sum _{A\in A_{n,m}}{e}^{-⟨A,△(x,y)⟩/\gamma })\end{array}$$
在前向传播过程中，定义 An,m⊂{0,1}n×mA_{n,m} \subset \{0,1\}^{n \times m}An,m​⊂{0,1}n×m 为两个长度分别为 $n,m$ 的序列上的校准矩阵集，即对于 路径$A=[a_{i,j}],A\in A_{n,m}$ ，若 $a_{i,j}=1$ ，此时该点 $(i,j)$ 在最优路径上，若 $a_{i,j}=0$，此时点 $(i,j)$ 不在最优路径中。

因此 $A_{n,m}$ 代表从 $(1,1)\to (n,m)$ 的所有路径，且受于DTW相同的约束限制；内积 $⟨A,△(x,y)⟩$ 为路径下的代价和。$△(x,y):=[\sigma (x_i,y_i){]}_{ij}\in {\mathbb{R}}^{n\ast n}$ 为代价矩阵

参考文献

• https://blog.youkuaiyun.com/qq_40206371/article/details/122568648
• https://zhuanlan.zhihu.com/p/346674141
• https://blog.youkuaiyun.com/qq_40206371/article/details/122616951

四、DILATE

文章：Shape and Time Distortion Loss for Training Deep Time Series Forecasting Models

贡献：本文在Soft-DTW算法的基础上，实现了对非平稳信号多步预测问题的处理，使得在面对发生区域剧变以及无法依赖过去信息进行推理的信号时，在兼顾了信号形状吻合的同时，显著降低了时延所带来的损失。

DILATE 算法简介

对于一组长度为x $m$ 的序列 $x$ ，存在 $x=(x_1,x_2,...,x_m)\in {\mathbb{R}}^{p\times m}$。其所对应未来 $n$ 步的预测值为 $\hat{y}=({\hat{y}}_1,{\hat{y}}_2,...,{\hat{y}}_n)\in {\mathbb{R}}^{p\times n}$，未来 $n$ 步序列的真实值为 $y^{\ast }=(y_1^{\ast },y_2^{\ast },...,y_n^{\ast })$，$DILATE$ 目标函数根据超参数 $\alpha \in [0,1]$ 实现对形状项 $L_s$ 和时延项 $L_t$ 的平衡，计算公式如下：
$$L_{DILATE}(\hat{y},y^{\ast })=\alpha L_{shape}(\hat{y},y^{\ast })+(1-\alpha )L_{temporal}(\hat{y},y^{\ast })$$
在上述公式中，$L_{shape}$ 本质为Soft-DTW Loss，即 $L_{shape}=DT{W}_{\gamma }(\hat{y},y^{\ast })$，代价矩阵和校准矩阵集的定义于Soft-DTW算法中定义相同。

$L_{temporal}$ 则是用来约束DTW算法中的时延，距离最优路径越远的点，时延越大。因此 $L_{temporal}$ 通过惩罚这些偏离最优路径较远的匹配，来减少预测过程中的时延。该时延项的损失函数来自时间损失指数 TDI ,实现对路径 $A^{\ast }$ 的时延惩罚，定义如下：
$$TDI(x,y)=⟨A^{\ast },\Omega ⟩=⟨\underset{A\in A_{n,m}}{\mathrm{argmin}}⟨A,△(x,y)⟩,\Omega ⟩$$
通常可以先得到 $DT{W}_{shape}(x,y)$ 中最优路径 $A^{\ast }$，再基于该形状惩罚计算 $⟨A^{\ast },\Omega ⟩$，以获取时延惩罚。

但也可采用类似Soft-DTW算法中对 $DT{W}_{\gamma }(x,y)$ 求解的方法，同时对形状和时延进行惩罚，公式如下：
$$\begin{array}{rl}{L}_{DILATE}(x,y)& =min⟨\alpha \cdot DT{W}_{shape}(x,y)+(1-\alpha )\cdot ⟨A^{\ast },\Omega ⟩⟩\\ & =-\gamma log(\sum _{A\in A_{n,m}}{e}^{-⟨\alpha \cdot ⟨A,△(x,y)⟩+(1-\alpha )\cdot ⟨A^{\ast },\Omega ⟩⟩})\\ & =-\gamma log(\sum _{A\in A_{n,m}}{e}^{-⟨A,\alpha \cdot △(x,y)+(1-\alpha )\cdot \Omega ⟩/r})\end{array}$$

ri,jDILATE=α⋅δi,j+(1−α)⋅ωi,j+minγ{ri,j−1,ri−1,j−1,ri−1,j} r_{i,j}^{DILATE}=\alpha \cdot \delta_{i,j}+(1-\alpha) \cdot \omega_{i,j} + min^\gamma\{r_{i,j-1},r_{i-1,j-1},r_{i-1,j}\} ri,jDILATE​=α⋅δi,j​+(1−α)⋅ωi,j​+minγ{ri,j−1​,ri−1,j−1​,ri−1,j​}

参考文献

• https://zhuanlan.zhihu.com/p/346674141

总结

综上，从DTW到DILATE，上述算法的本质均是对动态规划问题 $r_{i,j}$ 的不断改进
DTW: δi,j⋅wi,j+min{ri,j−1,ri−1,j−1,ri−1,j}Soft−DTW:δi,j⋅wi,j+minγ{ri,j−1,ri−1,j−1,ri−1,j}DILATE: α⋅δi,j+(1−α)⋅ωi,j+minγ{ri,j−1,ri−1,j−1,ri−1,j} \begin{aligned} &DTW: \ \delta_{i,j} \cdot w_{i,j} + min\{r_{i,j-1},r_{i-1,j-1},r_{i-1,j}\}\\ &Soft-DTW: \delta_{i,j} \cdot w_{i,j} + min^{\gamma}\{r_{i,j-1},r_{i-1,j-1},r_{i-1,j}\}\\ &DILATE: \ \alpha \cdot \delta_{i,j}+(1-\alpha) \cdot \omega_{i,j} + min^\gamma\{r_{i,j-1},r_{i-1,j-1},r_{i-1,j}\} \end{aligned} ​DTW: δi,j​⋅wi,j​+min{ri,j−1​,ri−1,j−1​,ri−1,j​}Soft−DTW:δi,j​⋅wi,j​+minγ{ri,j−1​,ri−1,j−1​,ri−1,j​}DILATE: α⋅δi,j​+(1−α)⋅ωi,j​+minγ{ri,j−1​,ri−1,j−1​,ri−1,j​}​
针对原始DTW算法，Soft-DTW通过修改 $mi{n}^{\gamma }$ 使其可微，用于loss函数计算；DILATE则是在Soft-DTW算法的基础上，添加了时延惩罚项，使其尽量减少时延带来的负面影响。