自回归模型的数学推导

#AR

自回归模型基础

自回归模型（Autoregressive Model，简称AR模型）是时间序列分析中的核心统计方法，通过历史数据预测当前值。其基本假设是当前时刻的值$X_t$可表示为前$p$个时刻值的线性组合加白噪声误差： $$X_t = c + \sum_{i=1}^p \phi_i X_{t-i} + \varepsilon_t$$ 其中$c$为常数项，$\phi_i$为自回归系数，$\varepsilon_t \sim \mathcal{N}(0,\sigma^2)$为白噪声误差项[citation:2][citation:7]。模型阶数$p$决定了历史窗口长度，直接影响模型预测能力。

数学推导与求解方法

Yule-Walker方程

基于自协方差函数建立参数与序列统计量的关系。定义自协方差函数： $$\gamma(k) = \text{Cov}(X_t, X_{t-k})$$ 推导得Yule-Walker方程组： $$ \begin{pmatrix} \gamma(0) & \gamma(1) & \cdots & \gamma(p-1) \ \gamma(1) & \gamma(0) & \cdots & \gamma(p-2) \ \vdots & \vdots & \ddots & \vdots \ \gamma(p-1) & \gamma(p-2) & \cdots & \gamma(0) \end{pmatrix} \begin{pmatrix} \phi_1 \ \phi_2 \ \vdots \ \phi_p \end{pmatrix} = \begin{pmatrix} \gamma(1) \ \gamma(2) \ \vdots \ \gamma(p) \end{pmatrix} $$ 该方程组将参数估计转化为自协方差矩阵求逆问题[citation:2][citation:7]。

Levinson-Durbin递推算法

为解决高阶模型矩阵求逆的计算效率问题，采用递归求解：

初始化： $$ \phi_{1,1} = \frac{\gamma(1)}{\gamma(0)},\quad \sigma_1^2 = \gamma(0)(1-\phi_{1,1}^2) $$
迭代计算（$k=2$至$p$）： $$ \kappa_k = \frac{\gamma(k) - \sum_{j=1}^{k-1}\phi_{k-1,j}\gamma(k-j)}{\sigma_{k-1}^2} $$ $$ \phi_{k,j} = \phi_{k-1,j} - \kappa_k\phi_{k-1,k-j},\quad 1\leq j\leq k-1 $$ $$ \phi_{k,k} = \kappa_k,\quad \sigma_k^2 = \sigma_{k-1}^2(1-\kappa_k^2) $$ 此算法时间复杂度$O(p^2)$，显著优于直接矩阵求逆的$O(p^3)$，且保证参数估计的稳定性[citation:3]。

平稳性条件

AR模型需满足平稳性要求：特征方程$1 - \phi_1 z - \cdots - \phi_p z^p = 0$的根全在复平面单位圆外。此时序列均值恒定： $$\mu = \frac{c}{1 - \sum_{i=1}^p \phi_i}$$ 方差有限且自协方差仅与时间间隔相关[citation:7]。

参数估计方法

最小二乘估计

通过最小化残差平方和求解： $$\min_{\phi} \sum_{t=p+1}^T \left( X_t - \hat{X}t \right)^2,\quad \hat{X}t = c + \sum{i=1}^p \phi_i X{t-i}$$ 转化为线性回归问题：

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
from sklearn.linear_model import LinearRegression
import numpy as np

# 生成滞后矩阵
X_lag = np.column_stack([X[i:-p+i] for i in range(p)])
y = X[p:]

# 拟合模型
model = LinearRegression().fit(X_lag, y)
params = model.coef_

此方法对异常值敏感但计算高效[citation:7]。

极大似然估计

假设误差服从高斯分布，最大化对数似然函数： $$\mathcal{L}(\phi,\sigma^2) = -\frac{T}{2}\ln(2\pi\sigma^2) - \frac{1}{2\sigma^2}\sum_{t=p+1}^T \varepsilon_t^2$$ 通过梯度下降等优化算法求解：

1
2
3
4
5
6
7
8
9
from scipy.optimize import minimize

def neg_log_likelihood(params, X, p):
    phi = params[:-1]
    sigma2 = params[-1]
    residuals = X[p:] - np.dot(X_lag, phi)
    return 0.5 * len(residuals)*np.log(2*np.pi*sigma2) + np.sum(residuals**2)/(2*sigma2)
    
result = minimize(neg_log_likelihood, initial_guess, args=(X, p))

此方法具有统计最优性但计算复杂[citation:7]。

模型选择与验证

阶数选择准则

信息准则法：
- AIC：$2p - 2\ln(L)$
- BIC：$p\ln(T) - 2\ln(L)$ 平衡模型复杂度与拟合优度[citation:2]
自相关函数分析：
- ACF（自相关函数）：$\rho(k) = \frac{\gamma(k)}{\gamma(0)}$
- PACF（偏自相关函数）：排除中间变量影响后的相关性通过PACF截尾位置确定$p$值[citation:7]

残差诊断

使用Ljung-Box检验验证残差独立性： $$Q = T(T+2)\sum_{k=1}^m \frac{\hat{\rho}_k^2}{T-k}$$ 原假设为残差是白噪声，$p$值>0.05接受原假设[citation:1]

应用场景与实例

视频编码

在AV1标准中，利用AR(24)模型模拟胶片颗粒噪声特性，通过Yule-Walker方程估计参数[citation:4]：

1
2
3
% 伪代码：Film Grain参数估计
autocorr = xcorr(film_grain_samples);
phi = levinson(autocorr, 24);

金融预测

S&P 500指数收益率预测案例：

1
2
3
4
5
6
from statsmodels.tsa.ar_model import AutoReg

# 拟合AR(4)模型
model = AutoReg(returns, lags=4, old_names=False)
results = model.fit()
print(results.params)  # 输出系数估计

融合语言模型

将时间序列转换为文本描述，利用BERT提取特征：

1
2
3
4
5
6
from transformers import AutoTokenizer, AutoModel

text_series = [f"t={t}, return={r:.4f}" for t,r in enumerate(returns)]
inputs = tokenizer(text_series, return_tensors='pt', padding=True)
outputs = model(**inputs)
features = outputs.last_hidden_state.mean(dim=1)  # 时间序列语义特征

结合LSTM进行长期预测，提升金融时间序列预测精度[citation:5]。