时间序列分析（浙大备考版）大纲

时间序列分析——总结笔记

本文基于《应用时间序列分析》教材整理，适用于备考浙江大学时间序列分析（3学分，课程号06121291）。基于不同教材的时间序列分析侧重点不同，考点也不同，请勿完全参考。

大纲如下：

1. 时间序列基础
   • 平稳序列与部分典型平稳序列
   • Hilbert空间
   • 谱函数与谱密度
2. 三大模型
   • AR(p)模型与Yule Walker方程
   • MA(q)模型
   • ARMA(p, q)模型
3. 参数估计与假设检验
   • 均值、自协方差函数的估计
   • 白噪声检验
   • 时间序列的预报
   • ARMA模型参数的估计
4. 三大模型的相关计算（另开一篇）

由于是复习笔记，所以内容较为简略。

目录

• 时间序列分析——总结笔记
  • Part 1：时间序列基础
    • 1-1 时间序列与平稳序列
    • 1-2 典型平稳序列
    • 1-3 Hilbert空间
    • 1-4 谱函数与谱密度
  • Part 2：三大模型
    • 2-1 常系数线性差分方程理论
    • 2-2 AR(p)模型
    • 2-3 Yule-Walker方程与偏相关系数
    • 2-4 MA(q)模型
    • 2-5 ARMA(p, q)模型
  • Part 3：参数估计与假设检验
    • 3-1 均值与自协方差函数的估计
    • 3-2 白噪声检验
    • 3-3 最佳线性预测与投影
    • 3-4 时间序列的线性预测
    • 3-5 AR模型参数估计
    • 3-6 MA模型的参数估计
    • 3-7 ARMA模型的参数估计

Part 1：时间序列基础

1-1 时间序列与平稳序列

时间序列是随机变量列，课程涉及的时间序列均为离散时间指标的。平稳序列是特殊的时间序列，满足的基本条件是：

1. \(\forall t,\mathbb E(X_t^2)<\infty\)，即存在二阶矩。
2. \(\forall t,\mathbb E(X_t)=\mu\)，即均值平稳。
3. \(\forall t,s,\mathbb E[(X_t-\mu)(X_s-\mu)]=\gamma_{t-s}\)，即相关性结构平稳。

平稳序列的自协方差函数仅与时间差有关，按照时间差排列构成的实数列称为自协方差函数：\(\gamma_0,\gamma_1,\cdots\)。自协方差函数具有如下性质：

1. 对称性，\(\gamma_k=\gamma_{-k}\)。
2. 非负定性，即自协方差矩阵非负定。
3. 有界性，\(\forall k,|\gamma_k|\le \gamma_0\)。

自协方差矩阵具有重要作用，它指的是

\[\Gamma_n=(\gamma_{|i-j|})_{n\times n}=\begin{bmatrix} \gamma_0 & \gamma_1 & \cdots & \gamma_{n-1} \\ \gamma_1 & \gamma_0 & \cdots & \gamma_{n-2} \\ \vdots & \vdots & \ddots & \vdots \\ \gamma_{n-1} & \gamma_{n-2} & \cdots & \gamma_0 \end{bmatrix}. \]

其下标为矩阵的阶数，它是一个非负定矩阵，且绝大多数时候是正定的。特别\(\det(\Gamma_n)=0\)时，\(X_0,X_1,\cdots,X_{n-1}\)是线性相关的。因此，如果\(\Gamma_{n}\)正定而\(\det(\Gamma_{n+1})=0\)，则\(X_n\)可以被\(X_0,\cdots,X_{n-1}\)线性表示。

• 如果\(k\to \infty\)时\(\gamma_k\to 0\)，则自协方差矩阵任意阶正定。（反证法）

自协方差函数的有界性，使得\(\gamma_k/\gamma_0\in[-1,1]\)，将所有自协方差函数除以方差\(\gamma_0\)得到自相关函数：\(\rho_1,\rho_2,\cdots\)，显然\(\rho_0=1\)。

最为简单的平稳序列是白噪声序列，满足：

1. \(\mathbb E(\varepsilon_t)=\mu\)，均值为常数。
2. \(\mathbb D(\varepsilon_t)=\mathbb E(\varepsilon_t^2)=\sigma^2\)，二阶矩为有限常数。
3. \(\mathbb{Cov}(\varepsilon_t,\varepsilon_s)=\delta_{t-s}\sigma^2\)，序列不相关。

平稳序列之间具有两种特殊的关系：正交、不相关。

• 正交：\(\forall s, t,\mathbb E(X_sY_t)=0\)。
• 不相关：\(\forall s,t,\mathbb E(X_sY_t)=\mu_X\mu_Y\)。
• 如果任一平稳序列是零均值的，则正交与不相关等价。

正交（不相关）平稳序列的和仍然是平稳序列。如果\(\{X_t\}\)和\(\{Y_s\}\)是正交的零均值平稳序列，则\(\mu_Z=0\)，\(\gamma_Z(k)=\gamma_X(k)+\gamma_X(k)\)。

1-2 典型平稳序列

白噪声的滑动和：无穷数列\(\{a_j\}\)平方可和；\(\{\varepsilon_t\}\)是零均值白噪声序列\({\rm WN}(0,\sigma^2)\)，定义白噪声的滑动和是

\[X_t=\sum_{j=-\infty}^\infty a_j\varepsilon_{t-j}. \]

当\(\{a_k\}\)平方可和时（保证了随机变量的收敛性），\(\{X_t\}\)是一个零平稳序列，自协方差函数为

\[\gamma_k=\sigma^2\sum_{j=-\infty}^\infty a_ja_{j+k}. \]

它随时间差\(k\)衰减的，即\(\gamma_k\to 0\)。

由滑动平均导出以下两种重要的平稳序列：有限滑动平均、单边无穷滑动和。它们是具有现实意义的，因为未来的白噪声不可观测。

有限滑动平均：\(\{a_j\}\)是一个有限数列，除了\(k\in[0,q]\)都有\(a_k=0\)，即

\[X_t=\sum_{j=0}^pa_j\varepsilon_{t-j}. \]

其自协方差函数与白噪声的滑动平均类似，有

\[\gamma_k=\sigma^2\sum_{j=0}^{q-k}a_ja_{j+k},\quad 0\le k\le q. \]

如果\(k>q\)，则\(\gamma_k=0\)。

单边无穷滑动和：\(k<0\)都有\(a_k=0\)，即

\[X_t=\sum_{j=0}^\infty a_j\varepsilon_{t-j}. \]

在实际应用中，我们可能经常假定白噪声服从正态分布，以便于求出随机变量的分布。

线性滤波是一般平稳序列的线性滑动和，要求\(\{h_j\}\)平方可和，\(X_t\)是零均值列，则

\[Y_t=\sum_{j=-\infty}^\infty h_jX_{t-j}. \]

如果\(\{X_t\}\)的自协方差函数是\(\gamma_k\)，则\(Y_t\)的自协方差函数为

\[\gamma_Y(k)=\sum_{i,j=-\infty}^\infty h_ih_j\gamma_{k+i-j}. \]

严平稳序列并不是平稳序列，它指的是联合分布函数与时间无关的时间序列，但二阶矩有限的严平稳序列一定是平稳的。严平稳序列具有遍历性，可以从时间序列的一次实现推出任意有限维分布的性质。

如果\(\{X_t\}\)严平稳遍历，则对任何有限维多元函数\(\varphi(\cdot)\)，\(\varphi(X_1,\cdots,X_m)\)严平稳遍历。特别地，白噪声的滑动和是严平稳遍历的。

1-3 Hilbert空间

Hilbert空间指完备的内积空间，平稳序列总是定义在Hilbert空间上的，这部分内容是理论准备。

内积空间首先是一个向量空间，在其中定义内积后成为内积空间，如果度量也成立，成为距离空间。而内积空间的完备性，指的是柯西序列总收敛于空间中一点。

定义\(L^2\)为所有二阶矩存在的随机变量构成的空间，内积\(\langle X,Y\rangle\)定义为\(\mathbb E(XY)\)，正交即\(\mathbb E(XY)=0\)。范数定义为\(\|X\|=\sqrt{{\mathbb E}(X^2)}\)，\(X\)和\(Y\)的距离定义为\(\|X-Y\|\)。在这样的定义下\(L^2\)是一个Hilbert空间，此时的收敛指均方收敛。

平稳序列\(\{X_t\}\)的所有有限线性组合构成的空间是\(L^2(X)\)，这个空间不是完备的，而所有二阶矩存在的随机变量构成的空间\(L^2\)是完备的。因此通常在线性闭包下研究平稳序列，即

\[\bar L^2(X), \]

它是\(L^2\)中包含\(L^2(X)\)的最小闭子空间。

Hilbert空间中存在如下的定理：

• 平行四边形公式：

  \[\|X+Y\|^2+\|X-Y\|^2=2(\|X\|^2+\|Y\|^2). \]

• Cauchy-Schwarz不等式：

  \[\langle X,Y\rangle^2\le\|X\|^2\|Y\|^2. \]

• 勾股定理：

  \[\langle X,Y\rangle=0\Leftrightarrow \|X\|^2+\|Y\|^2=\|X+Y\|^2. \]

• 投影定理：设\(M\)是Hilbert空间\(H\)的一个闭线性子空间，则\(\forall X\in H\)，\(\exists !Y\in M\)，使得\(Z=X-Y\)满足

  \[\langle X,Z\rangle =0, \]

  且\(\forall Y'\in M\)，\(\|X-Y\|<\|X-Y'\|\)。\(Y\)称为\(X\)在子空间\(M\)上的投影。

最小序列：设\(\{X_t\}\)是零均值平稳序列，\(H_x\)是其生成的Hilbert空间\(\bar L^2(X)\)，令

\[H_x(s)=\{X_t:t\ne s\}, \]

如果\(\exists s,H_x\ne H_x(s)\)，则称\(\{X_t\}\)是最小序列。由平稳性，\(\forall s,H_x\ne H_x(s)\)，也就是说最小序列里的每一项都是重要的。

• 最小序列比自协方差矩阵正定的要求更强。

• 如果\(\{X_t\}\)有谱密度，则\(\{X_t\}\)是最小序列等价于

  \[\int_{-\pi}^\pi\frac{{\rm d}\lambda}{f(\lambda)}<\infty. \]

1-4 谱函数与谱密度

复随机变量：\(Z=X+{\rm i}Y\)，\(X,Y\)是实随机变量。

• \(\mathbb E(Z)=\mathbb E(X)+{\rm i}\mathbb E(Y)\)。
• \(\mathbb E(Z^2)=\mathbb E(X^2)+\mathbb E(Y^2)\)。
• \(\mathbb{Cov}(X,Y)=\mathbb{E}[(X-\mu_X)\overline{(Y-\mu_Y)}]\)。

复时间序列：由复随机变量构成的随机变量列称为复时间序列。一个典型的复时间序列：\(\epsilon_t=e^{{\rm i}tY},Y\sim U(-\pi,\pi)\)，它的期望为\(\mathbb{E}(\epsilon_t)=\delta_t\)，自协方差函数为\(\mathbb{Cov}(\epsilon_t\epsilon_s)=\delta_{t-s}\)，其中关联着一个重要积分：

\[\int_{-\pi}^\pi e^{{\rm i}tx}{\rm d}x=2\pi\delta_t. \]

平稳序列的谱函数是唯一存在的（Herglotz定理），关联着它的自协方差函数。定义是某个函数\(F(\lambda)\)，满足

\[\gamma_k=\int_{-\pi}^\pi e^{{\rm i}k\lambda}{\rm d}F(\lambda),\quad F(-\pi)=0. \]

在谱函数可微的情况下，谱密度是谱函数的导函数。我们所感兴趣的谱函数大多是可微的，如果\(F'(\lambda)=f(\lambda)\)，则\(f(\lambda)\)就是平稳序列的谱密度，即满足下式：

\[\gamma_k=\int_{-\pi}^\pi e^{{\rm i}k\lambda}f(\lambda){\rm d}\lambda. \]

反演公式：如果平稳序列\(\{X_t\}\)的自协方差函数\(\{\gamma_k\}\)绝对可和，则谱密度可以表示成

\[f(\lambda)=\frac1{2\pi}\sum_{k=-\infty}^\infty \gamma_ke^{-{\rm i}k\lambda}. \]

• 如果\(X_t\)谱密度存在，则其自协方差矩阵任意阶正定。

无穷滑动和具有谱密度：白噪声\({\rm WN}(0,\sigma^2)\)关于数列\(\{a_j\}\)的无穷滑动和\(\{X_t\}\)谱密度为

\[f(\lambda)=\frac{\sigma^2}{2\pi}\left|\sum_{j=-\infty}^\infty a_je^{-{\rm i}j\lambda} \right|^2. \]

• 这说明无穷滑动和的自协方差矩阵总是任意阶正定的。

• 特别地，白噪声的谱密度为\(f(\lambda)=\sigma^2/2\pi\)。

如果平稳序列\(\{X_t\}\)关于数列\(\{h_j\}\)的线性滤波\(\{Y_t\}\)具有谱密度\(f_X(\lambda)\)，则\(Y_t\)的谱密度是

\[f_Y(\lambda)=\left|\sum_{j=-\infty}^\infty h_je^{-{\rm i}j\lambda} \right|^2f_X(\lambda). \]

正交零均值平稳序列的谱函数（谱密度）之和，就是和的谱函数（谱密度）。

如果谱函数是阶梯型的，此时的平稳序列称为离散谱序列（类比离散型随机变量与其分布函数）。

离散谱序列举例：调和平稳序列。设随机变量\(\xi\)和\(\eta\)满足\(\mathbb E(\xi)=\mathbb E(\eta)=0\)，\(\mathbb E(\xi\eta)=0\)，\(\mathbb E(\xi^2)=\mathbb E(\eta^2)=\sigma^2\)，\(\lambda_0\in(0,\pi]\)。定义时间序列为

\[Z_t=\xi\cos(t\lambda_0)+\eta\sin(t\lambda_0)\xlongequal{def}A\cos(t\lambda_0-\theta),\\ A\xlongequal{def}\sqrt{\xi^2+\eta^2},\cos\theta=\frac{\xi}{A},\sin\theta=\frac{\eta}{A}. \]

\(Z_t\)是一个平稳序列，因为\(\mathbb E(Z_t)=0\)，\(\mathbb {Cov}(Z_t,Z_s)=\sigma^2\cos((t-s)\lambda_0)\)。计算得\(Z_t\)的谱密度为

\[F(\lambda)=\frac{\sigma^2}2[I_{[-\lambda_0,\pi]}(\lambda)+I_{[\lambda_0,\pi]}(\lambda)],\quad (\lambda_0\ne \pi). \]

这个函数是右连续的，在\(-\lambda_0,\lambda_0\)处发生了两次阶梯跳跃。同时，这个函数的自协方差函数不收敛到0，是周期函数。若干个调和平稳序列的叠加仍然是一个离散谱序列。

• 对于离散谱序列，如果其谱函数恰有\(n\)个跳跃点，则\(\Gamma_n\)正定，\(\Gamma_{n+1}\)退化。如果\(F(\lambda)\)有无穷多个跳跃点，则自协方差矩阵任意阶正定。

调和模型：尽管离散谱序列是平稳序列，但其每一次观测是确定的三角函数的简单叠加，所以在实际工作中常常被人们视作非随机序列处理。

谱密度与角频率：如果平稳序列的谱密度存在，且在某一点\(\lambda_0\)有一个明显的峰值，则平稳序列应该有一个以\(\lambda_0\)为角频率的频率成分；这个峰值约陡峭，相应的频率成分就越重要。

Part 2：三大模型

2-1 常系数线性差分方程理论

引入推移算子以方便表达，即\(\mathscr B(X_t)=X_{t-1}\)，\(\mathscr B^n(X_t)=X_{t-n}\)。推移算子对非时间序列随机变量无效果，即\(\mathscr B(Y)=Y\)。需要注意的是推移算子与多项式的结合（它具有重要的作用）：

\[A(z)=\sum_{j=-\infty}^\infty a_jz^j\Rightarrow A(\mathscr B)X_t=\sum_{j=-\infty}^\infty a_jX_{t-j}. \]

多项式的无穷滑动和与谱密度可以用多项式来表达：

\[X_t=\sum_{j=-\infty}^\infty a_j\varepsilon_{t-j}=A(\mathscr B)\varepsilon_t,\\ f_X(\lambda)=\frac{\sigma^2}{2\pi}|A(e^{-{\rm i}\lambda})|^2. \]

常系数齐次线性差分方程的求解：设常系数线性差分方程为\(A(\mathscr B)X_t=0\)，这里限定

\[A(z)=1-\sum_{j=1}^pa_jz^j=\prod_{j=1}^k(1-z_j^{-1}z)^{r(j)}, \]

这里\(z_j\)是\(A(z)\)的根，\(r(j)\)是其重数，则其有通解

\[X_t=\sum_{j=1}^k\sum_{l=0}^{r(j)-1}U_{l,j}t^lz_j^{-t}. \]

通解的形式是我们日后讨论三大模型的基础，如果\(z_j\)全部（严格）大于1，则随着\(t\)增大，\(X_t\to 0\)，称为以负指数阶收敛到0。如果有\(|z_j|=1\)，则存在周期解；如果有\(|z_j|<1\)，则存在爆炸解。

由线性方程组的理论，如果\(A(\mathscr B)X_t=0\)的右端不是0，就称为非齐次常系数线性差分方程。非齐次形式下，方程的通解是对应齐次通解加上非齐次下的一个特解。如果对应的齐次方程通解以负指数阶收敛到0，那么方程的解随着时间推移，受到通解的影响就会越来越小。

平稳解指的是满足非齐次线性差分方程组的平稳序列解，其核心是平稳序列。注意到通解都不是平稳解，如果\(z_j\)全部严格大于1，那么差分方程的解随着时间推移，只与平稳解差一个无穷小量，可以忽略。

2-2 AR(p)模型

AR(p)模型：\(\{\varepsilon_t\}\)是白噪声\({\rm WN}(0,\sigma^2)\)，现有一个零点都在单位圆外的\(p\)阶多项式\(A(z)\)，常数项为1，AR(p)模型如下定义：

\[A(\mathscr B)X_t=\varepsilon_t.\\ A(z)=1-\sum_{j=1}^pa_jz^j,\\ \forall |z|\le 1,A(z)\ne 0. \]

又叫自回归模型，相当于用时间序列的历史信息，对将来作回归预测。由于是用时间序列自身的信息作回归，所以称为自回归模型。相应的平稳解叫做AR(p)序列。

模型的自回归系数为

\[\boldsymbol a=(a_1,\cdots,a_p)', \]

\(A(z)\)称为特征多项式，特征多项式的零点都在单位圆外，这一条件称为最小相位条件（稳定性条件）。只要满足稳定性条件，齐次通解就会很快收敛于0，所以只用考虑特解。

AR(p)模型平稳解的求解：AR(p)模型只有一个平稳解，为

\[X_t=A^{-1}(\mathscr B)\varepsilon_t=\sum_{j=0}^\infty\psi_j\varepsilon_{t-j}. \]

注意\(A^{-1}(z)\)在\([-1,1]\)上解析，故可以泰勒展开，得到的系数\(\{\psi_j\}\)称为Wold系数，Wold系数是以负指数阶收敛到0的。

另外，这是白噪声的单边无穷滑动和，自然得到其自协方差函数和谱密度为

\[\gamma_k=\sigma^2\sum_{j=0}^\infty\psi_j\psi_{j+k}, \\ f(\lambda)=\frac{\sigma^2}{2\pi|A(e^{-{\rm i}\lambda})|^2}. \]

其自协方差函数具有短记忆性，是以负指数阶收敛到0的，这比一般的无穷滑动和更快。

2-3 Yule-Walker方程与偏相关系数

AR(p)序列是单边无穷滑动和，定义\(\boldsymbol \gamma_n=(\gamma_1,\cdots,\gamma_n)'\)，\(\forall j>p,a_j=0\)，\(\boldsymbol a_n=(a_1,\cdots,a_p,a_{p+1},\cdots,a_n)'\)，则由

\[\begin{bmatrix} X_t \\ X_{t+1} \\ \vdots \\ X_{t+n-1} \end{bmatrix}=\begin{bmatrix} X_{t-1} & X_{t-2} & \cdots & X_{t-n} \\ X_{t} & X_{t-1} & \cdots & X_{t+1-n} \\ \vdots & \vdots & \ddots & \vdots \\ X_{t+n-2} & X_{t+n-3} & \cdots & X_{t-1} \end{bmatrix}\begin{bmatrix} a_1 \\ a_2 \\ \vdots \\ a_n \end{bmatrix}+\begin{bmatrix} \varepsilon_t \\ \varepsilon_{t+1} \\ \vdots \\ \varepsilon_{t+n-1} \end{bmatrix}. \]

两边同时乘以\(X_{t-1}\)并取期望，就得到

\[\boldsymbol \gamma_n=\Gamma_n\boldsymbol a_n,\quad n\ge p. \]

由于AR(p)序列是白噪声的单边滑动和，所以自协方差矩阵可逆，于是

\[\boldsymbol a_n=\Gamma_n^{-1}\boldsymbol \gamma_n. \]

这是Yule-Walker方程中的第一个式子，能够得到自回归系数；为了得到白噪声方差，考察\(\gamma_0\)，有

\[\gamma_0=\mathbb E(X_t^2)=\mathbb E[(A(\mathscr B)\varepsilon_t)^2]=\boldsymbol a_n'\Gamma_n\boldsymbol a_n+\sigma^2. \]

于是

\[\sigma^2=\gamma_0-\boldsymbol a_n'\boldsymbol \gamma_n=\gamma_0-\boldsymbol \gamma_n'\boldsymbol a_n. \]

以上两个式子为AR(p)模型的矩估计打下了基础。

对于一般的（不一定是AR(p)的）平稳序列\(X_t\)，如果其自协方差矩阵正定，则也可以类似地解出一组自回归系数：

\[\boldsymbol a_n=\Gamma_n^{-1}\boldsymbol \gamma_{n}. \]

当阶数为\(n\)时，记\(\boldsymbol a_n=(a_{n,1},a_{n,2},\cdots,a_{n,n})'\)，则称\(\{X_t\}\)的偏相关系数为

\[\{a_{n,n}\}. \]

偏相关系数的p后截尾：如果存在某个\(p\)，使得\(n>p\)时总有\(a_{n,n}=0\)，在\(n=p\)时的\(a_{p,p}\ne 0\)，则称这个平稳序列偏相关系数p后截尾。如果一个零均值平稳序列的偏相关系数p后截尾，则它一定是AR(p)序列。

随着问题规模的增大，偏相关系数不好求，所以有Levinson递推公式用于求解偏相关系数：

\[a_{1,1}=\gamma_1/\gamma_0,\\ \sigma_0^2=\gamma_0, \\ \sigma_k^2=\sigma_{k-1}^2(1-a_{k,k}^2), \\ a_{k+1,k+1}=\frac{\gamma_{k+1}-\gamma_ka_{k,1}-\gamma_{k-1}a_{k,2}-\cdots-\gamma_1\gamma_{k,k}}{\gamma_0-\gamma_1a_{k,1}-\gamma_2a_{k,2}-\cdots-\gamma_ka_{k,k}},\\ a_{k+1,j}=a_{k,j}-a_{k+1,k+1}a_{k,k+1-j},1\le j\le k. \]

这里\(\boldsymbol X_k=(X_k,X_{k+1},\cdots,X_1)\)，

\[\sigma^2_k\xlongequal{def}\mathbb E(X_{k+1}-\boldsymbol a_k'\boldsymbol X_k)^2. \]

定理4.1内容是：\(\boldsymbol a_n=(a_{n,1},\cdots,a_{n,n})\)满足最小相位条件，故把它们视为AR(n)模型的自回归系数，得到的模型确实是自回归模型。

对于AR(p)模型，其最佳线性预测是

\[\hat X_{p+1}=\boldsymbol a_p'\boldsymbol X_p. \]

如果不能确定平稳序列是AR(p)模型，用偏相关系数\(\boldsymbol a_{n}\)拟定一个AR(n)模型并类似用\(\boldsymbol a_n'\boldsymbol X_n\)是合理的，因为其满足最小相位条件。

2-4 MA(q)模型

MA(q)模型：\(\{\varepsilon_t\}\)是\({\rm WN}(0,\sigma^2)\)，\(B(z)\)在单位圆内没有根（注意与AR(p)模型不同），且\(B(z)\)的常数项为1，则MA(q)模型如下定义：

\[X_t=B(\mathscr B)\varepsilon_t,\\ B(z)=1+\sum_{j=1}^qb_jz^j,\\ \forall |z|<1,B(z)\ne 0. \]

又叫滑动平均模型，是白噪声的有限滑动和，故比AR(p)模型简单的多。相应的平稳解叫做MA(q)序列。

如果进一步要求\(B(z)\)在单位圆上也没有根，则相应的平稳解称为可逆的MA(q)序列。可逆的MA(q)序列，白噪声序列可以用\(X_t\)的历史信息表达。

MA(q)序列属于白噪声的滑动和，所以自协方差函数和谱密度为

\[\gamma_k=\sum_{j=0}^{q-k}a_ja_{j+k},0\le k\le q.\\ f(\lambda)=\frac{\sigma^2}{2\pi}|B(e^{-{\rm i}\lambda})|^2. \]

当\(k>q\)时\(\gamma_k=0\)，这个性质说明MA(q)序列的自协方差函数是q后截尾的，并有引理保证，自协方差函数q后截尾的平稳序列一定是MA(q)序列。

MA(q)序列与AR(p)序列的谱密度在除去常数项后互为倒数。MA(q)序列的谱密度倒数可积，所以MA(q)序列是最小序列。

MA(q)的模型参数可以如此计算：

\[A=\begin{bmatrix} 0 & 1 & 0 & \cdots & 0 & 0\\ 0 &0 & 1 & \cdots & 0&0 \\ \vdots& \vdots & \vdots &\ddots & \vdots &\vdots\\ 0&0&0&\cdots &0&1 \\ 0&0&0&\cdots &0&0 \end{bmatrix}_{q\times q},C=\begin{bmatrix} 1 \\ 0 \\ \vdots \\ 0 \end{bmatrix}_{q\times 1}, \\ \Omega_k=\begin{bmatrix} \gamma_1 & \gamma_2 & \cdots & \gamma_k \\ \gamma_2 & \gamma_3 & \cdots & \gamma_{k+1} \\ \vdots & \vdots & \ddots & \vdots \\ \gamma_q & \gamma_{q+1} & \cdots & \gamma_{k+q-1} \end{bmatrix}_{q\times k},\boldsymbol \gamma_q=\begin{bmatrix} \gamma_1 \\ \gamma_2 \\ \vdots \\\gamma_q \end{bmatrix}_{q\times 1},\\ \Pi=\lim_{\to \infty}\Omega_k\Gamma_k^{-1}\Omega_k', \]

有

\[\sigma^2=\gamma_0-C'\Pi C,\quad \boldsymbol b_q=\frac1{\sigma^2}(\boldsymbol \gamma_q-A\Pi C). \]

以上两个式子为基于观测样本估计MA(q)的模型系数打下了基础。

2-5 ARMA(p, q)模型

将AR(p)模型与MA(q)模型结合，就得到ARMA(p, q)模型：

\[A(\mathscr B)X_t=B(\mathscr B)\varepsilon_t,\\ A(z)=1-\sum_{j=1}^pa_jz^j,\quad \forall |z|\le 1,A(z)\ne 0,\\ B(z)=1+\sum_{j=1}^qb_jz^j,\quad \forall |z|<1,B(z)\ne 0. \]

另外还要求两个多项式没有公共根。由ARMA(p, q)模型决定的平稳解称为ARMA(p, q)序列，它是唯一的。

如果MA(q)模型可逆，即\(B(z)\)在单位圆上也没有根，则对应的ARMA序列称为可逆ARMA(p, q)序列。

同理，有

\[X_t=A^{-1}(\mathscr B)B(\mathscr B)\varepsilon_t\xlongequal{def}\Phi(\mathscr B)\varepsilon_t=\sum_{j=0}^\infty \psi_j\varepsilon_{t-j}. \]

这里的\(\psi_j\)也称为Wold系数，规定\(j>q\)时\(b_j=0\)，则递推计算公式为

\[\psi_j=0,\quad j<0,\\ \psi_0=1,\\ \psi_j=b_j+\sum_{k=1}^pa_j\psi_{j-k},\quad j>0. \]

ARMA(p, q)模型的可识别性：模型参数\((\boldsymbol a',\boldsymbol b',\sigma^2)\)可以由平稳解的自协方差函数唯一确定，这是要求\(A(z)\)和\(B(z)\)没有公共根的结果。

• 引理：如果\(A(\mathscr B)X_t=B(\mathscr B)\varepsilon_t\)的平稳解\(\{X_t\}\)同时对另一个组多项式\(C(z), D(z)\)和白噪声\(\{\eta_t\}\)，满足
  \[C(\mathscr{B})X_t=D(\mathscr{B})\varepsilon_t, \]
  则\(C(z)\)的阶数\(\ge p\)，\(D(z)\)的阶数\(\ge q\)。

在此引理的帮助下，如果要求模型参数，可以先用Yule-Walker方程，在确定阶数的情况下确定AR(p)部分的模型参数，扣掉这部分得到MA(q)部分，再用MA(q)的参数估计方式得到另一部分模型参数。

ARMA(p, q)序列也是白噪声的无穷滑动和，其自协方差函数与AR(p)模型类似：

\[\gamma_k=\sigma^2\sum_{j=0}^{\infty}\psi_j\psi_{j+k},\\ f(\lambda)=\frac{\sigma^2}{2\pi}\left|\frac{B(e^{-{\rm i}\lambda})}{A(e^{-{\rm i}\lambda})} \right|^2. \]

这种谱密度的形式称为有理谱密度，因为\(A(z)\)和\(B(z)\)不具有公因子，所以中间的分式是不可约的。

ARIMA(p, d, q)模型：求和ARMA(p, q)模型，\(Y_t\)是ARIMA(p, d, q)序列等价于其d阶差分是ARMA(p, q)模型，即

\[(1-\mathscr{B})^dY_t=X_t,\quad A(\mathscr{B})X_t=B(\mathscr{B})\varepsilon_t.\\ A(z)=1-\sum_{j=1}^pa_jz^j,\quad \forall |z|\le 1,A(z)\ne 0,\\ B(z)=1+\sum_{j=1}^qb_jz^j,\quad \forall |z|<1,B(z)\ne 0. \]

同时，也要求\(A(z)\)和\(B(z)\)没有公共根。要注意的是，ARIMA序列不是平稳的。

单位根过程：当\(d=1\)时，ARIMA(p, 1, d)模型又称为单位根模型，这是因为\(A(z)(1-z)\)具有一个单位根\(z=1\)，因而不是平稳的。

Part 3：参数估计与假设检验

3-1 均值与自协方差函数的估计

平稳序列\(\{X_t\}\)的观测值为\(x_1,\cdots,x_N\)，\(N\)是观测值数量，样本均值\(\mu\)的点估计是

\[\hat \mu=\bar x_N=\frac 1N\sum_{i=1}^Nx_i. \]

接下来用\(\bar X_N\)代替此统计量，研究均值估计的性质。\(\bar X_N\)具有以下性质：

1. \(\bar X_N\)是\(\mu\)的无偏估计。
2. 如果\(\gamma_k\to 0\)，则\(\bar X_N\)是\(\mu\)的相合估计，且至少是均方相合。
3. 如果\(\{X_t\}\)是严平稳遍历的，则\(\bar X_N\)还是\(\mu\)的强相合估计。

对于线性平稳序列，其样本均值还服从中心极限定理：即对于

\[X_t=\sum_{j=-\infty}^{\infty}\psi_jX_{t-j}+\mu, \\ f(\lambda) = \frac{\sigma^2}{2\pi}\left|\sum_{j=-\infty}^\infty\psi_je^{-{\rm i}j\lambda} \right|^2, \]

只要\(f(\lambda)\)在\(\lambda=0\)处联系续，且\(f(0)\ne 0\)，就有

\[\sqrt{N}(\bar X_N-\mu)\stackrel{d}\to N(0, 2\pi\sqrt{0}). \]

只要\(\{\psi_j\}\)绝对可和，就能满足谱密度在\(\lambda =0\)处连续的条件，如果再追加非零的条件，就满足此中心极限定理，依分布收敛于正态分布。另外，由于AR(p)、MA(q)、ARMA(p, q)都是系数绝对可和的，所以都可以依据CLT构造\(\bar X_N\)的区间估计。

\(\bar X_N\)的收敛速度在服从重对数律时，阶数一般是不能再被改进的

\[O\left(\sqrt{\frac{2\ln\ln N}{N}} \right). \]

这是因为重对数律有如下的约束：如果\(f(0)\ne 0\)的线性平稳序列满足\(\psi_{|k|}\)以负指数阶收敛于0，或者\(f(0)\)处连续且\(\mathbb{E}|\epsilon_t|^r<\infty\)对某个\(r>2\)成立，就有

\[\lim_{N\to \infty}\sup\sqrt{\frac{N}{2\ln\ln N}}(\bar X_N-\mu)=\sqrt{2\pi f(0)},{\rm a.s.} \]

由此，\(\{X_t\}\)和\(\{-X_t\}\)必定同时满足重对数律的成立条件，但是它们的上下界将收敛于不同的值\(\pm \sqrt{2\pi f(0)}\)，因此不可能再改进阶数。

对于自协方差函数\(\gamma_k\)的估计，一般采用

\[\hat\gamma_k=\frac1N\sum_{j=1}^{N-k}(x_j-\bar x_N)(x_{j+k}-\bar x_N),\quad 0\le k \le N-1,\\ \hat\gamma_{-k}=\hat\gamma_k,\quad 1\le k\le N-1. \]

这是因为我们希望使得自协方差函数收敛于0的速度更快，因而分母不用\(N-k\)。另一方面的原因是这样定义的协方差矩阵\(\hat \Gamma_k\)是正定的（只要\(x_1,\cdots, x_N\)不全相同）。

对于自相关函数的估计，一般采用

\[\hat\rho_k=\frac{\hat\gamma_k}{\hat\gamma_0},\quad k\le |N-1|. \]

关于自协方差函数估计量的性质，有

1. 如果\(\gamma_k\to 0\)，则对每一个确定的\(k\)，\(\hat\gamma_k\)是\(\gamma_k\)的渐进无偏估计。
2. 如果\(\{X_t\}\)是严平稳遍历的，则对每个确定的\(k\)，\(\hat\gamma_k,\hat\rho_k\)分别是\(\gamma_k, \rho_k\)的强相合估计。

可以看到，如果要使\(\bar X_N\)和\(\hat \gamma_k, \hat \rho_k\)都是强相合估计量，需要\(\{X_t\}\)是严平稳遍历序列，这是因为严平稳遍历序列的一次观测值可以推知任意维联合分布的信息。而线性平稳序列总是白噪声的无穷滑动和，因而AR、MA、ARMA模型总可以用\(\bar X_N\)和\(\hat \gamma_k,\hat\rho_k\)来估计\(\mu\)和\(\gamma_k,\rho_k\)。

对于线性平稳序列，其协方差函数还服从中心极限定理，对于

\[X_t=\sum_{j=-\infty}^\infty \psi_j\varepsilon_{t-j},\\ \gamma_k=\sigma^2\sum_{j=-\infty}^\infty \psi_j\psi_{j+k},\\ f(\lambda)=\frac{\sigma^2}{2\pi}\left|\sum_{j=-\infty}^\infty \psi_je^{-{\rm i}j\lambda} \right|^2. \]

只要谱密度\(f(\lambda)\)平方可积，即

\[\int_{-\pi}^\pi f^2(\lambda){\rm d}\lambda<\infty \]

这等价于自协方差函数平方可和，就有

\[\sqrt{N}(\hat\gamma_0-\gamma_0,\cdots,\hat\gamma_h-\gamma_h)\stackrel {d}\to (\xi_0,\cdots,\xi_h), \\ \sqrt{N}(\hat\rho_1-\rho_1,\cdots,\hat\rho_h-\rho_h)\stackrel{d}\to (R_1,\cdots,R_h). \]

其中，定义\(\{W_t\}\)是独立同分布标准正态白噪声，\(\mu_4=\mathbb E(\varepsilon_t^4)\)，有

\[\xi_j=(M_0\lambda_0)W_0+\sum_{t=1}^\infty (\gamma_{t-j}+\gamma_{t+j})W_t,\quad j\ge 0,\\ R_j=\sum_{t=1}^\infty (\rho_{t+j}+\rho_{t-j}-2\rho_t\rho_j)W_t,\quad j\ge 1.\\ M_0=\frac1{\sigma^2}(\mu_4-\sigma^4)^{1/2}. \]

由于AR、MA、ARMA序列的自协方差函数总是以负指数阶收敛到零，所以必定是平方可和的，于是\(\hat \gamma_k\)服从这样的中心极限定理，可以构造区间估计。但计算不便，这部分结论可以跳过。

3-2 白噪声检验

Box-Pierce白噪声检验基于自协方差函数的中心极限定理中，将线性平稳序列换成白噪声得到的结果：如果\(\{X_t\}\)是独立同分布的白噪声，样本自相关系数定义为

\[\hat\rho_k=\frac{\sum_{t=1}^{N-k}(x_t-\bar x_N)(x_{t+k}-\bar x_N)}{\sum_{t=1}^{N-k}(x_t-\bar x_N)^2}, \]

注意这里分母的项数，则对任意\(h\)，

\[\sqrt{N}(\hat\rho_1,\cdots,\hat\rho_h)\stackrel{d}\to N_h(0,I_h). \]

所以对充分大的\(N\)，\(\sqrt{N}\hat\rho_1,\cdots,\sqrt{N}\hat\rho_m\)近似独立服从\(m\)维标准正态分布，所以定义卡方统计量为

\[\chi^2(m)\xlongequal {def}N(\hat\rho_1^2+\cdots+\hat\rho_m^2)\stackrel {d}\to \chi^2_m.,\quad m\le \sqrt {N}. \]

如果\(\{X_t\}\)是白噪声序列，则理论上\(\rho_1^2+\cdots+\rho_m^2=0\)，所以当\(\chi^2(m)\)较大的时候应拒绝原假设，即找到\(\chi^2_m\)的上\(\alpha\)分位数\(\lambda_{m, \alpha}\)，也就是\(\mathbb{P}(\chi^2_m>\lambda_{m, \alpha})=\alpha\)，拒绝域为

\[\chi^2(m) \ge \lambda_{m, \alpha}. \]

白噪声的简化检验法：构造\(Q\)统计量

\[Q(m)=\frac1m\#\{j|\sqrt{N}|\hat\rho_j|\ge 1.96,1\le j\le m \}. \]

同样不能把\(m\)取得过大。如果\(Q(m)\ge 0.05\)，就拒绝\(\{X_t\}\)是白噪声这一假设。

3-3 最佳线性预测与投影

（零均值的）最佳线性预测：用\(\boldsymbol{X}=(X_1,\cdots,X_n)'\)的线性组合对\(Y\)进行预测，使得预测的均方误差最小，即

\[\arg\min_\boldsymbol{a=(a_1,\cdots,a_n)'}\mathbb{E}(Y-\boldsymbol a'\boldsymbol{X}). \]

如果\(\boldsymbol {X}\)的自协方差函数为\(\Gamma\)，则预测系数\(\boldsymbol a\)满足预测方程：

\[\Gamma\boldsymbol{a} = \mathbb{E}(\boldsymbol{X}Y), \]

预测方程总是有解的，并且如果\(\Gamma\)可逆，预测系数就是唯一的\(\boldsymbol a=\Gamma^{-1}\mathbb{E}(\boldsymbol{X}Y)\)。但不管\(\Gamma\)是否可逆，最佳线性预测总是\({\rm a.s.}\)唯一的。

最佳线性预测实质就是Hilbert空间上的投影，这里\(Y\in H\)，\(X_1,\cdots,X_n\in \bar L^2(X)=L^2(X)\)（由于\(\boldsymbol X\)是有限维的），将\(Y\)关于\(\boldsymbol{X}\)的最佳线性预测记作

\[L(Y|\boldsymbol{X})=P_{L^2(X)}({Y}). \]

结合Hilbert空间上的相关定理，得到最佳线性预测的如下性质：

1. 如果\(\mathbb{E}(\boldsymbol{X}Y)=\boldsymbol{0}\)，则\(L(Y|\boldsymbol{X})=0\)。这指的是\(Y\)与\(L^2(X)\)正交的情况下投影为0。

2. 如果\(Y\)是\(X_1,\cdots,X_n\)的线性组合，则\(L(Y|\boldsymbol{X})=Y\)。这指的是\(Y\in L^2(X)\)，故投影是其自身。

3. 投影是线性算子，这等价于

   \[Y=\sum_{j=1}^m b_jY_j\rightarrow L(Y|\boldsymbol{X})=\sum_{j=1}^m b_jL(Y_j|\boldsymbol{X}). \]

4. 如果\(\boldsymbol{X}=(X_1,\cdots,X_n)\)与\(\boldsymbol{Y}=(Y_1,\cdots,Y_m)\)满足\(\mathbb{E}(\boldsymbol{X}'\boldsymbol{Y})=\boldsymbol{O}\)，则

   \[L(Z|\boldsymbol{X},\boldsymbol{Y})=L(Z|\boldsymbol{X})+L(Z|\boldsymbol{Y}). \]

   这表明\(\boldsymbol{X}\)与\(\boldsymbol{Y}\)是直和。

5. \(X_1,\cdots, X_n\)的线性组合\(\boldsymbol{b}'\boldsymbol{X}\)是\(Y\)的最佳线性预测的等价条件是

   \[\mathbb{E}[X_j(Y-\boldsymbol{b}'\boldsymbol{X})]=0,\forall j=1,\cdots,n. \]

   这表明剩余部分与\(L^2(X)\)正交，因为\(X_1,\cdots,X_n\)是\(L^2(X)\)的基。

6. 如果\(\hat Y=L(Y|X_1,\cdots,X_n),\tilde Y=L(Y|X_1,\cdots,X_{n-1})\)，则

   \[L(\hat Y|X_1,\cdots,X_n)=\tilde{Y},\\ \mathbb{E}(Y-\hat{Y})^2 \le \mathbb{E}(Y-\tilde{Y})^2. \]

   这是投影的可分离性，可以理解成\(Y\)在\(L^2(X)\)子空间上的投影，在子空间的子空间上的投影\(\tilde Y\)与\(Y\)直接在子空间的子空间上投影效果一致。并且在父空间上的投影必定比子空间上更精准。

7. 如果\(\boldsymbol{X}=(X_1,\cdots,X_n)\)和\(\boldsymbol{Y}=(Y_1,\cdots,Y_n)\)满足存在矩阵\(A,B\)，使得\(A\boldsymbol{X}=\boldsymbol{Y}\)与\(B\boldsymbol{Y}=\boldsymbol{X}\)都成立，则

   \[L(Z|\boldsymbol{X}) = L(Z|\boldsymbol{Y}). \]

   这表明\(L^2(X)=L^2(Y)\)，也就是两个子空间中的基是等价向量组。

8. 勾股定理，即

   \[\mathbb{E}(Y^2)=\mathbb{E}[Y-L(Y|\boldsymbol{X})]^2+\mathbb{E}[L(Y|\boldsymbol{X})]^2. \]

3-4 时间序列的线性预测

时间序列的线性预测，实际上就是用时间序列的历史信息对未来进行线性预测的过程。由于此过程与回归分析类似（用一部分变量解释另一个变量），而解释变量出自时间序列自身，所以也可以称作自回归。

理论上，时间序列是没有源头的，即

\[\{X_t\}:\cdots,X_0,X_1,\cdots,X_n,\cdots \]

现实中，时间序列必定存在某个源头，记作\(X_1\)，即

\[\{X_t\}:X_1,\cdots,X_n,\cdots \]

因此对这两种情况的研究都有意义。但涉及到无穷时，我们要用线性闭包来表征Hilbert空间，即使用\(\bar L^2(X)\)来表征。

用无穷历史预测一步未来，即

\[\boldsymbol{X}_{t,m}=(X_t, X_{t-1},\cdots,X_{t-m+1}),\\ \hat{X}_{t+1, m}=L(X_{t+1}|\boldsymbol{X}_{t, m}),\\ \hat{X}_{t+1,\infty}=L(X_{t+1}|X_t, X_{t-1},\cdots) \]

其一步预测均方误差是

\[\sigma_1^2=\mathbb{E}(X_{t+1}-\hat{X}_{t+1, \infty})^2. \]

根据一步预测的均方误差，可以把时间序列分成两类：

1. \(\sigma^2=0\)，此时\(\{X_t\}\)是决定性平稳序列。
2. \(\sigma^2\ne 0\)，此时\(\{X_t\}\)是非决定性平稳序列。

用无穷历史预测\(k\)步未来也是一样的，有

\[\hat X_{t+k,\infty}=L(X_{t+k}|X_t,X_{t-1},\cdots),\\ \sigma^2_k=\mathbb{E}(X_{t+k}-\hat{X}_{t+k,\infty}),\\ \]

根据\(k\)步预测的均方误差，可以把时间序列分为两类：

1. \(\sigma_k^2\to\gamma_0\)，此时\(\{X_t\}\)是纯非决定性平稳序列。
2. \(\sigma_k^2\nrightarrow\gamma_0\)，此时\(\{X_t\}\)不是纯非决定性平稳序列。

这里\(\sigma_k^2\ge\sigma_{k-1}^2\)，且\(\sigma_k^2\le \gamma_0\)，根据单调有界定理，\(\{\sigma_k^2\}\)必有不大于\(\gamma_0\)的极限。

Wold表示定理：任一非决定性平稳序列（不一定零均值）\(\{X_t\}\)可以表示成

\[X_t=\sum_{j=0}^\infty a_j\varepsilon_{t-j}+V_{t},\\ \varepsilon_t = X_t - L(X_{t}|X_{t-1},X_{t-2},\cdots). \]

这里\(\varepsilon_t\)是零均值白噪声\({\rm WN}(0, \sigma^2)\)，\(a_0=1\)，且\(\{a_j\}\)平方可和；\(\{V_t\}\)是与白噪声序列正交的决定性平稳序列。令其白噪声部分为

\[U_t=\sum_{t=0}^{\infty} a_j\varepsilon_{t-j}, \]

则\(\{U_t\}\)是和\(\{V_t\}\)正交的纯非决定性平稳序列，并且\(\{\varepsilon_t\}\)与\(\{U_t\}\)有同样的历史张成空间。这指的是\(\forall t\)，\((\cdots,\varepsilon_{t-1},\varepsilon_t)\)张成的Hilbert空间与\((\cdots,U_{t-1},U_t)\)张成的Hilbert空间完全一致。

1. 称\(\{U_t\}\)是\(\{X_t\}\)的纯非决定性部分，\(\{V_t\}\)是\(\{X_t\}\)的决定性部分。
2. 称\(\{a_j\}\)是\(\{X_t\}\)的Wold系数。
3. 称一步预测误差\(\varepsilon_t\)是\(\{X_t\}\)的新息序列，\(\sigma^2=\mathbb{E}(\varepsilon_t)^2\)是一步预测的均方误差。
4. Wold表示定理说明，任何非决定性平稳序列可以表达为新息的单边滑动和。

而AR、MA、ARMA序列都是非决定性平稳序列，所以都可以表达为新息的单边滑动和。但并不意味着白噪声的单边滑动和中，白噪声就是新息，举例：

\[X_t=\varepsilon_t+2\varepsilon_{t-1}. \]

对于可逆的ARMA(p, q)序列来说，其Wold表示就是常规的表示，即

\[X_t=A^{-1}(\mathscr{B})B(\mathscr{B})\varepsilon_t, \]

这里的\(\varepsilon_t\)就是新息，Wold系数就是Wold表达中的Wold系数。

Kolmogorov公式：由谱密度可以求出均方误差（新息的二阶矩），为

\[\sigma^2=\mathbb{E}[X_{t}-L(X_t|{X}_{t-1},X_{t-2},\cdots)]^2=2\pi\exp\left\{\frac1{2\pi}\int_{-\pi}^\pi\ln f(\lambda){\rm d}\lambda \right\}. \]

用有限时间序列替代无限时间序列的合理性在于，设\(\boldsymbol{X}_{n,m}=(X_n,X_{n-1},\cdots,X_{n-m+1})'\)，\(Y\in L^2\)，则当\(m\to \infty\)时

\[L(Y|\boldsymbol{X}_{n,m})\stackrel{m.s.}{\to}L(Y|\boldsymbol{X}_{n,\infty}). \]

接下来转向有限时间序列的预测。设\(\{Y_t\}\)是方差有限的零均值时间序列（不一定平稳），接下来定义\(\boldsymbol{Y}_n=(Y_1,\cdots,Y_n)'\)，也用它表示由自己生成的线性闭包；\(W_1=Y_1,W_t=Y_t-L(Y_t|\boldsymbol{Y}_{t-1})\)，定义\(\boldsymbol{W}_n=(W_1,\cdots,W_n)'\)，也用它表示由它自己生成的线性闭包。可以论证\(\boldsymbol{Y}_n=\boldsymbol{W}_n\)：

\[W_n = Y_n - L(Y_n|\boldsymbol{Y}_{n-1}) \in \boldsymbol{Y_n}.\\ \]

采用归纳法，已知\(Y_1=W_1\)故\(\boldsymbol{Y}_1=\boldsymbol{W}_1\)，假设\(\boldsymbol{Y}_k=\boldsymbol{W}_k\)，则对\(k+1\)，有

\[Y_{n+1} = W_{n+1} + L(Y_{n+1}|\boldsymbol{W_n})\in \boldsymbol{W}_{n+1}, \]

结合\(W_{n+1}\in\boldsymbol{Y}_{n+1}\)就得到生成空间的等价性。这里的关键在于，对于有源头的时间序列有\(Y_1=W_1\)。

这个等价性表明，用\(\boldsymbol{Y_n}\)预测\(Y_{n+1}\)和用\(\boldsymbol{W}_n\)是等价的，而用\(\boldsymbol{W}_{n}\)进行预测的好处在于，\(\boldsymbol{W}_n\)中每一个随机变量都正交，从而是一个正交向量。

时间序列（不要求平稳）的递推预测公式：如果零均值时间序列\(\{Y_t\}\)的\(m\)阶协方差矩阵正定，则最佳线性预测为

\[\hat Y_{n+1} = \sum_{j=1}^n\theta_{n,j} W_{n+1-j},\quad \forall n \le m. \]

这里\(\theta_{n,j}\)是线性递推预测系数，定义预测的均方误差为\(\nu_n=\mathbb{E}(W_{n+1}^2)\)，则预测系数与均方误差满足如下的递推关系：

\[\nu_0=\mathbb{E}(Y_1^2), \\ \theta_{n,n}=\frac{\mathbb{E}(Y_{n+1}Y_1)}{\nu_0}, \\ \theta_{n,n-k}=\frac{\mathbb{E}(Y_{n+1}Y_{k+1})-\sum_{j=0}^{k-1}\theta_{k,k-j}\theta_{n,n-j}\nu_j}{\nu_k},1\le k\le n-1,\\ \nu_n = \mathbb{E}(Y_{n+1}^2)-\sum_{j=0}^{n-1}\theta_{n,n-j}^2\nu_j,1\le n\le m. \]

递推的顺序是一个三角形：

\[\begin{matrix} \nu_0 \\ \theta_{1,1} & \nu_1 \\ \theta_{2,2} & \theta_{2,1} & \nu_2 \\ \theta_{3,3} & \theta_{3,2} & \theta_{3,1} & \nu_3 \\ \vdots & \vdots & \vdots & \vdots & \ddots \end{matrix} \]

如果给定的时间序列\(\{X_t\}\)加上了平稳的条件，自协方差函数为\(\gamma_k\)，则其递推公式可以写成下面的形式：

\[\nu_0=\gamma_0, \\ \theta_{n,n}=\frac{\gamma_n}{\nu_0}, \\ \theta_{n,n-k}=\frac{\gamma_{n-k}-\sum_{j=0}^{k-1}\theta_{k,k-j}\theta_{n,n-j}\nu_j}{\nu_k},1\le k\le n-1,\\ \nu_n = \gamma_0-\sum_{j=0}^{n-1}\theta_{n,n-j}^2\nu_j,1\le n\le m. \]

三大模型线性递推预测的相关结论放在下一部分。

3-5 AR模型参数估计

ARMA(p, q)模型的参数估计指的是根据观测数据\(x_1,\cdots,x_N\)建立一个ARMA(p, q)模型，并且估计出模型系数。这包括定阶与求模型系数。AR(p)模型与MA(q)模型的参数估计也大概如此。建立模型的第一步是进行零均值化，这里直接假设\(x_1,\cdots,x_N\)是零均值化以后的观测数据。

现在假设AR(p)模型的阶数已知，自回归系数为\(\boldsymbol a=(a_1,\cdots,a_p)'\)，则AR(p)模型系数的矩估计是基于Yule-Walker方程的。首先根据观测样本计算\(\hat\gamma_0,\hat\gamma_1,\cdots,\hat{\gamma}_p\)，然后通过Yule-Walker方程计算：

\[\begin{bmatrix} \hat{\gamma_1} \\ \hat{\gamma}_2 \\ \vdots \\ \hat{\gamma}_p \end{bmatrix}=\begin{bmatrix} \hat{\gamma}_0 & \hat{\gamma}_1 & \cdots & \hat{\gamma}_{p-1} \\ \hat{\gamma}_1 & \hat{\gamma}_0 & \cdots & \hat{\gamma}_{p-2}\\ \vdots & \vdots & \ddots & \vdots \\ \hat{\gamma}_{p-1} & \hat{\gamma}_{p-2} & \cdots & \hat{\gamma}_0 \end{bmatrix}\begin{bmatrix} \hat{a}_1 \\ \hat{a}_2 \\ \vdots \\ \hat{a}_p \end{bmatrix}.\\ \hat{\sigma^2} = \hat{\gamma}_0-\sum_{j=1}^p\hat{a}_j\hat\gamma_j. \]

当\(n\)很大时，可以使用Levinson递推计算\(\hat a_1,\cdots,\hat a_p\)，这样得到的估计称为矩估计，也叫Yule-Walker估计。其优点是：

1. 计算简便；

2. 满足最小相位条件，即

   \[\hat{A}(z)=1-\sum_{j=1}^p\hat a_jz^j\ne 0,\quad \forall |z|\le 1. \]

3. 只要样本自协方差函数是强相合估计，那么\(\hat{\boldsymbol a}\)是强相合估计。

AR(p)模型系数的OLS估计，指的是用历史\(p\)个样本解释下一个样本，建立一个回归模型。接下来用\(d_1,\cdots,d_p\)表示AR(p)模型的最小二乘估计，则残差为

\[\hat{\epsilon}_j=x_j-(d_1x_{j-1}+d_2x_{j-2}+\cdots+d_px_{j-p}),\quad j>p. \]

OLS估计使得残差平方和最小，\(N\)个数据可以定义\(N-p\)个回归方程，所以

\[\underline{X}=\begin{bmatrix} x_p & x_{p-1} & \cdots & x_1 \\ x_{p+1} & x_p & \cdots & x_2 \\ \vdots & \vdots & \ddots & \vdots \\ x_{N-1} & x_{N-2} & \cdots & x_{N-p} \end{bmatrix},\underline{Y}=\begin{bmatrix} x_{p+1} \\ x_{p+2} \\ \vdots \\ x_{N} \end{bmatrix},\boldsymbol d=\begin{bmatrix} d_1 \\ d_2 \\ \vdots \\ d_p \end{bmatrix}, \\ \boldsymbol{d} = (\underline{X}'\underline{X})^{-1}(\underline{X}'\underline{Y}). \]

白噪声方差的估计就是残差平方和除以其自由度\(N-p\)（见下面的MLE估计）。OLS估计与Yule-Walker估计差别不大，它们的差是与\(1/N\)同阶的。

AR(p)模型的MLE估计给白噪声\({\rm WN}(0, \sigma^2)\)提出了更高的要求——服从正态分布。在这个假定下，给出了残差\(\varepsilon_{p+1},\cdots,\varepsilon_{N}\)的联合密度函数：

\[\left(\frac{1}{2\pi\sigma^2} \right)^{\frac{N-p}{2}}\exp\left\{-\frac{1}{2\sigma^2}\sum_{t=p+1}^N\left(x_{t}-\sum_{j=1}^p a_jx_{t-j} \right)^2 \right\}. \]

将中间复杂的这一块写成

\[S(\boldsymbol{a})=\sum_{t=p+1}^{N}\left(x_t-\sum_{j=1}^pa_jx_{t-j} \right)^2, \]

解此似然方程，得到的最大似然估计其实是使得\(S(\boldsymbol{a})\)最小的点，因而就是OLS估计，也与Yule-Walker估计相差不大，故\(\sigma^2\)的估计也自然是

\[\hat\sigma^2=\frac{S(\boldsymbol{a})}{N-p}. \]

前面假定了AR(p)模型已知，实际上由观测值不能推知AR(p)模型的阶数，所以需要一种定阶方式。定理表明，AR(p)模型的Yule-Walker估计在\(N\to \infty\)的情况下也是\(p\)后截尾的，因此可以根据Yule-Walker估计来定阶。

AIC定阶方法：如果根据问题的背景或数据的特性能够判定AR模型的阶数上界是\(P_0\)，则对于\(k=0,1,\cdots,P_0\)，分别假设AR模型的阶数是\(k\)，估计相应的AR(k)模型白噪声方差\(\hat\sigma^2_k\)。引入AIC函数为

\[{\rm AIC}(k)=\ln(\sigma^2_k)+\frac{2k}{N},\quad k = 0,1,\cdots,P_0, \]

使得AIC函数最小的点\(\hat p\)（如不唯一，应取最小值）称为AR(p)模型的AIC定阶。

1. AIC定阶并不是相合的，即使样本量增大，\(\hat p\)也不会依概率收敛到真实的阶数\(p\)。
2. AIC定阶通常会对阶数略有高估。但一般说来，高估AR(p)模型的阶数并不引起严重后果，而低估阶数则会带来较大的模型误差。而且实际数据中不存在阶数，高估阶数有利于多使用历史数据。
3. 如果样本量不是很大，人们乐于使用AIC定阶。

BIC定阶方法：在AIC定阶方法中，将AIC函数换成BIC函数为

\[{\rm BIC}(k)=\ln(\sigma^2_k)+\frac{k\ln N}{N},\quad k=0,1,\cdots,P_0. \]

使得BIC函数最小的点\(\hat{p}\)称为AR(p)模型的BIC定阶。

1. BIC定阶是强相合的。
2. 当\(N\)不是很大时，BIC定阶有时会低估模型阶数，造成模型的较大失真。
3. 在实际问题中，特别当样本量不是很大时，BIC定阶的效果不如AIC定阶。

AR(p)模型的拟合检验：使用估计的阶数与参数计算残差，对残差序列进行白噪声检验。

AR(p)序列谱密度的估计：使用\(\hat\sigma^2\)代替谱密度中的\(\sigma^2\)，得到谱密度的估计，称为AR谱估计或者极大熵谱估计。对于AIC定阶或BIC定阶\(\hat p\)，如果\(\{\varepsilon_t\}\)是独立同分布的\({\rm WN}(0, \sigma^2)\)，则\(\hat f(\lambda)\)一致收敛到\(f(\lambda)\)。

3-6 MA模型的参数估计

这里直接假设\(x_1,\cdots,x_N\)是零均值化以后的观测数据，与AR模型相同，首先假设MA模型的阶数\(q\)已知进行讨论。假设这一列数据满足的MA(q)模型是

\[X_t=\varepsilon_t+\sum_{j=1}^q b_j\varepsilon_{t-j}, \]

模型参数是\(\boldsymbol{b}=(b_1,\cdots,b_q)'\)与\(\sigma^2\)。
MA(q)模型系数的矩估计：由于MA(q)模型的自协方差函数\(\gamma_k\)是q后截尾的，存在如下的方程组：

\[\gamma_k=\sigma^2(b_0b_{k}+b_1b_{k+1}+\cdots+b_{q-k}b_q),\quad 0\le k\le q, \]

所以可以由这\(k\)个方程组求解\(k\)个系数，这是一个非线性方程组，使用线性迭代法或者Newton-Raphson迭代法可以得到，但迭代得到的解不一定使可逆条件成立，这时需要改变初值重新进行迭代。

MA(q)模型系数存在着矩阵表示，所以可以用矩阵表示反解模型参数，此时要令\(k>q\)时的\(\hat\gamma_k=0\)，然后构造矩阵

\[A=\begin{bmatrix} 0 & 1 & 0 & \cdots & 0 & 0\\ 0 &0 & 1 & \cdots & 0&0 \\ \vdots& \vdots & \vdots &\ddots & \vdots &\vdots\\ 0&0&0&\cdots &0&1 \\ 0&0&0&\cdots &0&0 \end{bmatrix}_{q\times q},C=\begin{bmatrix} 1 \\ 0 \\ \vdots \\ 0 \end{bmatrix}_{q\times 1}, \\ \hat\Omega_k=\begin{bmatrix} \hat\gamma_1 & \hat\gamma_2 & \cdots & \hat\gamma_k \\ \hat\gamma_2 & \hat\gamma_3 & \cdots & \hat\gamma_{k+1} \\ \vdots & \vdots & \ddots & \vdots \\ \hat\gamma_q & \hat\gamma_{q+1} & \cdots & \hat\gamma_{k+q-1} \end{bmatrix}_{q\times k},\boldsymbol {\hat\gamma_q}=\begin{bmatrix} \hat\gamma_1 \\ \hat\gamma_2 \\ \vdots \\ \hat\gamma_q \end{bmatrix}_{q\times 1},\\ \hat\Pi=\lim_{k\to \infty}\hat\Omega_k\hat\Gamma_k^{-1}\hat\Omega_k', \]

基于此计算

\[\hat\sigma^2=\hat\gamma_0-C'\hat\Pi C,\quad \hat{\boldsymbol b}_q=\frac1{\sigma^2}(\boldsymbol {\hat\gamma}_q-A\hat\Pi C). \]

基于此计算得到的矩估计量具有以下性质：如果\(\{\varepsilon_t\}\)是独立同分布的\({\rm WN}(0,\sigma^2)\)，则几乎必然地，\(N\to \infty\)时矩估计量满足可逆条件。

MA(q)模型的逆相关函数法是基于逆谱密度得到的，由于MA(q)序列谱密度取倒数后，可以成为AR(q)序列的谱密度，所以定义逆谱密度：设平稳序列\(\{X_t\}\)有恒正谱密度\(f(\lambda)\)，则称其逆谱密度为

\[f_y(\lambda)=\frac{1}{4\pi^2f(\lambda)}, \]

由它决定的自协方差函数称为\(\{X_t\}\)的逆相关函数即

\[\gamma_y(k)=\int_{-\pi}^\pi e^{{\rm i}k\lambda}f_y(\lambda){\rm d}\lambda. \]

MA(q)模型\(X_t=B(\mathscr{B})\varepsilon_t\)的逆谱密度是AR(p)模型

\[B(\mathscr{B})Y_t=\epsilon_t,\quad \mathbb{D}(\epsilon_t)=\frac{1}{\sigma^2} \]

的谱密度，也就是\(\{Y_t\}\)的自协方差函数是\(\{X_t\}\)的逆相关函数。

逆相关函数法的思想是，根据\(\{X_t\}\)的逆相关函数，用Yule-Walker估计确定MA(q)模型的系数。然而，MA(q)模型的逆相关函数是不好计算的，AR(p)模型的逆相关系数却容易计算：设

\[X_t=\sum_{j=1}^p a_jX_{t-j}+\varepsilon_t,\mathbb{D}(\varepsilon_t)=\sigma^2, \]

则\(\{X_t\}\)的逆相关函数为

\[\gamma_y(k)=\frac{1}{\sigma^2}\sum_{j=0}^{p-l}a_ja_{j+k},\quad 0\le k\le p,a_0=-1. \]

由于可逆MA(q)模型可以看成是一个AR(无穷)模型，且自协方差函数列以负指数阶收敛到0，所以可以将MA(q)模型近似为一个AR(p)模型，取较大的p即可。操作步骤如下：

1. 根据观测值，建立一个AR(p)序列，这里取阶数为AR模型的AIC定阶，或者取一个正数\(K\)，用\(K\ln N\)的正数部分作阶。
2. 对这个AR(p)模型，利用Yule-Walker方程解出AR模型系数。
3. 计算样本逆相关函数，将其视为MA(q)模型的逆相关函数。
4. 将这个逆相关函数视为另一个AR(q)模型的自协方差函数，再次运用Yule-Walker方程求解模型系数。
5. 将此时AR(q)模型的特征多项式，视为MA(q)模型的特征多项式；将此时AR(q)模型白噪声方差的倒数，视为MA(q)模型的白噪声方差。

MA(q)模型的新息估计方法思想是，根据递推预测系数来拟合模型参数，取\(m=O(N^{1/3})\)，计算样本自协方差函数\(\hat\gamma_0,\cdots,\hat\gamma_m\)，则新息估计为

\[(\hat{b}_1,\cdots,\hat{b}_q)=(\hat\theta_{m,1},\cdots,\hat\theta_{m,q}),\hat{\sigma}^2=\hat\nu_m. \]

这里样本递推系数与样本均方误差的线性递推公式为

\[\hat\nu_0=\hat\gamma_0, \\ \hat{\theta}_{n,n}=\frac{\hat{\gamma}_k}{\hat\nu_0}, \\ \hat{\theta}_{n,n-k}=\frac{\hat\gamma_{n-k}-\sum_{j=0}^{k-1}\hat\theta_{k,k-j}\hat\theta_{n,n-j}\hat\nu_j}{\nu_k^2}, \\ \hat{\nu}_n=\hat{\nu_0}-\sum_{j=0}^{n-1}\hat\nu_j\hat\theta^2_{n,n-j}. \]

MA(q)模型也有AIC和BIC两种定阶方法，类似于AR(p)模型的，即

\[{\rm AIC}(k)=\ln(\sigma^2_k)+\frac{2N}{k}, \\ {\rm BIC}(k)=\ln(\sigma^2_k)+\frac{k\ln N}{N}. \]

MA(q)模型的拟合检验也是用白噪声检验，令

\[\hat\varepsilon_{1-\hat q}=\hat\varepsilon_{2-\hat{q}}=\cdots=\hat\varepsilon_0=0,\\ \hat{\varepsilon_t}=x_t-\sum_{j=1}^q \hat{b}_j\hat\varepsilon_{t-j}. \]

取\(L=O(N^{1/3})\)，如果\(L\)之后的白噪声能够通过白噪声检验，就认为选择的模型合适。

3-7 ARMA模型的参数估计

ARMA模型是当AR模型和MA模型都不合适时才选择的模型，它的模型系数是\((\boldsymbol{a}',\boldsymbol{b}',\sigma^2)\)，模型形式是

\[X_t=\sum_{j=1}^p a_jX_{t-j}+\varepsilon_t+\sum_{j=1}^qb_j\varepsilon_{t-j}. \]

ARMA模型的矩估计，指的是使用延伸的Yule-Walker方程

\[\begin{bmatrix} \gamma_{q+1} \\ \gamma_{q+2} \\ \vdots \\ \gamma_{q+p} \end{bmatrix}=\begin{bmatrix} \gamma_q & \gamma_{q-1} & \cdots & \cdots_{q-p+1} \\ \gamma_{q+1} & \gamma_q & \cdots & \gamma_{q-p+2} \\ \vdots & \vdots & \ddots & \vdots \\ \gamma_{q+p-1} & \gamma_{q+p-2} & \cdots & \gamma_q \end{bmatrix}\begin{bmatrix} a_1 \\ a_2 \\ \vdots \\ a_p \end{bmatrix} \]

估计模型系数\(\boldsymbol{a}\)，然后扣掉AR部分得到一个MA(q)模型，继续估计MA(q)模型的参数。

ARMA模型的最小二乘估计基于自回归逼近思想，指的是先为数据建立一个AR模型。由于ARMA模型与AR模型不同，所以要对AR模型重新AIC估阶\(\hat p\)，计算得到残差序列\(\hat\varepsilon_t\)，然后用残差序列代替白噪声写出近似的ARMA模型，使用最小二乘估计同时得到\(\boldsymbol{a},\boldsymbol{b}\)，残差平方和为\(Q(\boldsymbol{a}, \boldsymbol{b})\)，最后得到\(\sigma^2\)的OLS估计为

\[\frac{1}{N-L}Q(\boldsymbol{a}, \boldsymbol{b}), \\ L=\max\{\hat p, p, q\}. \]

ARMA模型的极大似然估计首先要知道时间序列的联合分布，这里我们总假设预测误差是服从正态分布的。现在假设\(\{X_t\}\)是一个ARMA序列，假设\(\boldsymbol{X}_n=(X_1,\cdots,X_n)'\)且其自协方差函数正定，则由时间序列的递推预测，有

\[L(X_n|\boldsymbol{X}_{n-1})=\sum_{j=1}^{n-1}\theta_{n-1, n-j}Z_j,\\ X_n=L(X_n|\boldsymbol{X}_{n-1})+Z_n. \]

为了简化表达，令\(\boldsymbol{Z}_n=(Z_1,\cdots,Z_n)'\)，则有

\[\begin{bmatrix} X_1 \\ X_2 \\ \vdots \\ X_n \end{bmatrix}=\begin{bmatrix} 1 \\ \theta_{1,1} & 1 \\ \theta_{2,2} & \theta_{2,1} & 1 \\ \theta_{3,3} & \theta_{3,2} & \theta_{3,1} & 1 \\ \vdots & \vdots & \vdots & \vdots & \ddots \\ \theta_{n-1,n-1} & \theta_{n-1,n-2} & \theta_{n-1,n-3} & \cdots & \theta_{n-1,1} & 1 \end{bmatrix}\begin{bmatrix} Z_1 \\ Z_2 \\ \vdots \\ Z_n \end{bmatrix}, \\ \boldsymbol{X}_n=C\boldsymbol{Z}_n. \]

这样表达是因为\(Z_1,\cdots,Z_n,\cdots\)相互正交，方差定义为\(\mathbb{D}(\nu_k)=r_{k-1}\)，则

\[D\xlongequal{def}\mathbb{E}(\boldsymbol{Z}_n'\boldsymbol{Z}_n)={\rm diag}(\nu_0,\cdots,\nu_{n-1}). \]

计算\(\boldsymbol{X}_n\)的分布，它应该服从多元正态分布，且

\[\boldsymbol{Z}_n\sim N_n(0,D),\\ \boldsymbol{X}_n\sim N_n(0,CDC')\xlongequal{def}N_n(0,\Gamma_n).\\ \]

所以基于\(\Gamma_n\)可以得到似然函数为

\[L(\Gamma_n)=\frac{1}{(2\pi)^{n/2}\det(\Gamma_n)^{1/2}}\exp\left(-\frac12\boldsymbol{X}_n'\Gamma_n^{-1}\boldsymbol{X}_n \right). \]

并且

\[\det(\Gamma_n)=\nu_0\nu_1\cdots\nu_{n-1}, \\ \boldsymbol{X}_n'\Gamma^{-1}\boldsymbol{X}_n=\sum_{j=1}^{n} \frac{Z_j^2}{\nu_{j-1}},\\ L(\Gamma_n)=\frac{1}{(2\pi)^{n/2}(\nu_0\nu_1\cdots\nu_{n-1})^{1/2}}\exp\left(-\frac12\sum_{j=1}^n\frac{Z_j^2}{\nu_{j-1}} \right). \]

此时令\(\boldsymbol{\beta}=(a_1,\cdots,a_p,b_1,\cdots,b_q)'\)，只要找到模型参数\(\boldsymbol{\beta},\sigma^2\)使\(L(\Gamma_n)\)最小。事实上递推系数仅与\(\boldsymbol \beta\)有关，而均方误差则与\(\sigma^2\)和\(\boldsymbol \beta\)都有关，因此这里可以作方差提取处理，就得到

\[L(\boldsymbol {\beta},\sigma^2)=\frac{1}{(2\pi)^{n/2}(\sigma^{2n}\nu_0\nu_1\cdots\nu_{n-1})^{1/2}}\exp\left(-\frac{1}{2\sigma^2}\sum_{j=1}^n\frac{Z_j^2}{\nu_{j-1}} \right), \]

这里的\(\nu\)与上式中的\(\nu\)不一样，这也使得\(\nu\)也只与\(\boldsymbol \beta\)有关（参考ARMA模型最佳线性预测的证明）。令

\[\sum_{j=1}^{n} \frac{Z_j^2}{\nu_{j-1}}\xlongequal{def}S(\boldsymbol{\beta}), \]

让\(L(\Gamma_n)\)最大化，得到MLE估计\(\hat{\boldsymbol{\beta}}\)，进而得到MLE估计的\(\sigma^2\)为

\[\hat\sigma^2=\frac{1}{N}S(\hat{\boldsymbol{\beta}}). \]

这个\(\hat{\boldsymbol{\beta}}\)难以解析，所以通常使用搜索，而且往往用矩估计或者OLS估计得到的估计量作为搜索的初值。
将这个估计值\(\hat\sigma^2\)代回对数似然函数\(l(\boldsymbol{\beta},\sigma^2)\)就得到约化似然函数\(l(\boldsymbol{\beta})\)。对于较大的\(N\)，\(l(\boldsymbol{\beta})\)和\(S(\boldsymbol{\beta})\)的最小值点近似相等，所以也可以用\(S(\boldsymbol{\beta})\)的最小值点\(\tilde{\boldsymbol \beta}\)作为\(\boldsymbol{\beta}\)的估计，此时对应的\(\sigma^2\)估计为

\[\tilde{\sigma}^2=\frac{1}{N-p-q}S(\tilde{\boldsymbol{\beta}}). \]

ARMA模型的矩估计在自回归部分的参数表现较好，在滑动平均部分的参数估计比较差，极大似然估计能改进滑动平均部分的参数估计。

ARMA模型的检验：取

\[x_0=x_{-1}=\cdots=x_{-p+1} = 0; \\ \hat\varepsilon_0 = \hat{\varepsilon}_{-1} = \cdots=\hat{\varepsilon}_{-q+1}=0; \\ \varepsilon_t=x_t-\sum_{j=1}^p\hat a_jx_{t-j}-\sum_{j=1}^q\hat b_j\varepsilon_{t-j},\quad t=1,2,\cdots,N. \]

然后取\(m=O(N^{1/3})\)，对\(m\)以后的残差执行白噪声检验。
ARMA模型的AIC定阶与BIC定阶方法，与之前的AR模型、MA模型类似，即计算

\[{\rm AIC}(k,j)=\ln(\hat\sigma^2_{k,j})+\frac{2(k+j)}{N}, \\ {\rm BIC}(k,j)=\ln(\hat{\sigma}^2_{k,j})+\frac{(k+j)\ln N}{N}. \]

如果最小值不唯一，优先选择\(k+j\)最小的，然后选\(j\)最小的（滑动平均部分）。

自动建模：根据数据的性质有时能知道\(p\)和\(q\)的上界\(P_0\)、\(Q_0\)，可以在范围内把所有模型建出来，一共建立\(P_0Q_0\)个模型，将其中通过假设检验的\(p+q\)最小的模型作为实际模型。如果不能确定上界，则从\(p+q=1,2,\cdots\)开始一步步建模。