动态模式分解（DMD）

本文地址：https://goodgoodstudy.blog.youkuaiyun.com/article/details/111873026

动态模式分解

动态模式分解（Dynamic Mode Decomposition，DMD）使用随时间增长、衰减和振荡的相干结构来求解或近似动力学系统。将相干结构称为DMD模式。换句话说，DMD将动力学系统转换为模式的叠加，每个模式的强度由特征值控制。

令人惊讶的是，尽管识别DMD模式和特征值的数学过程是纯线性的，但系统本身却可以是非线性的！

有一个合理的理论基础可以证明非线性系统可以用一组模式和特征值对来描述。在这里不多做讨论，感兴趣的话，可以阅读 Koopman 算子及其与DMD的联系¹²³。

DMD不仅是分析系统内部运行情况的有用诊断工具，而且还可以用于预测系统的未来状态。所需要的只是模式和特征值。

正式描述

考虑一组 $n$ 个数据对 { ( x 0 , y 0 ) , ( x 1 , y 1 ) , … , ( x n , y n ) } \{(x_0, y_0), (x_1, y_1), \ldots, (x_n, y_n)\} {(x0​,y0​),(x1​,y1​),…,(xn​,yn​)}，其中 $x_i$ 和 $y_i$ 是长度为 $m$ 的列向量。现在，我们定义两个 $m\times n$ 矩阵：
$$X=\left[\begin{array}{cccc}{x}_0& x_1& \dots & x_n\end{array}\right],Y=\left[\begin{array}{cccc}{y}_0& y_1& \dots & y_n\end{array}\right]$$

如果我们将符号 $A$ 定义为
$$A=Y{X}^{†}$$
其中$X^{†}$是$X$的伪逆⁴，则对 $(X,Y)$ 的动态模式分解由 $A$ 的特征分解给出：即DMD模式和特征值是A的特征向量和特征值。

上面的定义来自 (Tu et al)² ，被称为 exact DMD。它是当前最通用的定义，可以应用于满足给定要求的任何数据集。

在这篇文章中，我们对A满足（或近似满足）方程 $y_i=A{x}_i,\forall i$的情况最感兴趣:
$$Y=AX$$
$X$ 是一组输入向量，$Y$是相应的输出向量集合。 DMD的这种表述非常强大，因为它为分析（和预测）控制方程未知的动力学系统提供了一种便捷的方法。

DMD算法

乍一看，寻找 $A=Y{X}^{†}$ 的特征分解的任务似乎没什么大不了的。确实，当$X$ 和 $Y$ 的大小合理时，从 Numpy 或 MATLAB 调用pinv和eig方法可以解决问题。当 $A\in {\mathbb{R}}^{m\times m}$ 确实很大时，特征分解可能变得麻烦。

幸运的是，借助 exact DMD算法，可以将问题分解为如下几个子任务。

1. 对 $X$ 奇异值分解，可以选择同时执行低阶截断：
   $$X=U\Sigma V^{\ast }$$易知$X^{†}=V{\Sigma }^{-1}{U}^{\ast },A=YV{\Sigma }^{-1}{U}^{\ast }$.

2. 计算 $\tilde{A}$，即将完整矩阵 $A$ 投影到 $U$ 上：
   $$\tilde{A}=U^{\ast }AU=U^{\ast }Y{X}^{†}U=U^{\ast }YV{\Sigma }^{-1}$$

3. 计算 $\tilde{A}$ 的特征值${\lambda }_i$ 和特征向量 $w_i$：
   $$\tilde{A}W=W\Lambda$$

4. 从$W$和 $\Lambda$ 重构 $A$ 的特征分解。 $A$的特征值等于$\tilde{A}$ 的特征值。 $A$ 的特征向量由 $\Phi$ 的列给出：
   $$A\Phi =\Phi \Lambda ，\Phi =YV{\Sigma }^{-1}W$$
   因为
   $$\begin{array}{rl}A\Phi & =YV{\Sigma }^{-1}{U}^{\ast }YV{\Sigma }^{-1}W\\ & =YV{\Sigma }^{-1}\tilde{A}W\\ & =YV{\Sigma }^{-1}W\Lambda \\ & =\Phi A\end{array}$$

在参考文献¹²中可以找到对该算法推导的更深入的解释。以上是 exact DMD模式。

从理论上讲，${\Phi }^{\circ }=UW$ 是 $\Phi$ 的另一推导，也称为 projected DMD 模式：
$$\begin{array}{rl}A{\Phi }^{\circ }& =YV{\Sigma }^{-1}{U}^{\ast }UW\\ & =YV{\Sigma }^{-1}W\\ & =U{U}^{\ast }YV{\Sigma }^{-1}W\\ & =U\tilde{A}W\\ & =UW\Lambda \\ & ={\Phi }^{\circ }\Lambda \end{array}$$

动态系统

本文仅考虑方程式 $Y=AX$ 的两种解释。

第一种解释是 $A$定义了一个差分方程：
$$x_{i+1}=A{x}_i$$
在这种情况下，算子 $A$ 动态系统状态 $x_i$ 在时间上向前迈出了一步。

假设我们有一个时间序列 $D$：
$$D=\left[\begin{array}{cccc}{x}_0& x_1& \dots & x_{n+1}\end{array}\right]$$
其中 $x_i$ 是在时间步 $i$ 处定义系统的 $m$ 维状态的列向量。可以这样定义 $X$ 和 $Y$：
$$X=\left[\begin{array}{ccc}{x}_0& \dots & x_n\end{array}\right],Y=\left[\begin{array}{ccc}{x}_1& \dots & x_{n+1}\end{array}\right]$$
这样，来自$X$和 $Y$ 的每对列向量对应于差分方程式的单次迭代：
$$Y=AX$$
使用DMD，计算特征分解 $A\Phi =\Phi \Lambda$。掌握了DMD模式和特征值，动态系统的状态方程可写成：
$$x_{k+1}=A{x}_k=\Phi \Lambda {\Phi }^{†}{x}_k$$
因此系统状态关于时间步 $k$ 的表达式为：
$$x_k=\Phi {\Lambda }^k{\Phi }^{†}{x}_0$$
若离散时间步长为 $\Delta t$，连续时间 $t$ 中的对应函数为
$$x(t)=\Phi {\Lambda }^{t/\Delta t}{\Phi }^{†}x(0)$$

真正令人惊讶的是，仅使用数据就及时定义了显式函数！这是无方程式建模的一个很好的例子。

本文讨论的 $Y=AX$ 的第二种解释是$A$定义了一个微分方程组：
$$\dot{x}=Ax$$
在这种情况下，矩阵 $X$ 和 $Y$ 将由向量场的 $n$ 个样本组成：$X$ 的第 $i$ 列是位置向量$x_i$； $Y$ 的第 $i$ 列是速度向量 $\dot{x}$。

在计算DMD之后，时间函数与之前差分方程的情形非常相似。如果 $x(0)$ 是任意初始条件，而 $t$ 是连续时间，则
$$x(t)=\Phi \exp (\Lambda t){\Phi }^{†}x(0)$$

参考文献

---

1. Kutz, J. Nathan. Data-driven modeling & scientific computation: methods for complex systems & big data. OUP Oxford, 2013. ↩︎ ↩︎

2. Tu, Jonathan H., et al. “On dynamic mode decomposition: theory and applications.” arXiv preprint arXiv:1312.0041 (2013). ↩︎ ↩︎ ↩︎

3. Wikipedia contributors. “Composition operator.” Wikipedia, The Free Encyclopedia. Wikipedia, The Free Encyclopedia, 28 Dec. 2015. Web. 24 Jul. 2016. ↩︎

4. Wikipedia contributors. “Moore–Penrose pseudoinverse.” Wikipedia, The Free Encyclopedia. Wikipedia, The Free Encyclopedia, 12 May. 2016. Web. 24 Jul. 2016. ↩︎