跳到主内容
研究主线时间序列与统计套利 · 第43–50讲7/8模块导览 →
讲义

第49讲:卡尔曼滤波 (Kalman Filters)

从状态空间模型的预测与更新步骤出发,演示目标跟踪、动态均值和时变 Beta 估计及参数敏感性。

完成状态保存在当前浏览器,可随时取消。

一、什么是卡尔曼滤波器?#

卡尔曼滤波器是一种通过噪声观测数据来估计系统状态并进行预测的算法。它通过迭代的预测-更新过程,动态调整对系统状态的估计。该算法特别适用于动态系统中存在噪声的场景。

算法流程

  1. 预测阶段:根据系统模型预测当前状态

  2. 更新阶段:使用新观测值修正预测

核心公式

  • 状态预测方程:xk=Fkx^k1+Bkuk{x}_k = F_k \hat{x}_{k-1} + B_k u_k

  • 协方差预测方程:P_k = F_k P_{k-1} F_k^T + Q_k

  • 卡尔曼增益计算:Kk=PkHkT(HkPkHkT+Rk)1K_k = P_k H_k^T (H_k P_k H_k^T + R_k)^{-1}

  • 状态更新方程:x^k=x^k+Kk(zkHkx^k)\hat{x}_k = \hat{x}_k + K_k (z_k - H_k \hat{x}_k)

  • 协方差更新方程:Pk=(IKkHk)PkP_k = (I - K_k H_k) P_k

核心参数说明#

参数符号参数名称说明
F状态转移矩阵描述系统状态如何随时间变化
H观测矩阵连接系统状态与观测值的矩阵
Q过程噪声协方差矩阵系统模型的不确定性
R观测噪声协方差矩阵测量设备的误差
P估计误差协方差矩阵当前状态估计的可信度

二、实例分析#

2.1 实例1:下落小球跟踪#

物理模型设定

  • 状态变量:位置(x)和速度(v)

  • 状态转移方程:

xk=xk1+vk1τ+0.5gτ2vk=vk1+gτx_{k} = x_{k-1} + v_{k-1}*τ + 0.5*g*τ² v_{k} = v_{k-1} + g*τ
code
import yfinance as yf
import numpy as np
import matplotlib.pyplot as plt
from pykalman import KalmanFilter

# 参数设置
tau = 0.1  # 时间间隔
g = -9.8   # 重力加速度

# 初始化卡尔曼滤波器
kf = KalmanFilter(
    n_dim_obs=1,
    n_dim_state=2,
    initial_state_mean=[30, 10],      # 初始状态估计(位置,速度)
    initial_state_covariance=np.eye(2),
    transition_matrices=[[1, tau], [0, 1]],  # 状态转移矩阵
    observation_matrices=[[1, 0]],           # 观测矩阵(只能观测位置)
    transition_offsets=[0.5*g*tau**2, g*tau], # 重力影响
    observation_covariance=3,        # 观测噪声方差
    transition_covariance=np.zeros((2,2))  # 过程噪声设为0(完美模型)
)

# 生成模拟数据
t = np.arange(40)
real_positions = 30 + 10*t + 0.5*g*tau**2*t**2
noisy_observations = real_positions + np.random.randn(len(t))*3

# 运行滤波器
filtered_means, _ = kf.filter(noisy_observations)

# 结果可视化
plt.figure(figsize=(12,6))
plt.plot(t, real_positions, 'g--', label='真实轨迹')
plt.plot(t, noisy_observations, 'r.', label='噪声观测')
plt.plot(t, filtered_means[:,0], 'b-', label='滤波估计')
plt.legend()
plt.title('下落小球跟踪示例')
plt.xlabel('时间步')
plt.ylabel('高度');

练习:尝试修改观测噪声参数(observation_covariance),观察滤波效果的变化


2.2 实例2:动态移动平均估计#

与传统移动平均对比

卡尔曼滤波无需固定窗口长度,自适应调整权重

code
# 获取股票数据
data = yf.download('AAPL', start='2015-01-01', end='2020-01-01')
prices = data['Close']

# 初始化卡尔曼滤波器
kf = KalmanFilter(
    transition_matrices=[1],
    observation_matrices=[1],
    initial_state_mean=prices[0],
    initial_state_covariance=1,
    observation_covariance=1,
    transition_covariance=0.01
)

# 运行滤波器
state_means, _ = kf.filter(prices.values)
kalman_ma = pd.Series(state_means.flatten(), index=prices.index)

# 计算传统移动平均
ma_30 = prices.rolling(30).mean()
ma_60 = prices.rolling(60).mean()

# 可视化对比
plt.figure(figsize=(12,6))
prices.plot(label='原始价格', alpha=0.5)
kalman_ma.plot(label='卡尔曼滤波', lw=2)
ma_30.plot(label='30日MA', ls='--')
ma_60.plot(label='60日MA', ls='--')
plt.title('移动平均估计对比')
plt.legend();

关键优势

  • 自动调整响应速度

  • 无需人工选择窗口长度

  • 实时更新估计值


2.3 实例3:动态Beta估计#

股票与市场的关系建模

使用线性回归模型:rstock=α+βrmarket+ϵr_{stock} = \alpha + \beta \cdot r_{market} + \epsilon

code
# 获取市场数据
market = yf.download('SPY', start='2015-01-01', end='2020-01-01')['Close'].pct_change().dropna()
stock = yf.download('AAPL', start='2015-01-01', end='2020-01-01')['Close'].pct_change().dropna()

# 对齐数据
common_dates = market.index.intersection(stock.index)
market = market.loc[common_dates]
stock = stock.loc[common_dates]

# 构建观测矩阵(包含市场收益和截距项)
obs_matrix = np.ones((len(market), 1, 2))
obs_matrix[:, 0, 0] = market.values  # 市场收益作为x变量

# 配置卡尔曼滤波器
kf = KalmanFilter(
    n_dim_state=2,
    n_dim_obs=1,
    initial_state_mean=[0, 0],       # 初始alpha和beta
    initial_state_covariance=np.eye(2),
    transition_matrices=np.eye(2),    # 状态保持稳定
    observation_matrices=obs_matrix,
    observation_covariance=0.1,
    transition_covariance=0.01*np.eye(2)
)

# 运行滤波器
state_means, _ = kf.filter(stock.values)

# 可视化时变参数
plt.figure(figsize=(12,6))
plt.subplot(2,1,1)
plt.plot(state_means[:,0], label='动态Alpha')
plt.legend()
plt.subplot(2,1,2)
plt.plot(state_means[:,1], label='动态Beta', color='orange')
plt.legend();

应用场景

  • 动态对冲策略

  • 风险敞口实时监控

  • 因子投资组合管理


三、进阶主题#

3.1 处理非线性系统#

当系统存在非线性关系时,可以使用:

  • 扩展卡尔曼滤波(EKF):通过泰勒展开局部线性化

  • 无迹卡尔曼滤波(UKF):使用sigma点传播技术

code
from pykalman import AdditiveUnscentedKalmanFilter

# 示例:非线性观测模型
def nonlinear_observation(state):
    return np.sin(state[0]) + state[1]**2

ukf = AdditiveUnscentedKalmanFilter(
    transition_functions=lambda x: x,
    observation_functions=nonlinear_observation,
    initial_state_mean=[0, 1],
    initial_state_covariance=np.eye(2)
)

3.2 参数调优#

使用EM算法自动估计参数:

code
kf = kf.em(observations, n_iter=10)  # 进行10次迭代优化

常见问题解答#

Q1:如何选择初始参数?

A1:初始状态可以设为第一个观测值,协方差矩阵初始值较大表示初始不确定性高

Q2:过程噪声和观测噪声如何确定?

A2:可以通过以下方法估算:

  • 过程噪声Q:根据系统物理特性估算

  • 观测噪声R:测量设备的技术参数或历史数据统计

Q3:卡尔曼滤波适用于高频交易吗?

A3:是的,其O(n)时间复杂度适合实时处理,但要注意数值稳定性问题


练习:尝试用卡尔曼滤波实现波动率估计(提示:使用GARCH模型结构)

通过本教程,您应该已经掌握卡尔曼滤波的基本原理和实际应用方法。接下来可以尝试:

  1. 将示例扩展到多资产场景

  2. 结合其他技术指标构建交易策略

  3. 研究不同噪声假设下的滤波器表现