模型定阶：计算自相关函数ACF 和 偏自相关函数PACF

若你检验方法选择方法二，则你已经得到了ACF和PACF，但检验方法最好还是选择单位根检验

流程：

1.根据ACF图和PACF图判断，序列平稳化后，使用AR、MA、ARMA中哪种模型

2.判断p，q阶数是多少

（前两步都是看ACF和PACF图来判断）

3.可能ARMA有多个p，q决定的模型（因为可能ACF和PACF图可以看出多个p和q），这时要通过信息标准AIC和BIC来协助选择模型

选择AIC、BIC值最小的，因为AIC、BIC越小，模型越好

但要注意的是，这些准则不能说明某一个模型的精确度，也即是说，对于三个模型Ａ，Ｂ，Ｃ，我们能够判断出Ｃ模型是最好的，但不能保证Ｃ模型能够很好地刻画数据，因为有可能三个模型都是糟糕的。

赤池信息准则 (akaike information criterion) :AIC鼓励数据拟合的优良性但是尽量避免出现过度拟合(Overfitting)的情况。所以优先考虑的模型应是AIC值最小的那一个

贝叶斯信息准则 (bayesian information criterion) :

其中 L 是该模型下的最大似然, n 是数据数量， k 是模型的变量个数。

检查平稳时间序列的自相关图和偏自相关图。

通过sm.graphics.tsa.plot_acf和sm.graphics.tsa.plot_pacf得到图形其中 lags 表示滞后的阶数，以上分别得到 ACF 图和 PACF 图

例子1

1) 观察ACF和PACF图，判断是哪种模型

dta= dta.diff(1)#我们已经知道要使用一阶差分的时间序列，之前判断差分的程序可以注释掉
fig = plt.figure(figsize=(12,8))
ax1=fig.add_subplot(211)
fig = sm.graphics.tsa.plot_acf(dta,lags=40,ax=ax1)
ax2 = fig.add_subplot(212)
fig = sm.graphics.tsa.plot_pacf(dta,lags=40,ax=ax2)

通过两图观察得到：

自相关图显示滞后有三个阶超出了置信边界；

偏相关图显示在滞后1至7阶（lags 1,2,…，7）时的偏自相关系数超出了置信边界，从lag 7之后偏自相关系数值缩小至0 则有以下模型可以供选择：

ARMA(0,1)模型：即自相关图在滞后1阶之后缩小为0，且偏自相关缩小至0，则是一个阶数q=1的移动平均模型；

ARMA(7,0)模型：即偏自相关图在滞后7阶之后缩小为0，且自相关缩小至0，则是一个阶层p=3的自回归模型；

ARMA(7,1)模型：即使得自相关和偏自相关都缩小至零。则是一个混合模型。

2) AIC BIC选择其中最佳模型

dta是差分后的序列因为使用的是ARMA模型

arma_mod20 = sm.tsa.ARMA(dta,(7,0)).fit() # dta是差分后的序列
print(arma_mod20.aic,arma_mod20.bic,arma_mod20.hqic)
arma_mod30 = sm.tsa.ARMA(dta,(0,1)).fit()
print(arma_mod30.aic,arma_mod30.bic,arma_mod30.hqic)
arma_mod40 = sm.tsa.ARMA(dta,(7,1)).fit()
print(arma_mod40.aic,arma_mod40.bic,arma_mod40.hqic)
arma_mod50 = sm.tsa.ARMA(dta,(8,0)).fit()
print(arma_mod50.aic,arma_mod50.bic,arma_mod50.hqic)

可以看到ARMA(7,0)的AIC，BIC，HQIC均最小，因此是最佳模型。

例子2:

1) 观察ACF和PACF图，判断是哪种模型

从一阶差分序列的自相关图和偏自相关图可以发现:

自相关图拖尾或一阶截尾

偏自相关图一阶截尾,

所以我们可以建立ARIMA(1,1,0)、ARIMA(1,1,1)、ARIMA(0,1,1)模型。

2) AIC BIC选择其中最佳模型

data[“xt”]是没有差分的序列因为使用的是ARIMA模型

arma_mod20 = sm.tsa.ARIMA(data["xt"],(1,1,0)).fit() # data["xt"]是没有差分的数据
arma_mod30 = sm.tsa.ARIMA(data["xt"],(0,1,1)).fit()
arma_mod40 = sm.tsa.ARIMA(data["xt"],(1,1,1)).fit()
values = [[arma_mod20.aic,arma_mod20.bic,arma_mod20.hqic],[arma_mod30.aic,arma_mod30.bic,arma_mod30.hqic],[arma_mod40.aic,arma_mod40.bic,arma_mod40.hqic]]
df = pd.DataFrame(values,index=["AR(1,1,0)","MA(0,1,1)","ARMA(1,1,1)"],columns=["AIC","BIC","hqic"

选择模型MA(0, 1, 1),即ARIMA(0, 1, 1)

参数估计

from statsmodels.tsa.arima_model import ARIMA
model = ARIMA(data["xt"], order=(0,1,1))
result = model.fit()
print(result.summary())

                             ARIMA Model Results                              
==============================================================================
Dep. Variable:                   D.xt   No. Observations:                   36
Model:                 ARIMA(0, 1, 1)   Log Likelihood                -122.987
Method:                       css-mle   S.D. of innovations              7.309
Date:                Tue, 22 Dec 2020   AIC                            251.973
Time:                        09:11:55   BIC                            256.724
Sample:                    01-01-1953   HQIC                           253.631
                         - 01-01-1988                                         
==============================================================================
                 coef    std err          z      P>|z|      [0.025      0.975]
------------------------------------------------------------------------------
const          4.9956      2.014      2.481      0.013       1.048       8.943
ma.L1.D.xt     0.6710      0.165      4.071      0.000       0.348       0.994
                                    Roots                                    
=============================================================================
                  Real          Imaginary           Modulus         Frequency
-----------------------------------------------------------------------------
MA.1           -1.4902           +0.0000j            1.4902            0.5000
-----------------------------------------------------------------------------

模型检验

1) 参数的显著性检验

P<α（α通常为0.05），拒绝原假设，，认为该参数显著非零MA(2)模型拟合该序列，残差序列已实现白噪声

2）模型的显著性检验

resid = result.resid#残差
fig = plt.figure(figsize=(12,8))
ax = fig.add_subplot(111)
fig = qqplot(resid, line='q', ax=ax, fit=True)

qq图显示，我们看到红色的KDE线与N(0,1)平行，这是残留物正太分布的良好指标，说明残差序列是白噪声序列，模型的信息的提取充分，当让大家也可以使用前面介绍的检验白噪声的方法LB统计量来检验

ARIMA(0,1,1)模型拟合该序列，残差序列已实现白噪声，且参数均显著非零。说明ARIMA(0,1,1)模型是该序列的有效拟合模型

模型预测

pred = result.predict('1988', '1990',dynamic=True, typ='levels')
print (pred)

1988-01-01   278.35527
1989-01-01   283.35088
1990-01-01   288.34649
Freq: AS-JAN, dtype: float64

plt.figure(figsize=(12, 8))
plt.xticks(rotation=45)
plt.plot(pred)
plt.plot(data.xt)
plt.show()

预测结果分析