Koopman Operators and DMD: Do Market Modes Survive Out-of-Sample?

动态模式分解为您提供了一个频谱：一些复杂的特征值，每个特征值都有一个增长率和一个频率，每个特征值都附加到资产横截面的空间模式。看起来像结构。整个问题在于它是否是结构，或者是否是一个忠实地记住一个噪声窗口的线性算子。

这个问题有两个可测试的部分，本文是围绕它们构建的：

1. 模式持久性。 将 DMD 安装在窗口 $t$ 和窗口 $t+1$ 上。主要模式是否跨越相同的子空间，或者它们是否会在每次改装时重新洗牌？如果他们重新洗牌，DMD 就是样本内分解，仅此而已 - 显然这比其他教程更有用。
2. 谱半径作为领先指标。 最大特征值模量 $\rho = \max_j |\lambda_j|$ 是一个单一标量，总结了拟合动力学的爆炸性。它是“领先”已实现的波动性、滞后于它还是仅仅重述了它？任一答案均可发布；只有第一个是可以交易的。

这里假设市场是非平稳非线性系统，而不是争论——博客已经争论过这一点，在algotrading中的吸引子中使用相空间几何，并在用HMMs进行政权检测中测量每个政权的BTC统计数据]。接下来是一个更狭隘的问题，即非平稳性对拟合的库普曼算子有何影响，以及如何测量它。

1. 核心思想：线性化非线性动力学

考虑状态空间 $\mathcal{M} \subseteq \mathbb{R}^n$ 上的离散时间动力系统：

$x_{k+1} = F(x_k)$

其中 $F: \mathcal{M} \to \mathcal{M}$ 可能是非线性映射。对于市场而言，$x_k$ 是时间步长 $k$ 的资产回报、波动性或订单簿失衡的向量。

Koopman 运算符 $\mathcal{K}$ 不直接作用于状态 $x$，而是作用于标量值可观察函数 $g: \mathcal{M} \to \mathbb{C}$：

$[\mathcal{K} g](x) = g(F(x))$

关键属性：$\mathcal{K}$ 是线性的，即使 $F$ 不是线性的。代价是维度——$\mathcal{K}$ 作用于无限维函数空间。良好的有限维近似可以为您带来非线性动力学的表现力和线性代数的易处理性：预测变成矩阵求幂，并且每种模式都是可检查的，而不是埋藏在网络权重中。

2. Koopman 算子的谱分解

如果 $\mathcal{K}$ 具有特征值 $\lambda_j$ 和特征函数 $\varphi_j$，则 $\mathcal{K} \varphi_j = \lambda_j \varphi_j$ 以及这些特征函数范围内任何可观测的 $g$ 分解为：

$g(x_k) = \sum_{j=1}^{\infty} \varphi_j(x_0) \, \lambda_j^{k} \, v_j$

其中 $v_j$ 是 Koopman 模式 - 向量值系数，描述每个特征函数如何对完整可观测向量做出贡献。

每个特征值 $\lambda_j = |\lambda_j| e^{i\omega_j}$ 都编码增长或衰减率 ($|\lambda_j|$) 和振荡频率 ($\omega_j$)：

这张桌子就是承诺。第 4 部分是根据数据检查 Promise 的地方。

3. 动态模式分解（DMD）

DMD 是从数据近似 $\mathcal{K}$ 的主力算法。给定排列成矩阵的快照：

$X = \begin{bmatrix} x_0 & x_1 & \cdots & x_{m-1} \end{bmatrix}, \quad X' = \begin{bmatrix} x_1 & x_2 & \cdots & x_m \end{bmatrix}$

DMD 寻求最佳拟合线性算子 $A$ 和 $X' \approx AX$：

1. 计算 SVD：$X = U \Sigma V^*$ 2.项目：$\tilde{A} = U^* X' V \Sigma^{-1}$ 3.特征分解：$\tilde{A} W = W \Lambda$ 4.恢复全空间模式：$\Phi = X' V \Sigma^{-1} W$

$\Phi$的列是DMD模式； $\Lambda$ 的对角线包含 DMD 特征值。

import numpy as np
from numpy.linalg import svd, eig, lstsq

def dmd(X: np.ndarray, rank: int | None = None) -> tuple:
    """
    Dynamic Mode Decomposition.

    Parameters
    ----------
    X : np.ndarray, shape (n_features, n_snapshots)
        Data matrix where each column is a state snapshot.
    rank : int or None
        Truncation rank for the SVD. None = no truncation.

    Returns
    -------
    eigenvalues : np.ndarray, shape (r,)
        DMD eigenvalues (approximating Koopman eigenvalues).
    modes : np.ndarray, shape (n_features, r)
        DMD modes (columns), L2-normalised.
    amplitudes : np.ndarray, shape (r,)
        Mode amplitudes fitted to the FINAL snapshot, so that a one-step
        forecast is simply modes @ (eigenvalues * amplitudes).
    """
    X0 = X[:, :-1]
    X1 = X[:, 1:]

    U, S, Vh = svd(X0, full_matrices=False)

    if rank is not None:
        U = U[:, :rank]
        S = S[:rank]
        Vh = Vh[:rank, :]

    S_inv = np.diag(1.0 / S)

    A_tilde = U.conj().T @ X1 @ Vh.conj().T @ S_inv
    eigenvalues, W = eig(A_tilde)

    modes = X1 @ Vh.conj().T @ S_inv @ W

    norms = np.linalg.norm(modes, axis=0)
    norms[norms == 0] = 1.0
    modes = modes / norms

    amplitudes = lstsq(modes, X[:, -1].astype(complex), rcond=None)[0]

    return eigenvalues, modes, amplitudes

这里有两个经过深思熟虑的选择。模式是 L2 归一化，因为第 4.2 节比较跨窗口的模式子空间，并且非归一化幅度会淹没比较。幅度适合最后快照而不是第一个，这意味着预测永远不会将特征值提高到很大的幂 - 这是数值爆炸的来源，使原始 DMD 信号在浮点溢出时看起来像动态。

4. 测量

这是文章中非教科书的部分。以上均为拟合程序；以下是用于查明适合度是否有意义的协议。

数据。 使用项目自己的交易数据——在一致的分钟或刻度网格上的 BTC、ETH 和流动性替代品的横截面，而不是股票 ETF 的每日下载。该博客是加密优先的，微观结构的论点不会转移。构建 $X$，其中资产位于行、时间位于列、日志返回、每个窗口内的每个资产均被贬低。

4.1 特征值谱

报告，针对代表性窗口：有多少个 $r$ 特征值落在单位圆的公差范围内，每个振荡周期对应于**以小时为单位（周期 $= 2\pi / |\omega_j|$ 条形图，已转换），以及顶部 $r$ 模式重建的返回方差的分数。

def spectrum_report(eigenvalues: np.ndarray, bar_minutes: float,
                    tol: float = 0.05) -> list[dict]:
    """
    Turn a DMD spectrum into human-readable rows: modulus, period in hours,
    and whether the eigenvalue sits on the unit circle within `tol`.
    """
    rows = []
    for lam in eigenvalues:
        modulus = float(np.abs(lam))
        omega = float(np.angle(lam))
        period_hours = (2 * np.pi / abs(omega)) * bar_minutes / 60 if omega else np.inf
        rows.append({
            "modulus": modulus,
            "period_hours": period_hours,
            "on_unit_circle": abs(modulus - 1.0) < tol,
            "regime": "unstable" if modulus > 1 + tol
                      else "persistent" if abs(modulus - 1.0) <= tol
                      else "decaying",
        })
    return sorted(rows, key=lambda r: -r["modulus"])


def reconstruction_r2(X: np.ndarray, modes: np.ndarray,
                      eigenvalues: np.ndarray, amplitudes: np.ndarray) -> float:
    """Fraction of in-window return variance captured by the truncated modes."""
    n_steps = X.shape[1]
    powers = eigenvalues[:, None] ** np.arange(-(n_steps - 1), 1)
    X_hat = (modes @ (amplitudes[:, None] * powers)).real
    resid = np.var(X - X_hat)
    return 1.0 - resid / np.var(X)

一份诚实的频谱报告已经值得这篇文章了。如果单位圆附近没有任何东西，则不存在持续的交易周期，并且股票文献中的“年度季节性”故事根本不会延续到 24/7 加密货币。

4.2 相邻窗口之间的模式稳定性

决定性的考验。在窗口 $t$ 上安装 DMD，然后在窗口 $t+1$ 上安装 DMD，并测量有多少主模子空间幸存。正确的统计量不是模式向量的简单相关性（模式排序和复相位是任意的），而是两个子空间之间的主角度。

def subspace_stability(modes_a: np.ndarray, modes_b: np.ndarray,
                       k: int = 3) -> float:
    """
    Overlap between the leading-k DMD mode subspaces of two adjacent windows.

    Returns the mean cosine of the principal angles: 1.0 = identical subspace,
    0.0 = orthogonal. Immune to mode reordering and complex phase, both of
    which are arbitrary in a DMD fit.
    """
    Qa, _ = np.linalg.qr(modes_a[:, :k])
    Qb, _ = np.linalg.qr(modes_b[:, :k])
    sing = np.linalg.svd(Qa.conj().T @ Qb, compute_uv=False)
    return float(np.mean(np.clip(sing, 0.0, 1.0)))


def stability_curve(returns: np.ndarray, window: int, step: int,
                    rank: int, k: int = 3) -> np.ndarray:
    """Subspace overlap between every pair of adjacent windows."""
    fits = []
    for t_end in range(window, returns.shape[1], step):
        evals, modes, _ = dmd(returns[:, t_end - window:t_end], rank=rank)
        order = np.argsort(-np.abs(evals))
        fits.append(modes[:, order])
    return np.array([subspace_stability(fits[i], fits[i + 1], k=k)
                     for i in range(len(fits) - 1)])

报告这种重叠的分布，并针对空值报告它：对相同回报的相位随机替代项计算的相同统计数据。重叠度较高但不高于代理零值意味着模式正在跟踪协方差结构，而不是动态。

4.3 针对实际波动率的滚动谱半径

要测试的声明：$\rho_t = \max_j |\lambda_j|$，在滚动窗口上重新计算，在“之前”实现波动性而不是随波动性移动。标量在机械上是新的——该博客之前已经实时监控过几何标量，特别是算法交易的复杂流形中的小林曲率——但是库普曼特征值模量是具有不同故障模式的不同量，并且它值得有自己的超前滞后测试，而不是继承的音调。

def rolling_spectral_radius(returns: np.ndarray, window: int = 1440,
                            step: int = 60, rank: int = 5) -> dict:
    """
    Rolling DMD spectrum for regime monitoring.

    returns : np.ndarray, shape (n_assets, n_timesteps)
    window  : rolling window length in bars
    step    : bars between refits
    """
    idx, radii, dom_freq = [], [], []

    for t_end in range(window, returns.shape[1], step):
        X_win = returns[:, t_end - window:t_end]
        try:
            evals, _, _ = dmd(X_win, rank=rank)
        except np.linalg.LinAlgError:
            continue

        idx.append(t_end)
        radii.append(float(np.max(np.abs(evals))))

        on_circle = np.abs(np.abs(evals) - 1.0) < 0.1
        if on_circle.any():
            sel = evals[on_circle]
            dom_freq.append(float(np.abs(np.angle(sel[np.argmax(np.abs(sel))])) / (2 * np.pi)))
        else:
            dom_freq.append(0.0)

    return {"index": np.array(idx),
            "spectral_radius": np.array(radii),
            "dominant_frequency": np.array(dom_freq)}


def lead_lag(signal: np.ndarray, target: np.ndarray, max_lag: int = 24) -> dict:
    """
    Cross-correlation of `signal` against `target` over +/- max_lag steps.
    A peak at negative lag means the signal LEADS the target.
    """
    s = (signal - signal.mean()) / (signal.std() + 1e-12)
    y = (target - target.mean()) / (target.std() + 1e-12)
    lags = np.arange(-max_lag, max_lag + 1)
    corrs = []
    for L in lags:
        if L < 0:
            corrs.append(float(np.corrcoef(s[:L], y[-L:])[0, 1]))
        elif L > 0:
            corrs.append(float(np.corrcoef(s[L:], y[:-L])[0, 1]))
        else:
            corrs.append(float(np.corrcoef(s, y)[0, 1]))
    corrs = np.array(corrs)
    return {"lags": lags, "corr": corrs, "peak_lag": int(lags[np.argmax(np.abs(corrs))])}

将 spectral_radius 与在同一网格上计算的已实现波动率对齐并读取 peak_lag。滞后 0 处的峰值意味着 $\rho_t$ 是具有额外步骤的波动率重述。负滞后峰值在整个样本和各个等级中保持稳定，是本文中唯一具有可交易主张的版本。

关于阳性结果后会发生什么的说明

如果 $\rho_t$ 确实领先，那么明显的下一步行动是横截面构建：根据 DMD 预测的下一步回报对资产进行排名，做多预测的赢家并做空预测的输家。这种结构在这里并不新鲜——它是crypto中的统计套利和配对交易以及向量和矩阵的复杂套利;的复杂套利]第4节中已经涵盖的因子残差交易。唯一真正针对库普曼的变化是特征组合携带时变特征值而不是静态 PCA 载荷。

本文特意省略了策略部分，因为这里尚未对费用和滑点进行回溯测试。当它出现时，它必须清除博客自己的栏：缩减夏普比率和多重测试中的显着性测试，与诚实的否定的常设反例]。 DMD 频谱图不是结果。

对于一步预测本身，正确的实现是根据最后一个快照进行预测，而不是将特征值传播到窗口长度的幂：

def dmd_one_step(returns: np.ndarray, rank: int = 4) -> np.ndarray:
    """
    One-step-ahead prediction from the final snapshot of the window.

    Never raise eigenvalues to the window length: any |lambda| != 1 then
    overflows or underflows and the "signal" becomes numerical garbage.
    """
    evals, modes, amplitudes = dmd(returns, rank=rank)
    return (modes @ (evals * amplitudes)).real

5. 扩展 DMD (EDMD)：非线性可观测量

标准 DMD 对原始状态向量进行操作。 EDMD 首先通过非线性基函数字典提升数据。

给定标量函数 $d_i: \mathbb{R}^n \to \mathbb{R}$ 的字典 $\mathbf{D} = \{d_1, \ldots, d_p\}$，定义提升状态：

$z_k = \begin{bmatrix} d_1(x_k) \\ d_2(x_k) \\ \vdots \\ d_p(x_k) \end{bmatrix} \in \mathbb{R}^p$

EDMD 寻找 $K \in \mathbb{R}^{p \times p}$ 和 $z_{k+1} \approx K z_k$。在下面的代码使用的约定中 - 表示为列，$K$ 作用在左侧：

$K = A G^{-1}, \quad G = \sum_{k=0}^{m-1} z_k z_k^{T}, \quad A = \sum_{k=0}^{m-1} z_{k+1} z_k^{T}$

字典是领域知识进入的地方，而时间延迟行是库普曼关心延迟坐标的特定原因：它们不是附加的单独技术，它们是提升图的又一个块。嵌入本身——延迟 $\tau$ 的选择、嵌入维度 $d$ 及其背后的重构定理——已经在算法交易的复杂流形；从那里获取延迟向量并将它们作为额外行直接输入到 build_financial_dictionary 中。

import numpy as np
from itertools import combinations_with_replacement

def build_financial_dictionary(X: np.ndarray, max_poly_degree: int = 2,
                               include_volatility: bool = True,
                               delay_steps: int = 0) -> np.ndarray:
    """
    Build a dictionary of nonlinear observables for EDMD.

    X : np.ndarray, shape (n_features, n_snapshots)
    Returns Z of shape (n_dict, n_snapshots - delay_steps).
    """
    n_features, n_snapshots = X.shape
    offset = max(delay_steps, 0)
    X_eff = X[:, offset:]
    n_eff = X_eff.shape[1]

    lifted = [X_eff]  # degree-1 terms (identity)

    if max_poly_degree >= 2:
        for deg in range(2, max_poly_degree + 1):
            for combo in combinations_with_replacement(range(n_features), deg):
                term = np.ones(n_eff)
                for idx in combo:
                    term *= X_eff[idx]
                lifted.append(term.reshape(1, -1))

    if include_volatility:
        lifted.append(np.abs(X_eff))   # absolute returns
        lifted.append(X_eff ** 2)      # squared returns

    for d in range(1, delay_steps + 1):
        lifted.append(X[:, offset - d : n_snapshots - d])

    return np.vstack(lifted)


def edmd(X: np.ndarray, dictionary_fn=None, reg: float = 1e-8,
         **dict_kwargs) -> tuple:
    """
    Extended Dynamic Mode Decomposition.

    Solves Z1 ~= K @ Z0 in the least-squares sense. Uses a least-squares
    solve rather than an explicit Gram inverse: `inv` on a near-singular
    dictionary Gram matrix is how EDMD spectra get silently corrupted.
    """
    if dictionary_fn is None:
        dictionary_fn = lambda x: build_financial_dictionary(x, **dict_kwargs)

    Z = dictionary_fn(X)
    Z0, Z1 = Z[:, :-1], Z[:, 1:]

    p = Z0.shape[0]
    G = Z0 @ Z0.T + reg * np.eye(p)   # regularised Gram matrix
    A = Z1 @ Z0.T

    K = np.linalg.solve(G, A.T).T

    eigenvalues, eigenvectors = np.linalg.eig(K)
    return K, eigenvalues, eigenvectors

6. 深度 Koopman 网络

EDMD 字典是手工制作的，当库普曼不变子空间未知时，这是一个真正的限制。深度库普曼网络共同学习提升和算子。

该架构是一个自动编码器——编码器、潜在代码、解码器、重建损失，所有这些都如algotrading中的异常检测]中所介绍的——还有一个补充，这是本节的全部要点：线性损失通过单个学习矩阵 $K$ 强制潜在动态。

x_k  -->  [Encoder φ]  -->  z_k  -->  [Linear K]  -->  z_{k+1}  -->  [Decoder ψ]  -->  x̂_{k+1}
                              |                            |
                              +--- Linearity loss: ‖z_{k+1} - K z_k‖ ---+

$\mathcal{L} = \underbrace{\lVert x_{k+1} - \psi(K \varphi(x_k)) \rVert^2}_{\text{prediction}} + \alpha \underbrace{\lVert \varphi(x_{k+1}) - K \varphi(x_k) \rVert^2}_{\text{linearity}} + \beta \underbrace{\lVert x_k - \psi(\varphi(x_k)) \rVert^2}_{\text{reconstruction}}$

如果没有 $\alpha$ 项，您将拥有一个普通的自动编码器，其潜在空间恰好后面跟着一个矩阵乘法。有了它，网络就会因为任何演化不是线性的潜在表示而受到惩罚——这使得学习到的 $K$ 成为库普曼近似，其特征值与第 4 节的 DMD 谱相当。

import torch
import torch.nn as nn

class DeepKoopman(nn.Module):
    """Deep Koopman autoencoder: learned lifting + linear latent dynamics."""

    def __init__(self, input_dim: int, latent_dim: int, hidden_dim: int = 128):
        super().__init__()

        self.encoder = nn.Sequential(
            nn.Linear(input_dim, hidden_dim), nn.ReLU(),
            nn.Linear(hidden_dim, hidden_dim), nn.ReLU(),
            nn.Linear(hidden_dim, latent_dim),
        )
        self.decoder = nn.Sequential(
            nn.Linear(latent_dim, hidden_dim), nn.ReLU(),
            nn.Linear(hidden_dim, hidden_dim), nn.ReLU(),
            nn.Linear(hidden_dim, input_dim),
        )
        self.K = nn.Linear(latent_dim, latent_dim, bias=False)

    def encode(self, x: torch.Tensor) -> torch.Tensor:
        return self.encoder(x)

    def decode(self, z: torch.Tensor) -> torch.Tensor:
        return self.decoder(z)

    def forward(self, x_k: torch.Tensor) -> dict:
        z_k = self.encode(x_k)
        z_k1_pred = self.K(z_k)
        return {
            "z_k": z_k,
            "z_k1_pred": z_k1_pred,
            "x_k1_pred": self.decode(z_k1_pred),
            "x_k_recon": self.decode(z_k),
        }

    def multi_step_predict(self, x_0: torch.Tensor, n_steps: int) -> torch.Tensor:
        """Roll out by repeated application of the linear operator."""
        z = self.encode(x_0)
        preds = []
        for _ in range(n_steps):
            z = self.K(z)
            preds.append(self.decode(z))
        return torch.stack(preds, dim=1)

    def latent_spectrum(self) -> np.ndarray:
        """Eigenvalues of the learned K — directly comparable to DMD's."""
        return np.linalg.eigvals(self.K.weight.detach().cpu().numpy())


def koopman_loss(model: DeepKoopman, x_k: torch.Tensor, x_k1: torch.Tensor,
                 alpha: float = 1.0, beta: float = 0.5) -> torch.Tensor:
    out = model(x_k)
    z_k1_true = model.encode(x_k1)

    prediction = nn.functional.mse_loss(out["x_k1_pred"], x_k1)
    linearity = nn.functional.mse_loss(out["z_k1_pred"], z_k1_true)
    reconstruction = nn.functional.mse_loss(out["x_k_recon"], x_k)

    return prediction + alpha * linearity + beta * reconstruction

latent_spectrum 使得这个值得训练而不是达到序列模型：学习的算子仍然是一个矩阵，因此第 4.2 节稳定性测试和第 4.3 节超前滞后测试不变地适用于深度模型。

7. 实际考虑因素和陷阱

排名选择。 截断排名 $r$ 是一个偏差-方差刻度盘 - 太低会错过动态，太高会适合噪声。不要关注奇异值弯头；该博客已经正确回答了“在拟合噪声之前有多少个分量”，其中使用向量和矩阵的复杂套利中的随机矩阵理论中的Marchenko-Pastur界限。保留奇异值超过您的窗口形状的 Marchenko-Pastur 边缘的分量，并在您发布的每个频谱旁边明确说明生成的 $r$。

窗口长度。 Koopman 理论假设 $F$ 固定；市场不提供这种产品。滚动改装是强制性的，第 4.2 节的稳定性曲线精确地诊断了所选窗口是否足够长以进行估计并且足够短以保持在一种状态内。

噪声敏感性。 金融数据的信噪比较低，标准 DMD 会受到 $X$ 中噪声的偏差。在得出模式不稳定的结论之前值得尝试的补救措施：

• 总 DMD (TDMD) — 通过总最小二乘法将 $X$ 和 $X'$ 视为噪声。
• 优化 DMD — 直接针对剩余 Frobenius 范数优化特征值模式分解。
• 内核 EDMD — 在高维特征空间中隐式工作，无需构建字典。

如果模式稳定性在 TDMD 下显着提高，则不稳定因素是测量噪声。如果没有，那就是市场了。

8. DMD 的位置

|方法|线性|可解释|多步预测| | ------------ | ----------------------- | ------------------------ | | -------------------- | | DMD |状态空间中的线性 |是（众数 + 特征值）|稳定（矩阵电源）| | EDMD |提升空间中的线性 |是的，根据字典 |稳定（矩阵电源）| |深库夫曼 |学习空间中的线性 |中等（检查潜在 K）|稳定（矩阵电源）|

对于特定波动率的预测，比较点是 GARCH 系列 — 请参阅 crypto 的 GARCH 波动率预测。有关完全非线性序列模型及其可解释性和误差累积权衡，请参阅交易中的时间融合变压器 。

DMD 所占据的利基市场虽然狭窄，但却是真实的：由单一矩阵幂生成的多步预测，而不是自回归推出，并且每种模式都可以检查。该利基市场是否包含 alpha 是第 4 节的问题，而不是此表的问题。

＃＃ 结论

库普曼理论是一种真正优雅的看待市场动态的方式，而优雅正是它需要最严格的测试的原因。要点：

1. 拟合 DMD，然后立即测试拟合情况。 相邻窗口之间的子空间重叠（根据相位随机零值进行测量）会在一个下午内告诉您是否已找到结构或记住了窗口。
2. 滚动谱半径是值得监控的标量，其值完全取决于超前-滞后结果。滞后 0 时，它是波动率代理；负滞后是政权警告。
3. 将振幅锚定在最后一个快照，并且永远不要将特征值提高到窗口长度。野外 DMD“信号”的很大一部分是浮点伪影。
4. 根据 Marchenko-Pastur 选择排名，而不是通过目视肘部，并发布每个频谱的排名。
5. 谱图不是结果。 在此基础上建立的任何策略都必须经受住费用、滑点和收缩的夏普测试，然后才能发挥作用。

对于这些之外的实现，PyDMD 库全面涵盖了 DMD 变体，Mallen 等人的参考代码涵盖了深层的 Koopman 方面。

市场将保持混乱、不稳定和部分可观察的状态。库普曼理论为从混乱中提取结构提供了一个原则性的视角——前提是你下周检查该结构是否仍然存在。

参考文献和进一步阅读：

• B. O. Koopman，“希尔伯特空间中的哈密尔顿系统和变换”，美国国家科学院院刊，1931 年。
• J. H. Tu 等人，“动态模式分解：理论与应用”，计算动力学杂志，2014 年。
• M. O. Williams、I. G. Kevrekidis、C. W. Rowley，“库普曼算子的数据驱动近似：扩展动态模式分解”，非线性科学杂志，2015 年。
• B. Lusch、J. N. Kutz、S. L. Brunton，“非线性动力学通用线性嵌入的深度学习”，《自然通讯》，2018 年。
• J. Mann 和 J. N. Kutz，“金融交易策略的动态模式分解”，定量金融，2016 年。
• A. Mallen 等人，“具有时间分布变化的时间序列的 Koopman 神经预测器”，ICML，2023 年。
• E. Gonzalez 和 M. Generelo，“通过 Koopman 算子、EDMD、Takens 定理和机器学习分析混沌经济模型”，金融和经济中的数据科学，2022 年。