Mean Reversion（OU过程与Kalman滤波）

日期: 2026-07-30 方向: 量化 / 统计套利 / Alpha 阶段: Phase 2 - 统计套利与Alpha Research (Day 89-102) 标签: #量化策略 #均值回归 #OU过程 #Kalman滤波 #Avellaneda

---

今日目标

---

一、理论与模型

1.1 Ornstein-Uhlenbeck 过程

OU 是均值回归的连续时间随机过程：

$$ dX_t = \theta (\mu - X_t) dt + \sigma dW_t $$

• $\theta > 0$：回归速度（mean reversion rate）
• $\mu$：长期均值
• $\sigma$：波动率
• $W_t$：标准布朗运动

离散化（Euler-Maruyama）：

$$ X_{t+\Delta t} = X_t + \theta(\mu - X_t)\Delta t + \sigma \sqrt{\Delta t}, Z_t, \quad Z_t \sim N(0, 1) $$

或更精确的精确离散化：

$$ X_{t+1} = \mu(1 - e^{-\theta \Delta t}) + e^{-\theta \Delta t} X_t + \sigma \sqrt{\frac{1 - e^{-2\theta \Delta t}}{2\theta}} Z_t $$

关键导出量：

• 半衰期：$t_{1/2} = \ln 2 / \theta$
• 长期方差：$\text{Var}_\infty = \sigma^2 / (2\theta)$
• 长期分布：$X_\infty \sim N(\mu, \sigma^2/(2\theta))$

1.2 OU 最大似然估计

观察样本 ${X_0, X_1, \dots, X_n}$（等间隔 $\Delta t$），AR(1) 形式：

$$ X_{t+1} = a + b X_t + \epsilon_t, \quad \epsilon \sim N(0, \tilde\sigma^2) $$

通过 OLS 估计 $\hat a, \hat b, \hat{\tilde\sigma}$，再反推 OU 参数：

$$ \hat\theta = -\frac{\ln \hat b}{\Delta t}, \quad \hat\mu = \frac{\hat a}{1 - \hat b}, \quad \hat\sigma = \hat{\tilde\sigma} \sqrt{\frac{2 \hat\theta}{1 - \hat b^2}} $$

1.3 OU 下的最优交易边界（Avellaneda-Stoikov 思想）

最简单的 Z-score 阈值不是最优的。Bertram (2010) 推导了 OU 下的最优交易边界。考虑无成本无折现，最大化每单位时间收益：

$$ \max_{a, b} \frac{(b - a) - 2c}{T(a, b)} $$

其中 $a, b$ 是入场/出场水平（关于 $\mu$ 对称），$c$ 是单边交易成本，$T(a, b)$ 是平均交易周期。

经验近似（OU 标准化后 $X \to (X-\mu)/\sigma_\infty$）：

• 入场阈值 $\approx 1.0 - 1.5 \sigma_\infty$（远小于经典 2σ）
• 出场阈值 $\approx 0$
• 高成本环境（DEX）入场阈值 ↑ 至 2σ+

1.4 Kalman 滤波动态 hedge ratio

固定 $\beta$ 假设过强。真实市场中 $\beta_t$ 随时间漂移。State-space 模型：

$$ \begin{aligned} \beta_t &= \beta_{t-1} + w_t, \quad w_t \sim N(0, Q) & \text{(状态方程)} \ y_t &= \beta_t x_t + v_t, \quad v_t \sim N(0, R) & \text{(观测方程)} \end{aligned} $$

Kalman 递归（标量版本）：

$$ \begin{aligned} \hat\beta_{t|t-1} &= \hat\beta_{t-1|t-1} \ P_{t|t-1} &= P_{t-1|t-1} + Q \ K_t &= P_{t|t-1} x_t / (x_t^2 P_{t|t-1} + R) \ \hat\beta_{t|t} &= \hat\beta_{t|t-1} + K_t (y_t - \hat\beta_{t|t-1} x_t) \ P_{t|t} &= (1 - K_t x_t) P_{t|t-1} \end{aligned} $$

调参技巧：

• 大 $Q$（process noise）→ $\beta$ 适应快但噪声大
• 大 $R$（observation noise）→ $\beta$ 平滑但滞后
• 实战 $Q/R$ 比通常在 $10^{-5}$ 到 $10^{-3}$

1.5 自适应 Z-score

替代固定 lookback：

$$ Z_t = \frac{X_t - \hat\mu_t}{\hat\sigma_{t,\infty}} $$

其中 $\hat\mu_t, \hat\theta_t, \hat\sigma_t$ 由滚动 OU 校准得到，$\hat\sigma_{t,\infty} = \hat\sigma_t / \sqrt{2\hat\theta_t}$ 是稳态标准差。

---

二、直觉与陷阱

陷阱 1：OU 假设的退化

OU 假设方差稳态有限。但加密在牛熊切换时 $\sigma$ 跳跃，半衰期估计不稳定。 解法：分 regime 校准（HMM 或 GARCH-OU）；或滚动重校准 + 稳定性检验。

陷阱 2：Kalman 的 Q/R 调参陷阱

错误调参 → "动态 β" 退化为对噪声的过拟合。 诊断：

• 看 $\beta_t$ 时序：每天波动 > 5% 必然过拟合
• innovation 序列（$y_t - \hat\beta_{t|t-1} x_t$）应近似白噪声，做 Ljung-Box 检验

陷阱 3：MLE 下偏（small-sample bias）

短样本下 $\hat\theta$ 系统性高估、半衰期低估。 解法：Tang & Chen (2009) 偏差修正；或 Bayesian 估计（先验对 θ 加小约束）。

陷阱 4：Look-ahead bias 的隐形来源

• 用全样本估计 OU 参数生成 in-sample 信号 → 严重 look-ahead
• 用滚动估计但 lookback < 半衰期 × 5 → 估计不稳定
• 阈值用全样本分位数 → look-ahead

解法：必须 expanding 或 rolling，且 lookback ≥ 5× half-life。

陷阱 5：跳跃过程的 OU 误用

加密资产有显著跳跃（清算、闪崩）。OU 是连续过程，跳跃会污染参数估计。 解法：

1. 跳跃过滤：去除 |return| > 5σ 的样本
2. 用 jump-diffusion OU 模型
3. 鲁棒估计：M-estimator、Huber loss

陷阱 6：协整的"假平稳"

短样本可能伪平稳（spurious stationarity）。 解法：变结构断点检验 + 多窗口稳定性检验。

---

三、代码实现

3.1 OU 校准与策略

# mr.py
"""
Mean Reversion Engine
- OU process MLE calibration
- Kalman dynamic hedge ratio
- Adaptive entry/exit thresholds
- Bertram optimal levels
"""

import numpy as np
import pandas as pd
import statsmodels.api as sm
from scipy.optimize import minimize_scalar
from typing import Tuple


# ----------------------------- OU 校准 -----------------------------
def calibrate_ou_mle(x: np.ndarray, dt: float = 1.0) -> dict:
    """
    OU 参数 MLE 估计
    通过 AR(1) OLS 然后转换
    """
    x = np.asarray(x).astype(float)
    x_lag = x[:-1]
    x_next = x[1:]

    # AR(1): x_{t+1} = a + b * x_t + ε
    X = sm.add_constant(x_lag)
    model = sm.OLS(x_next, X).fit()
    a, b = model.params
    sigma_eps = np.sqrt(model.mse_resid)

    # 转换到 OU 参数
    if b >= 1 or b <= 0:
        return {'theta': np.nan, 'mu': np.nan, 'sigma': np.nan,
                'half_life': np.inf, 'sigma_inf': np.nan}

    theta = -np.log(b) / dt
    mu = a / (1 - b)
    sigma = sigma_eps * np.sqrt(2 * theta / (1 - b**2))
    half_life = np.log(2) / theta
    sigma_inf = sigma / np.sqrt(2 * theta)

    return {'theta': theta, 'mu': mu, 'sigma': sigma,
            'half_life': half_life, 'sigma_inf': sigma_inf,
            'r2': model.rsquared}


# ----------------------------- 最优入场（Bertram 近似） -----------------------------
def bertram_optimal_levels(theta: float, mu: float, sigma_inf: float,
                            cost: float) -> Tuple[float, float]:
    """
    Bertram (2010) 最优交易边界（数值近似）
    cost: 单边成本（绝对值，与 spread 同量纲）
    返回 (entry, exit) - entry 是远离 mu 的对称水平，exit = mu
    """
    # 标准化：u = (x - mu) / sigma_inf
    # 目标：max [(b - a) - 2c'] / T(a, b)，其中 c' = cost / sigma_inf
    c_norm = cost / sigma_inf

    def neg_profit_per_time(u: float) -> float:
        # u 是入场水平（标准差倍数），出场水平 = 0
        if u <= c_norm:
            return 1e10
        # 平均交易周期近似：T ≈ (π / θ) * ψ(u)
        # 简化：T ≈ ln(1 + u^2) / θ（粗略近似）
        T = np.log(1 + u**2) / theta + 1e-6
        profit = (u - 2 * c_norm) * sigma_inf  # 来回收益
        return -profit / T

    res = minimize_scalar(neg_profit_per_time, bounds=(c_norm + 0.1, 5.0),
                          method='bounded')
    u_opt = res.x
    entry_level = mu + u_opt * sigma_inf
    return entry_level, mu


# ----------------------------- Kalman 滤波动态 β -----------------------------
class KalmanHedgeRatio:
    """
    动态估计 y_t = β_t * x_t + v_t
    β_t = β_{t-1} + w_t  (随机游走)
    """
    def __init__(self, q: float = 1e-4, r: float = 1e-2,
                 beta_init: float = 1.0, p_init: float = 1.0):
        self.q = q          # process noise variance
        self.r = r          # observation noise variance
        self.beta = beta_init
        self.p = p_init
        self.history = []

    def update(self, y: float, x: float):
        # Predict
        beta_pred = self.beta
        p_pred = self.p + self.q

        # Update
        innovation = y - beta_pred * x
        s = x * x * p_pred + self.r  # innovation variance
        k = p_pred * x / s
        self.beta = beta_pred + k * innovation
        self.p = (1 - k * x) * p_pred

        self.history.append({'beta': self.beta, 'p': self.p,
                             'innovation': innovation, 's': s})
        return self.beta

    def fit_series(self, y: pd.Series, x: pd.Series) -> pd.DataFrame:
        betas, innovs, residuals = [], [], []
        for yi, xi in zip(y.values, x.values):
            self.update(yi, xi)
            betas.append(self.beta)
            innovs.append(self.history[-1]['innovation'])
            residuals.append(yi - self.beta * xi)
        return pd.DataFrame({
            'beta': betas, 'innovation': innovs, 'residual': residuals
        }, index=y.index)


# ----------------------------- 自适应均值回归策略 -----------------------------
def adaptive_mr_strategy(spread: pd.Series, lookback: int = 90,
                          fee_bps: float = 5,
                          min_half_life: float = 1,
                          max_half_life: float = 30) -> pd.DataFrame:
    """
    每天滚动校准 OU，自适应阈值
    """
    fee = fee_bps / 10_000

    out = pd.DataFrame(index=spread.index, columns=[
        'theta', 'mu', 'sigma_inf', 'half_life',
        'z', 'entry', 'signal'])

    position = 0
    for i in range(lookback, len(spread)):
        window = spread.iloc[i - lookback:i].values
        params = calibrate_ou_mle(window)

        if (np.isnan(params['theta']) or
            params['half_life'] < min_half_life or
            params['half_life'] > max_half_life):
            out.iloc[i] = [np.nan]*4 + [np.nan, np.nan, position]
            continue

        # Bertram 最优入场（用 sigma_inf 单位）
        try:
            entry_level, _ = bertram_optimal_levels(
                params['theta'], params['mu'], params['sigma_inf'],
                cost=2 * fee * abs(spread.iloc[i]))
            entry_z = (entry_level - params['mu']) / params['sigma_inf']
        except Exception:
            entry_z = 1.5

        z = (spread.iloc[i] - params['mu']) / params['sigma_inf']

        # 信号
        if position == 0:
            if z > entry_z:
                position = -1
            elif z < -entry_z:
                position = 1
        else:
            if abs(z) < 0.1:
                position = 0
            elif abs(z) > max(entry_z * 2, 4.0):
                position = 0  # 止损

        out.iloc[i] = [params['theta'], params['mu'], params['sigma_inf'],
                       params['half_life'], z, entry_z, position]
    return out


# ----------------------------- Kalman + OU 联合 -----------------------------
def kalman_ou_strategy(y: pd.Series, x: pd.Series,
                        q: float = 1e-4, r: float = 1e-2,
                        ou_lookback: int = 90) -> pd.DataFrame:
    """
    Step 1: Kalman 估动态 β_t
    Step 2: 残差序列上 rolling OU 校准
    Step 3: 自适应交易
    """
    # Step 1
    kf = KalmanHedgeRatio(q=q, r=r, beta_init=1.0)
    kf_result = kf.fit_series(y, x)
    spread = kf_result['residual']

    # Step 2 + 3
    strat = adaptive_mr_strategy(spread, lookback=ou_lookback)
    strat['beta'] = kf_result['beta']
    strat['spread'] = spread
    return strat


# ----------------------------- 回测 -----------------------------
def backtest_mr(y: pd.Series, x: pd.Series, strat: pd.DataFrame,
                fee_bps: float = 5) -> pd.DataFrame:
    rets_y = y.pct_change().fillna(0)
    rets_x = x.pct_change().fillna(0)

    # 用动态 beta（如果有）
    beta_series = strat['beta'] if 'beta' in strat.columns else 1.0
    pair_ret = strat['signal'].shift(1) * (rets_y - beta_series.shift(1) * rets_x)
    turnover = strat['signal'].diff().abs().fillna(0)
    cost = turnover * (fee_bps / 10_000) * 2
    net_ret = pair_ret - cost
    equity = (1 + net_ret).cumprod()
    return pd.DataFrame({'gross': pair_ret, 'cost': cost,
                         'net': net_ret, 'equity': equity})


# ----------------------------- 主流程 -----------------------------
if __name__ == '__main__':
    # 用合成 OU 数据测试
    np.random.seed(42)
    T = 1000
    theta_true = 0.05
    mu_true = 0.0
    sigma_true = 0.02
    x = np.zeros(T)
    for t in range(1, T):
        x[t] = x[t-1] + theta_true * (mu_true - x[t-1]) + sigma_true * np.random.randn()
    x_series = pd.Series(x, index=pd.date_range('2024-01-01', periods=T))

    print("=" * 60)
    print("OU 校准测试（合成数据）")
    params = calibrate_ou_mle(x)
    print(f"  True: theta={theta_true}, mu={mu_true}, sigma={sigma_true}")
    print(f"  Estimated: theta={params['theta']:.4f}, mu={params['mu']:.4f}, "
          f"sigma={params['sigma']:.4f}")
    print(f"  Half-life: {params['half_life']:.2f}")
    print(f"  Sigma_inf: {params['sigma_inf']:.4f}")

    print("\n自适应策略回测")
    strat = adaptive_mr_strategy(x_series, lookback=200)
    print(f"  非空信号占比: {strat['signal'].notna().sum() / len(strat):.2%}")
    print(f"  入场阈值范围: [{strat['entry'].min():.2f}, {strat['entry'].max():.2f}]")

---

四、真实数据/案例

案例 1：USDC 脱钩（2023-03-10 至 2023-03-13）

SVB 倒闭后 USDC 持有 $3.3B 受困，脱钩到 $0.87。

OU 模型分析：

• 历史 USDC/USDT half-life ≈ 4 小时（极快）
• 2023-03-11 偏离 13σ，远超 OU 模型预期
• 入场策略：Z = -3 时买 USDC，但 -10σ 时仍未止损 → 浮亏 8%
• 教训：OU 假设连续，不能处理跳跃；必须有结构性止损

案例 2：3AC 的 GBTC 套利失败（OU 假设破灭）

3AC 在 GBTC 折价 -10% 时大量买入，假设回归到长期均值（历史 +5% 溢价）。 但 GBTC 折价持续扩大到 -45%，OU 校准的 mu 一直在漂移，但 3AC 用的是历史固定均值。

正确做法：

1. 滚动校准 mu（动态 mean）
2. 检验 mu 的结构性漂移（Chow test）
3. 当 mu 漂移幅度 > 2σ_inf 时停止策略

案例 3：基差均值回归（BTC 季度合约）

CME BTC 季度合约相对现货的基差通常 $\theta \approx 0.1$（半衰期约 7 天）。

OU 模型在不同 regime 下参数差异巨大，必须按 regime 校准。

---

五、CEX vs DEX 策略差异

DeFi 原生 MR 机会：

1. Pendle PT 折价收敛：到期前 PT 必须收敛到 1。可用 OU 模型校准折价的回归速度
2. stETH 流动性危机：非线性折价，但事后回归可建模
3. Curve 池失衡：3pool 中某腿失衡时偏离 1:1，回归速度由套利者决定

---

六、风险管理

6.1 OU 模型特有风险

6.2 持仓时间约束

经验法则：最大持仓 = 3 × half_life

例：half_life = 5 days → 5 天内 spread 没回归则强制平仓（即使没到止损）

6.3 凯利公式仓位

OU 下凯利仓位（无成本）：

$$ f^* = \frac{|Z_t| \cdot \theta \cdot \Delta t}{\sigma_\infty^2 \cdot \Delta t} = \frac{|Z_t| \cdot \theta}{\sigma_\infty^2} $$

实际用 0.25-0.5 倍凯利。

---

七、关键速查

OU 参数估计

import statsmodels.api as sm
def ou_params(x, dt=1):
    X = sm.add_constant(x[:-1])
    res = sm.OLS(x[1:], X).fit()
    a, b = res.params
    sigma_eps = np.sqrt(res.mse_resid)
    theta = -np.log(b) / dt
    mu = a / (1 - b)
    sigma = sigma_eps * np.sqrt(2*theta / (1 - b**2))
    return theta, mu, sigma

Kalman Q/R 推荐值

入场阈值建议

---

八、面试题

Q1：OU 模型在加密市场的局限性？如何改进？

答：

1. 跳跃：OU 是连续过程，无法处理清算闪崩。改进：jump-diffusion OU
2. 异方差：σ 随波动率聚集变化。改进：GARCH-OU 或 stochastic volatility OU
3. Regime change：mu 和 theta 在牛熊切换中变化。改进：HMM-OU 或 RSDP（regime-switching）
4. 小样本偏差：theta 高估。改进：Bayesian 估计或偏差校正
5. 静态阈值：忽略时变成本。改进：自适应阈值（Bertram + 实时成本）

Q2：Kalman 滤波的 Q/R 怎么调？

答：

• 理论：Q 是状态过程噪声方差，R 是观测噪声方差。Q/R 比决定平滑度
• 实战调参：
  1. EM 算法极大化 likelihood
  2. 网格搜索 + 样本外 Sharpe
  3. 经验：先用 R = sample_var(residual)，Q = R / 100，再微调
• 诊断：innovation 序列应是白噪声（Ljung-Box test）。若有自相关，Q 调大；若过激（β 跳跃），Q 调小

Q3：Bertram 最优入场为什么不是 2σ？

答：

• 经典 2σ 来自正态分布尾部概率，但忽略了机会成本
• Bertram 框架最大化"每单位时间收益"，需要平衡：
  • 入场太接近 → 单笔利润小
  • 入场太远 → 等待时间长，机会少
• 解析解显示无成本下最优入场约 1σ（无量纲化），有成本时上升至 1.5-2.5σ
• 对加密：CEX 低成本场景 1σ 即可入场；DEX 高成本场景需 2σ+

Q4：你怎么处理 OU 模型的 look-ahead bias？

答：

1. 滚动校准：在 t 时刻只用 [t-L, t] 数据估计参数，禁用全样本
2. Lookback 长度：≥ 5 × 半衰期（半衰期 5 天 → lookback ≥ 25 天，建议 90 天）
3. 阈值滚动：分位数也必须滚动计算，不能用全样本分位数
4. walk-forward：每月重新优化超参数，避免 in-sample 过拟合
5. Out-of-sample 验证：留 30% 数据严格不参与任何调参

Q5：均值回归策略和动量策略可以结合吗？

答：

• 时间尺度分离：短期均值回归 + 中长期动量
• 经典框架：Mean reversion in the short, momentum in the long
• 实施：
  • 用 5min OU 信号做短期反向
  • 用 30day TSMOM 做中期方向
  • 当两者冲突，TSMOM 优先（趋势主导）
• 加密实例：BTC/USD 5min 经常 OU 但日线动量；FTX 暴雷时短期反弹（MR）但中期趋势向下（TSMOM 做空）

---

明日预告

Day 91: Momentum — 横截面动量、TSMOM、AQR 经典框架。今天的 MR 是反向，明天的 momentum 是顺向。两者结合就是完整的多策略组合基础。