时间序列中的动态模态分解

动态模态分解 (dynamic mode decomposition) 最早是被用来分析流体（例如水流）的动态过程，它可以把复杂的流动过程分解为低秩的时空特征 (low-rank spatiotemporal features)，这种方法的强大之处在于它不依赖于动态系统中的任何主方程。作为衍生，动态模态分解可以被用来分析多元时间序列 (multivariate time series)，进行短期未来状态预测。

动态模态分解是一种数据驱动的方法，其在描述一些动态过程时具有很多优势，包括：

• 动态模态分解不依赖于任何给定的动态系统表达式；
• 不同于奇异值分解，动态模态分解可以做短期状态预测，即模型本身具备预测能力。

本期文章为大家简要介绍该模型的原理。

测试源码：https://github.com/ekmmrs/DMD_sector_rotation

模型表达式

实际上，动态模态分解和我们比较熟悉的向量自回归一样，他们拥有完全一样的数学表达式。具体而言，若多元时间序列是由 M 条时间长度为 T 的时间序列组成，则对于时刻 t , 动态模态分解的表达式为：

其中，A 表示 Koopman 矩阵，大小为 M x M，当然，在向量自回归里面，我们会称矩阵 A 为系数矩阵 (coefficient matrix)。

在这里，如果令

则动态模态分解的表达式可以写成：

不过与向量自回归不同的是，A 作为动态模态分解中的 Koopman 矩阵时，它可以用一个低秩结构进行逼近。在向量自回归中，如果求解系数矩阵 A ，我们需要对向量自回归的残差平方和做最小化处理，即

模型求解

在动态模态分解中，如果求解 Koopman 矩阵，我们可以采用如下两步：

• 对矩阵 X1 进行奇异值分解，即

• 取矩阵 X1 的截断奇异值分解，截断的秩为 r，则可用如下矩阵：

对 Koopman 矩阵 A 进行近似，其中，矩阵

、

、

分别为 U, V, ∑ 的截断矩阵。

到这里，我们实际上已经完成了对 Koopman 矩阵 A 的近似求解。如果想借助动态模态分解进行时空特征分析，我们可以对矩阵

进行特征值分解，即

，其中，

是一个对角矩阵，对角线上面的元素为特征值，而矩阵

则由特征向量构成。通常来说，我们可以用特征值和特征向量来分析复杂流动过程的时空特征。

实际上，不管是向量自回归还是动态模态分解，它们都具备一定的预测能力。在动态模态分解中，定义

便可以根据

进行短期预测。

实验

01

案例代码

import numpy as np

def dmd(X1, X2, rank):
    """Dynamic Mode Decomposition, DMD."""
    
    u, s, v = np.linalg.svd(X1, full_matrices = 0)
    A_tilde = u[:, : rank].conj().T @ X2 @ v[: rank, :].conj().T @ np.linalg.inv(np.diag(s[: rank]))
    eigval, eigvec = np.linalg.eig(A_tilde)
    Phi = X2 @ v[: rank, :].conj().T @ np.linalg.inv(np.diag(s[: rank])) @ eigvec
    temp = Phi @ np.diag(eigval) @ (np.linalg.pinv(Phi) @ X1)
    
    return temp.real, eigval, Phi

02

实例分析

预测如下矩阵的最后一行，其他5行数据作为训练样本：

X = np.array([[-2,6,1,1,-1],
               [-1,5,1,2,-1],
               [0,4,2,1,-1],
               [1,3,2,2,-1],
               [2,2,3,1,-1],
               [3,1,3,2,-1]])
X1 = X[: 4, :].T
X2 = X[1 : 5, :].T
rank = 3
temp, eigval, Phi = dmd(X1, X2, rank)
print("矩阵X最后一行的预测值为：")
print((Phi @ np.diag(eigval) @ (np.linalg.pinv(Phi) @ X[4, :])).real)

输出结果为：

矩阵X最后一行的预测值为：
[ 3.          1.          3.00000002  1.99999999 -1.        ]

对比结果：