Koopman Operators and DMD: Do Market Modes Survive Out-of-Sample?
動態模式分解為您提供了一個頻譜:一些複雜的特徵值,每個特徵值都有一個成長率和一個頻率,每個特徵值都附加到資產橫斷面的空間模式。看起來像結構。整個問題在於它是否是一個結構,或者是否是一個忠實地記住一個雜訊視窗的線性算子。
這個問題有兩個可測試的部分,本文是圍繞它們建構的:
- **模式持久性。 ** 將 DMD 安裝在視窗 和視窗 上。主要模式是否跨越相同的子空間,或者它們是否會在每次改裝時重新洗牌?如果他們重新洗牌,DMD 就是樣本內分解,僅此而已 - 顯然這比其他教程更有用。
- **譜半徑作為領先指標。 ** 最大特徵值模量 是一個單一標量,總結了擬合動力學的爆炸性。它是「領先」已實現的波動性、落後於它還是僅僅重述了它?任一答案均可發布;只有第一個是可以交易的。
這裡假設市場是非平穩非線性系統,而不是爭論——博客已經爭論過這一點,在algotrading中的吸引子中使用相空間幾何,並在用HMMs進行政權檢測中測量每個政權的BTC統計數據]。接下來是一個更狹隘的問題,即非平穩性對擬合的庫普曼算子有何影響,以及如何測量它。
1. 核心思想:線性化非線性動力學

考慮狀態空間 上的離散時間動力系統:
其中 可能是非線性映射。對市場而言, 是時間步長 的資產回報、波動性或訂單簿失衡的向量。
Koopman 運算子 不直接作用於狀態 ,而是作用於標量值可觀察函數 :
關鍵屬性: 是線性的,即使 不是線性的。代價是維度- 作用於無限維函數空間。好的有限維近似可以為您帶來非線性動力學的表現力和線性代數的易處理性:預測變成矩陣求冪,並且每種模式都是可檢查的,而不是埋藏在網路權重中。
2. Koopman 算子的譜分解

如果 具有特徵值 和特徵函數 ,則 以及這些特徵函數範圍內任何可觀測的 分解為:
其中 是 Koopman 模式 - 向量值係數,描述每個特徵函數如何對完整可觀測向量做出貢獻。
每個特徵值 都編碼成長或衰減率 () 和振盪頻率 ():
| 組件 | 特徵值性質 | 財經解讀 |
|---|---|---|
| 趨勢 | 、 | 緩慢漂移,氣勢 |
| 週期 | 、 | 振盪,季節性 |
| 瞬態 | 衰退的衝擊,短暫的走勢 | |
| 不穩定模式 | 不斷增長的爆炸性動態 |
這張桌子就是承諾。第 4 部分是根據數據檢查 Promise 的地方。
3. 動態模式分解(DMD)

DMD 是從資料近似 的主力演算法。給定排列成矩陣的快照:
DMD 尋求最佳擬合線性算子 和 :
- 計算 SVD: 2.項目: 3.特徵分解: 4.恢復全空間模式:
的欄位是DMD模式; 的對角線包含 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 的每日下載。該部落格是加密優先的,微觀結構的論點不會轉移。建構 ,其中資產位於行、時間位於列、日誌返回、每個視窗內的每個資產均被貶低。
4.1 特徵值譜
報告,針對代表性視窗:有多少 特徵值落在單位圓的公差範圍內,每個振盪週期對應於**以小時為單位(週期 條形圖,已轉換),以及頂部 模式重建的返回方差的分數。
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 相鄰視窗之間的模式穩定性
決定性的考驗。在視窗 上安裝 DMD,然後在視窗 上安裝 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 針對實際波動率的滾動譜半徑
要測試的聲明:,在滾動視窗上重新計算,在「之前」實現波動性而不是隨波動性移動。標量在機械上是新的——該博客之前已經實時監控過幾何標量,特別是算法交易的複雜流形中的小林曲率——但是庫普曼特徵值模量是具有不同故障模式的不同量,並且它值得有自己的超前滯後測試,而不是繼承的音調。
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 處的峰值表示 是具有額外步驟的波動率重述。負滯後峰值在整個樣本和各個等級中保持穩定,是本文中唯一具有可交易主張的版本。
關於陽性結果後會發生什麼事的說明
如果 確實領先,那麼明顯的下一步是橫斷面構建:根據 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 首先透過非線性基底函數字典提升資料。
給定標量函數 的字典 ,定義提升狀態:
EDMD 尋找 和 。在下面的程式碼使用的約定中 - 表示為列, 作用在左側:
| 字典類型 | 功能 | 捕捉 |
|---|---|---|
| 多項式 | 非線性跨資產交互作用 | |
| 徑向基 (RBF) | 局部相似性、政權聚類 | |
| 延時嵌入 | 記憶/自迴歸結構 | |
| 傅立葉 | 已知週期(日內、每週) | |
| 波動特徵 | 異方差性,卷聚類 |
字典是領域知識進入的地方,而時間延遲行是庫普曼關心延遲座標的特定原因:它們不是附加的單獨技術,它們是提升圖的另一個區塊。嵌入本身——延遲 的選擇、嵌入維度 及其背後的重構定理——已經在演算法交易的複雜流形;從那裡獲取延遲向量並將它們作為額外行直接輸入到 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中的異常檢測]中所介紹的——還有一個補充,這是本節的全部要點:線性損失通過單個學習矩陣 強制潛在動態。
x_k --> [Encoder φ] --> z_k --> [Linear K] --> z_{k+1} --> [Decoder ψ] --> x̂_{k+1}
| |
+--- Linearity loss: ‖z_{k+1} - K z_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. 實際考慮因素與陷阱

**排名選擇。 ** 截斷排名 是偏差-方差刻度盤 - 太低會錯過動態,太高會適合噪音。不要關注奇異值彎頭;該部落格已經正確回答了“在擬合噪聲之前有多少個分量”,其中使用向量和矩陣的複雜套利中的隨機矩陣理論中的Marchenko-Pastur界限。保留奇異值超過您的視窗形狀的 Marchenko-Pastur 邊緣的分量,並在您發布的每個頻譜旁邊明確說明生成的 。
**視窗長度。 ** Koopman 理論假設 固定;市場不提供這種產品。滾動改裝是強制性的,第 4.2 節的穩定性曲線精確地診斷了所選視窗是否足夠長以進行估計並且足夠短以保持在一種狀態內。
**噪音敏感性。 ** 金融資料的訊噪比較低,標準 DMD 會受到 中雜訊的偏差。在得出模式不穩定的結論之前值得嘗試的補救措施:
- 總 DMD (TDMD) — 透過總最小平方法將 和 視為雜訊。
- 最佳化 DMD — 直接針對剩餘 Frobenius 範數最佳化特徵值模式分解。
- 核心 EDMD — 在高維度特徵空間中隱式工作,無需建立字典。
如果模式穩定性在 TDMD 下顯著提高,則不穩定因素是測量雜訊。如果沒有,那就是市場了。
8. DMD 的位置

|方法|線性|可解釋|多步驟預測| | ------------ | ----------------------- | ------------------------ | | -------------------- | | DMD |狀態空間中的線性 |是(眾數 + 特徵值)|穩定(矩陣電源)| | EDMD |提升空間中的線性 |是的,根據字典 |穩定(矩陣電源)| |深庫夫曼 |學習空間中的線性 |中(檢查潛在 K)|穩定(矩陣電源)|
特定波動率的預測,比較點是 GARCH 系列 — 請參閱 crypto 的 GARCH 波動率預測。有關完全非線性序列模型及其可解釋性和誤差累積權衡,請參閱交易中的時間融合變壓器 。
DMD 所佔據的利基市場雖然狹窄,但卻是真實的:由單一矩陣冪產生的多步驟預測,而不是自回歸推出,並且每種模式都可以檢查。該利基市場是否包含 alpha 是第 4 節的問題,而不是此表的問題。
## 結論

庫普曼理論是一種真正優雅的看待市場動態的方式,而優雅正是它需要最嚴格的測試的原因。要點:
- **擬合 DMD,然後立即測試擬合情況。 ** 相鄰視窗之間的子空間重疊(根據相位隨機零值進行測量)會在一個下午內告訴您是否已找到結構或記住了視窗。
- 滾動譜半徑是值得監控的標量,其值完全取決於超前-滯後結果。滯後 0 時,它是波動率代理;負滯後是政權警告。
- 將振幅錨定在最後一個快照,並且永遠不要將特徵值提高到視窗長度。野外 DMD「訊號」的很大一部分是浮點偽影。
- 根據 Marchenko-Pastur 選擇排名,而不是透過目視肘部,並發布每個頻譜的排名。
- **譜圖不是結果。 ** 在此基礎上建立的任何策略都必須經受住費用、滑點和收縮的夏普測試,然後才能發揮作用。
對於這些之外的實現,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 年。
Authors
Trading-systems engineer
Trading-systems engineer building bots since 2017: cross-exchange arbitrage (connected up to 30 venues), cointegration-based pairs arbitrage across spot and futures, scalping, news and sentiment-driven strategies, trend algorithms, and portfolio management and balancing algorithms. Also builds sub-millisecond order execution, big-data warehouses, backtesting engines, AI agents, and trading interfaces (incl. open-source profitmaker.cc). Stack: JS/TS, Python, Rust/Zig/Go, DevOps, backend, frontend, architecture.