← 返回模块
2.3.1.3beta 可读 · 未来付费校验通过内容版本 2026-05-26

ARMA 模型的识别、估计与预测

2.3.1 · 平稳性与 ARMA 模型 · 数学与统计能力

周一早盘，某私募的时间序列研究员把过去 200 个交易日的对冲组合超额收益丢进 statsmodels。她想确认这条曲线是不是一个干净的 ARMA 过程——若是，残差就是一组白噪声，可以挂上下一阶段的 GARCH；若不是，她得回去重做特征工程。问题是：用 AR(1)、MA(1)、ARMA(1, 1) 还是 ARMA(2, 1)？拟合完之后怎么知道这一支模型确实把信息榨干了？再要她出一份「未来三步预测均值 ± 95% 区间」的纪要，她敢报多少？

上一节把 ARMA(p, q) 的结构搭好了，但没有告诉你给定一条数据怎么挑、怎么估、怎么算预测。这一节把 Box-Jenkins 建模法（Box-Jenkins methodology）拆成四步：识别 → 估计 → 诊断 → 预测。每一步都给一个可执行的判据，最后用一段模拟序列把整条流水线走完。

步骤 1：识别（identification）

给定 {xt}t=1T\{x_t\}_{t=1}^T，肉眼无明显趋势、无显著方差变化——视作平稳（stationary，正式的非平稳检验留到下一节）。先画样本自相关函数 ρ^(k)\hat{\rho}(k) 与样本偏自相关函数 ϕ^kk\hat{\phi}_{kk}，k=1,…,mk = 1, \dots, m，mm 一般取 T/4T/4 或 20–40 中较小者。叠加 ±1.96/T\pm 1.96 / \sqrt{T} 的 Bartlett 带：lesson 1 已经证明，在白噪声原假设下样本 ACF 渐近服从 N(0,1/T)N(0, 1/T)，所以单个滞后落在带内即视作「与零不可区分」。日交易频率下 T=200T = 200 给出的带宽约 ±0.139\pm 0.139——也即任何一阶绝对值低于此都没必要建模。

对照上一节的识别表读图：

  • PACF 在 pp 阶截尾、ACF 拖尾 → AR(p)
  • ACF 在 qq 阶截尾、PACF 拖尾 → MA(q)
  • 两者都拖尾 → ARMA(p, q)，且 (p, q) 较小

识别只是「提案」，不是「定案」。常见做法是先列出 2–4 个候选阶数，把最终决定权交给后面的估计与 AIC / BIC。

步骤 2：估计（estimation）

纯自回归（autoregressive）AR(p)：Yule-Walker 矩估计

把 AR(p) 定义式两端同乘 Xt−kX_{t-k} 再取期望——关键观察是 Xt−kX_{t-k} 只依赖于 t−kt-k 时刻之前的创新，而 ϵt\epsilon_t 与那段历史独立，所以 E[ϵtXt−k]=0E[\epsilon_t X_{t-k}] = 0（k≥1k \geq 1），噪声项干净地落地：

  1. 起点：Xt=ϕ1Xt−1+⋯+ϕpXt−p+ϵtX_t = \phi_1 X_{t-1} + \cdots + \phi_p X_{t-p} + \epsilon_t。
  2. 两端乘 Xt−kX_{t-k}、取期望：γ(k)=ϕ1γ(k−1)+⋯+ϕpγ(k−p)\gamma(k) = \phi_1 \gamma(k-1) + \cdots + \phi_p \gamma(k-p)，k=1,…,pk = 1, \dots, p。
  3. 把 pp 条方程堆成矩阵，其中 Γp\Gamma_p 为 p×pp \times p Toeplitz 自协方差矩阵，γp=(γ(1),…,γ(p))T\gamma_p = (\gamma(1), \dots, \gamma(p))^T；再用样本量替换总体量得到 Yule-Walker 估计：
Γp ϕ=γp,ϕ^YW=Γ^p−1 γ^p.\Gamma_p\,\phi = \gamma_p, \qquad \hat{\phi}_{\mathrm{YW}} = \hat{\Gamma}_p^{-1}\,\hat{\gamma}_p.

Toeplitz 结构允许 Durbin-Levinson 递推 O(p2)O(p^2) 求解（按名字引用，不再展开）。ϕ^YW\hat{\phi}_{\mathrm{YW}} 相合且渐近正态分布，渐近协方差为 σ2Γp−1/T\sigma^2 \Gamma_p^{-1} / T，与上模块 2.2.1 lesson 2 给出的 Fisher 信息矩阵之逆一致——所以对 AR(p) 而言 Yule-Walker 与 MLE 渐近等价。

ARMA(p, q)：极大似然估计

当 q≥1q \geq 1 时 Yule-Walker 不再干净（噪声项 θ(L)ϵt\theta(L)\epsilon_t 进入右端）。换用极大似然估计（maximum likelihood estimation, MLE）。设创新 ϵt∼i.i.d.\epsilon_t \overset{\text{i.i.d.}}{\sim} 正态分布（Gaussian distribution）N(0,σ2)N(0, \sigma^2)，把 ϵ1,…,ϵq\epsilon_1, \dots, \epsilon_q 置零（条件似然），可写出条件对数似然：

ℓ(ϕ,θ,σ2;x)=−T2log⁡(2πσ2)−12σ2∑t=max⁡(p,q)+1Tϵt(ϕ,θ)2.\ell(\phi, \theta, \sigma^2; x) = -\tfrac{T}{2}\log(2\pi\sigma^2) - \dfrac{1}{2\sigma^2}\sum_{t = \max(p, q) + 1}^{T} \epsilon_t(\phi, \theta)^2.

ϵt(ϕ,θ)\epsilon_t(\phi, \theta) 由观测值通过逆 ARMA 表示递归重构。一般 MLE 的渐近结果在 2.2.1 lesson 2 已经证过，这里直接引用 Brockwell-Davis 定理 8.8.1：θ^MLE\hat{\theta}_{\mathrm{MLE}} 相合、渐近正态、渐近有效，协方差为 Fisher 信息矩阵之逆。实务上 statsmodels.tsa.arima.model.ARIMA 与 R 的 stats::arima 都给出条件 / 精确似然两种实现；状态空间形式 + 卡尔曼滤波是精确版的标准引擎（仅按名字引用）。

模型选择：AIC / BIC

对每个候选 (p,q)(p, q)，记最大化的对数似然为 ℓ^\hat{\ell}，定义信息准则：

AIC(p,q)=−2ℓ^+2K,BIC(p,q)=−2ℓ^+Klog⁡T,K=p+q+1.\mathrm{AIC}(p, q) = -2\hat{\ell} + 2K, \qquad \mathrm{BIC}(p, q) = -2\hat{\ell} + K\log T, \qquad K = p + q + 1.

挑使准则最小的 (p,q)(p, q)。BIC 的惩罚项 log⁡T\log T 随 TT 增长，所以在中大样本下 BIC 更偏好简约模型；AIC 在某些区间有更好的预测表现，BIC 在真实模型属于候选集时是相合的（T→∞T \to \infty 时以概率 1 选中真实阶数）。两者都是常规操作，工作流就是「拟一个小网格，挑最小值」。实务里常见的做法是同时报告 AIC 与 BIC 两个排序，如果两者选出同一阶数就直接用；不一致时通常以 BIC 为准并把 AIC 的次选作为备选稳健性检验。

步骤 3：诊断（diagnostic checking）

拟合完成后计算标准化残差 e^t=ϵ^t/σ^\hat{e}_t = \hat{\epsilon}_t / \hat{\sigma}。若模型设定正确，{e^t}\{\hat{e}_t\} 应近似白噪声。最常用的整体检验是 Ljung-Box 检验（portmanteau test）：

Q(m)=T(T+2)∑k=1mρ^e(k)2T−k∼χm−p−q2.Q(m) = T(T + 2)\sum_{k = 1}^{m} \dfrac{\hat{\rho}_e(k)^2}{T - k} \sim \chi^2_{m - p - q}.

在显著性水平 α\alpha 下，当且仅当 Q(m)>χm−p−q, α2Q(m) > \chi^2_{m - p - q,\,\alpha} 时拒绝模型设定正确这一原假设。mm 一般取 10 或 20，α=0.05\alpha = 0.05。直观上 Q(m)Q(m) 是前 mm 个残差自相关平方的加权和——若残差里还残留任何被模型遗漏的自相关结构，这一项会膨胀。配合画残差 ACF：若仅有单个滞后超出 Bartlett 带，把 pp 或 qq 加 1 再拟；若多个滞后同向偏离，往往是把信号留在了均值之外的某个维度（譬如条件方差），是切到 GARCH 的信号。

步骤 4：预测（forecasting）

hh 步预测定义为给定 FT=σ(X1,…,XT)\mathcal{F}_T = \sigma(X_1, \dots, X_T) 的条件期望，由 2.1.2 lesson 4 的最小均方误差结论它就是最佳预测：x^T+h∣T=E[XT+h∣FT]\hat{x}_{T+h \mid T} = E[X_{T+h} \mid \mathcal{F}_T]。对中心化常数 μ\mu 的 AR(1)：

x^T+h∣T=μ+ϕ1h(xT−μ),Var(XT+h−x^T+h∣T)=σ2 1−ϕ12h1−ϕ12.\hat{x}_{T + h \mid T} = \mu + \phi_1^h (x_T - \mu), \qquad \mathrm{Var}(X_{T + h} - \hat{x}_{T + h \mid T}) = \sigma^2\,\dfrac{1 - \phi_1^{2 h}}{1 - \phi_1^2}.

h→∞h \to \infty 时 x^T+h∣T→μ\hat{x}_{T+h \mid T} \to \mu，方差 →γ(0)\to \gamma(0)——预测向无条件均值回归，区间宽到无条件标准差。一般 ARMA(p, q) 的 hh 步预测误差方差由 Wold 表示 Xt=∑j≥0ψjϵt−jX_t = \sum_{j \geq 0} \psi_j \epsilon_{t-j}（上一节 MA(∞) 展开）给出：Var=σ2∑j=0h−1ψj2\mathrm{Var} = \sigma^2 \sum_{j=0}^{h-1} \psi_j^2。在正态创新下，95% 预测区间为 x^T+h∣T±1.96Var\hat{x}_{T+h \mid T} \pm 1.96 \sqrt{\mathrm{Var}}。

Formula Explorer

phi^h * x_T

把 ϕ1\phi_1 拨到 0.95、0.6、0.2，再看 h=1,5,20h = 1, 5, 20 的预测衰减——平稳性的可视化在这里最直接：∣ϕ1∣|\phi_1| 越接近 1，回归到 0 越慢。

端到端示例（ARMA(1, 1) 模拟序列）

固定参数 ϕ1=0.6\phi_1 = 0.6、θ1=0.3\theta_1 = 0.3、σ=1\sigma = 1、T=200T = 200，按定义模拟一条序列。

​(a) 识别​​：样本 ACF 在前 3 阶都明显出 Bartlett 带（±0.139\pm 0.139），且缓慢衰减；样本 PACF 在前 2 阶超出后塌回零附近。两者都拖尾，按识别表提案 ARMA(1, 1) 作为主候选，同时把 AR(1)、MA(1)、ARMA(2, 1) 列入对比。

​(b) 估计 + 选阶​​：四个候选拟合后的 AIC 如下表，最小者为 ARMA(1, 1)；MLE 给出 ϕ^1=0.58\hat{\phi}_1 = 0.58、θ^1=0.32\hat{\theta}_1 = 0.32、σ^=1.02\hat{\sigma} = 1.02，与真实参数一致。

候选 (p,q)(p, q)ℓ^\hat{\ell}KKAIC
(1, 0)−289.3-289.32582.6
(0, 1)−291.1-291.12586.2
(1, 1)−285.4-285.43576.8
(2, 1)−285.2-285.24578.4

​(c) 预测​​：取 xT=0.5x_T = 0.5、ϵ^T=0\hat{\epsilon}_T = 0。1 步预测 x^T+1∣T=ϕ^1xT+θ^1ϵ^T=0.58×0.5=0.29\hat{x}_{T+1 \mid T} = \hat{\phi}_1 x_T + \hat{\theta}_1 \hat{\epsilon}_T = 0.58 \times 0.5 = 0.29。MA(∞) 系数 ψ0=1\psi_0 = 1，1 步预测误差方差 =σ^2ψ02=1.04= \hat{\sigma}^2 \psi_0^2 = 1.04，标准差 1.02。95% 预测区间 0.29±1.96×1.02=(−1.71, 2.29)0.29 \pm 1.96 \times 1.02 = (-1.71,\,2.29)。

​(d) 诊断​​：Q(10)=7.8Q(10) = 7.8，自由度 10−1−1=810 - 1 - 1 = 8，χ8, 0.052=15.51\chi^2_{8,\,0.05} = 15.51，p-value≈0.45p\text{-value} \approx 0.45——不拒绝白噪声原假设，残差通过检验。再画残差 ACF 复核：所有滞后均落在 ±0.139\pm 0.139 的 Bartlett 带内，与 Q(10)Q(10) 给出的结论一致，模型已经把可解释的自相关结构榨干。

整条流水线在 statsmodels 里就是 ARIMA(x, order=(1, 0, 1)).fit() + .aic + .forecast(1) + acorr_ljungbox(...) 四行调用，但每一行背后都对应上面的某一步——任何一步的判据没看懂，就不知道该相信报告里的数字到几位小数。

练习

Exercise

你对 T=300T = 300 个月度收益拟合 ARMA(1, 1) 模型，得到 ϕ^1=0.4\hat{\phi}_1 = 0.4、θ^1=−0.2\hat{\theta}_1 = -0.2、σ^=0.018\hat{\sigma} = 0.018，最后一个观测 xT=0.025x_T = 0.025，取 μ=0\mu = 0。

(a) 算出 1 步与 3 步预测 x^T+1∣T\hat{x}_{T+1 \mid T}、x^T+3∣T\hat{x}_{T+3 \mid T}。
(b) 用 MA(∞) 展开 ψ0=1\psi_0 = 1、ψ1=ϕ1+θ1\psi_1 = \phi_1 + \theta_1、ψ2=ϕ1(ϕ1+θ1)\psi_2 = \phi_1(\phi_1 + \theta_1) 近似算出 3 步预测误差标准差。
(c) 报出 95% 的 3 步预测区间。

提示

1 步预测忽略未观测的 ϵT+1\epsilon_{T+1}：x^T+1∣T=ϕ^1xT+θ^1ϵ^T\hat{x}_{T+1 \mid T} = \hat{\phi}_1 x_T + \hat{\theta}_1 \hat{\epsilon}_T。把 ϵ^T\hat{\epsilon}_T 视为 0 只剩 ϕ^1xT\hat{\phi}_1 x_T；3 步预测再迭代两次。

提示

hh 步预测误差方差 =σ^2∑j=0h−1ψj2= \hat{\sigma}^2 \sum_{j=0}^{h-1} \psi_j^2。代入题给 ψ0,ψ1,ψ2\psi_0, \psi_1, \psi_2 算 3 步标准差，再乘 1.96 给出对称的 95% 预测区间。

下一节预告

整条 Box-Jenkins 流水线的前提是「序列已经平稳」。但在量化研究里最常见的反例——价格水平、利率、宏观存量——几乎都不平稳。下一节正面处理这一情况：单位根（unit root）、ADF 检验、差分（differencing）阶数 dd，以及把 ARMA 推广成 ARIMA(p, d, q)。本节的 ARMA 拟合工作流也是模块 2.3.2 中 ARMA + GARCH 联合建模的第一阶段。