时间序列丨九种经典的基于「自回归+移动平均」序列预测模型及其 Python 实践

背景

这篇文章中主要学习 & 总结下经典的基于「自回归」和「滑动平均」的时间序列预测模型及其变体形式，毕竟经典的才是最实用的！🤔

主要包含但不局限于以下九种：

• 单变量时间序列预测

  1. 「AR」Autoregression，自回归模型；
  2. 「MA」Moving Average，滑动平均模型；
  3. 「ARMA」Autoregressive Moving Average，自回归滑动平均模型；
  4. 「ARIMA」Autoregressive Integrated Moving Average，差分自回归移动平均模型；
  5. 「SARIMA」Seasonal Autoregressive Integrated Moving-Average，季节性差分自回归移动平均模型；
  6. 「SARIMAX」Seasonal Autoregressive Integrated Moving-Average with Exogenous Regressors，外生变量的季节性差分自回归移动平均模型；

• 多变量时间序列预测

  1. 「VAR」Vector Autoregression，向量自回归模型；
  2. 「VARMA」Vector Autoregression Moving-Average，向量自回归滑动平均模型；
  3. 「VARMAX」Vector Autoregression Moving-Average with Exogenous Regressors，外生变量的向量自回归滑动平均模型；

可能包含其它 ARMA 的变体形式，后续遇到再进行学习。

注：本文使用 Python statsmodels 库测试不同模型性能，时间序列相关概念参考另一篇博文：时间序列丨基础概念理论 & 异常检测算法概述

AR

自回归模型 AR 描述当前值与历史值间的关系，利用序列历史时间数据和线性函数对下一步进行预测。

一般 p 阶自回归模型 AR（P） 定义如下：

$$X_{t}=\alpha_{1} X_{t-1}+\alpha_{2} X_{t-2}+\ldots+\alpha_{p} X_{t-p}+u_{t}$$

其中$u_t$ 是随机扰动项，如果其为白噪声^[1] ($u_t=\varepsilon_t$) 那么称为纯 AR（p） 过程，表达形式如下：

$$X_{t}=\alpha_{1} X_{t-1}+\alpha_{2} X_{t-2}+\ldots+\alpha_{p} X_{t-p}+\varepsilon_{t}$$

局限性：(1) 适用于没有趋势和季节性的时间序列。(2) 适用于预测与历史数据相关的时间序列。(3) 时间序列数据必须是平稳的。

Python 代码示例

# AR example
from statsmodels.tsa.ar_model import AutoReg
from random import random
# contrived dataset
data = [x + random() for x in range(1, 100)]
# fit model
model = AutoReg(data, lags=1)
model_fit = model.fit()
# make prediction
yhat = model_fit.predict(len(data), len(data))
print(yhat)

参考：statsmodels.tsa.ar_model.AutoReg

MA

在 AR 模型中，如果$u_t$ 不是白噪声，那么通常认为其可以通过$q$ 阶移动平均值估计，$q$ 阶 MA（q） 模型表达形式如下：

$$u_{t}=\varepsilon_{t}+\beta_{1} \varepsilon_{t-1}+\ldots+\beta_{q} \varepsilon_{t-q}$$

其中$\varepsilon_t$ 为白噪声序列。特殊情况下如果$X_t=u_t$ ，表示时间序列与历史值没有任何关系，只依赖于历史白噪声，即$X_{t}=\varepsilon_{t}+\beta_{1} \varepsilon_{t-1}+\ldots+\beta_{q} \varepsilon_{t-q}$

注意两点：

• AR 模型中历史白噪声是通过影响历史值来间接影响当前预测值的。
• MA 模型与计算时间序列的滑动平均是不同的。

Python 代码示例

# MA example
from statsmodels.tsa.arima.model import ARIMA
from random import random
# contrived dataset
data = [x + random() for x in range(1, 100)]
# fit model
model = ARIMA(data, order=(0, 0, 1)) # order=(p,d,q)
model_fit = model.fit()
# make prediction
yhat = model_fit.predict(len(data), len(data))
print(yhat)

参考：statsmodels.tsa.arima.model.ARIMA

ARMA

将 AR（p） 与 MA（q）结合可得一般的自回归移动平均模型 ARMA（p, q）^[2] ，形式化定义如下：

$$X_{t}=\alpha_{1} X_{t-1}+\alpha_{2} X_{t-2}+\ldots+\alpha_{p} X_{t-p}+\varepsilon_{t}+\beta_{1} \varepsilon_{t-1}+\ldots+\beta_{q} \varepsilon_{t-q}$$

上式表示 ARMA 方法将序列中的下一步预测结果为历史的观测值和残差的线性函数。

局限性：适用于无周期性和趋势性的平稳时间序列。

Python 代码示例：

# ARMA example
from statsmodels.tsa.arima.model import ARIMA
from random import random
# contrived dataset
data = [random() for x in range(1, 100)]
# fit model
model = ARIMA(data, order=(2, 0, 1))
model_fit = model.fit()
# make prediction
yhat = model_fit.predict(len(data), len(data))
print(yhat)

参考：statsmodels.tsa.arima.model.ARIMA

ARIMA

将 AR（p）、MA（q）和差分法 I（d）结合可得差分自回归移动平均模型 ARIMA（p, d, q），其中 d 是需要对数据进行差分的阶数。

差分的主要目的是将非平稳序列变平稳的方法，即计算相邻观测值间的差值$y_{t}^{\prime}=y_{t}-y_{t-1}$

在现实应用中，通常没有必要进行二阶以上的差分，ARIMA 模型构建的一般步骤为：

1. 对序列进行 ADF 平稳性检验；对于非平稳时间序列要先进行 d 阶差分，转化为平稳时间序列；
2. 对平稳时间序列分别求得其自相关系数（ACF）和偏自相关系数（PACF），通过对自相关图和偏自相关图的分析，得到最佳的阶数p、q；
3. 基于以上得到的参数 d、q、p 构建 ARIMA 模型，最后对模型进行效果检验。

Python 代码示例如下：

# ARIMA example
from statsmodels.tsa.arima.model import ARIMA
from random import random
# contrived dataset
data = [x + random() for x in range(1, 100)]
# fit model
model = ARIMA(data, order=(1, 1, 1))
model_fit = model.fit()
# make prediction
yhat = model_fit.predict(len(data), len(data), typ='levels')
print(yhat)

参考：statsmodels.tsa.arima.model.ARIMA

SARIMA

SARIMA 方法将序列中的下一步预测值为先前时间步长的差异观测值、误差、差异季节性观测值和季节性误差的线性函数。

表达形式如下：

$$\text { SARIMA } \quad(p, d, q) \text{季节项} + (P, D, Q, s) \text{非季节项}$$

注意模型的季节性部分中含有和非季节性非常类似的项，但包含了季节性时段的回溯。

Python 代码示例如下

# SARIMA example
from statsmodels.tsa.statespace.sarimax import SARIMAX
from random import random
# contrived dataset
data = [x + random() for x in range(1, 100)]
# fit model
model = SARIMAX(data, order=(1, 1, 1), seasonal_order=(0, 0, 0, 0))
model_fit = model.fit(disp=False)
# make prediction
yhat = model_fit.predict(len(data), len(data))
print(yhat)

参考：statsmodels.tsa.statespace.sarimax.SARIMAX

SARIMAX

SARIMAX 带有外生变量的 SARIMA 的扩展模型。^[3]

外生变量定义^[3:1]：外生变量也称为协变量，可以被认为是并行输入序列，其观测值与原始序列的时间步长相同。 初级序列可被称为内源数据以将其与外源序列进行对比。 外生变量的观察结果在每个时间步直接包含在模型中，并且不以与主要内生序列相同的方式建模（例如作为 AR、MA 等过程）。

Python 代码示例如下

# SARIMAX example
from statsmodels.tsa.statespace.sarimax import SARIMAX
from random import random
# contrived dataset
data1 = [x + random() for x in range(1, 100)]
data2 = [x + random() for x in range(101, 200)]
# fit model
model = SARIMAX(data1, exog=data2, order=(1, 1, 1), seasonal_order=(0, 0, 0, 0))
model_fit = model.fit(disp=False)
# make prediction
exog2 = [200 + random()]
yhat = model_fit.predict(len(data1), len(data1), exog=[exog2])
print(yhat)

参考：statsmodels.tsa.statespace.sarimax.SARIMAX

VAR

VAR 方法使用 AR 模型对多个时间序列的进行建模，是 AR 在多变量时间序列上的推广。^[4]

表达形式如下：

$$y_{t}=c+A_{1} y_{t-1}+A_{2} y_{t-2}+\cdots+A_{p} y_{t-p}+e_{t}$$

其中$c, e_t \in R^{n\times 1}$，$A_i\in R^{n\times n}$。误差项$e_t$ 满足均值为 0 且不存在自相关。

Python 代码示例如下

# VAR example
from statsmodels.tsa.vector_ar.var_model import VAR
from random import random
# contrived dataset with dependency
data = list()
for i in range(100):
    v1 = i + random()
    v2 = v1 + random()
    row = [v1, v2]
    data.append(row)
# fit model
model = VAR(data)
model_fit = model.fit()
# make prediction
yhat = model_fit.forecast(model_fit.y, steps=1)
print(yhat)

参考：statsmodels.tsa.vector_ar.var_model.VAR

VARMA

VARMA 方法使用 ARMA 模型对多个时间序列的进行建模。

Python 代码示例如下

# VARMA example
from statsmodels.tsa.statespace.varmax import VARMAX
from random import random
# contrived dataset with dependency
data = list()
for i in range(100):
    v1 = random()
    v2 = v1 + random()
    row = [v1, v2]
    data.append(row)
# fit model
model = VARMAX(data, order=(1, 1))
model_fit = model.fit(disp=False)
# make prediction
yhat = model_fit.forecast()
print(yhat)

参考：statsmodels.tsa.statespace.varmax.VARMAX

VARMAX

带有外生回归量的 VARMAX 是 VARMA 模型的扩展，可以处理没有趋势和季节性的多变量时间序列。

Python 代码示例如下

# VARMAX example
from statsmodels.tsa.statespace.varmax import VARMAX
from random import random
# contrived dataset with dependency
data = list()
for i in range(100):
    v1 = random()
    v2 = v1 + random()
    row = [v1, v2]
    data.append(row)
data_exog = [x + random() for x in range(100)]
# fit model
model = VARMAX(data, exog=data_exog, order=(1, 1))
model_fit = model.fit(disp=False)
# make prediction
data_exog2 = [[100]]
yhat = model_fit.forecast(exog=data_exog2)
print(yhat)

参考：statsmodels.tsa.statespace.varmax.VARMAX

---

1. 时间序列丨基础概念理论 & 异常检测算法概述 ↩︎

2. https://zhuanlan.zhihu.com/p/60648709 ↩︎

3. Autoregressive integrated moving average on Wikipedia ↩︎ ↩︎

4. Vector autoregression on Wikipedia ↩︎