时间序列析步骤及程序详解（python）

时间序列析步骤及程序详解（python）

---

---

前言

透过实例《应用时间序列分析（第6版）》王燕 时间序列的两道例题4.4、4.7，以下均为python编程，其他语言的时间序列代码可私聊我。
注意：代码不把图叉掉是不会运行后续代码的，要是没有出现图，可以加一个plt.show()

---

城市未来的人口死亡率情况

先导入第三方库

clear;clc
from __future__ import print_function
import pandas as pd
import numpy as np
from scipy import  stats
import matplotlib.pyplot as plt
import statsmodels.api as sm
from statsmodels.graphics.api import qqplot
from statsmodels.stats.diagnostic import acorr_ljungbox

导入数据

data=pd.Series([3.665,4.247,4.674,3.669,4.752,4.785,5.929,4.468,5.102,4.831,6.899,5.337,5.086,
      5.603,4.153,4.945,5.726,4.965,1.82,3.723,5.663,4.739,4.845,4.535,4.774,5.962,
      6.614,5.255,5.355,6.144,5.59,4.388,3.447,4.615,6.032,5.74,4.391,3.128,3.436,
      4.964,6.332,7.665,5.277,4.904,4.83,])

1、绘制该序列的时序图

data=pd.Series(data)
data.index = pd.Index(sm.tsa.datetools.dates_from_range('1978','2022'))
data.plot(figsize=(12,8))

2、判断该序列的平稳性与纯随机性

(i)平稳性检验

# ADF单位根平稳型检验
ADF = sm.tsa.stattools.adfuller(data)
print(ADF[0:2]) #输出统计量和对应的p值
if ADF[4]['1%'] > ADF[0]:
      print('ADF单位根平稳型检验通过，显著的拒绝原假设，数据是平稳的')
else:
      print('ADF单位根平稳型检验不通过，接受原假设，数据是不平稳的')

显著性水平1%、%5、10%不同程度拒绝原假设的统计值和ADF Test result的比较，ADF Test result同时小于1%、5%、10%即说明非常好地拒绝该假设，本数据中，adf结果为-5.18696， 小于三个level的统计值。小于显著性水平0.01对应的统计量临界值-3.59250。对应的P值非常接近0，本数据中，P值为9e-6,接近0。
要是没有通过平稳性检验就需要进行差分，差分代码如下

#一阶差分
fig = plt.figure(figsize=(12,8))
ax1= fig.add_subplot(111)
diff1 = data.diff(1)
diff1.plot(ax=ax1)

#二阶差分
fig = plt.figure(figsize=(12,8))
ax2= fig.add_subplot(111)
diff2 = data.diff(2)
diff2.plot(ax=ax2)

(ii)纯随机性检验

# 白噪声检验
bzs = acorr_ljungbox(data.dropna(), lags = [i for i in range(1,12)],boxpierce=True)
print(bzs)

这里的第二行为p值，对应于每一个自回归系数对应的p值，如果全都大于0.05，那么这个序列就是白噪声序列，不具备任何价值。

3、考察该序列的自相关系数和偏自相关系数的性质

# 检查平稳时间序列的自相关图和偏自相关图
dta= data.diff(1)#我们已经知道要使用一阶差分的时间序列，之前判断差分的程序可以注释掉
fig = plt.figure(figsize=(12,8))
ax1=fig.add_subplot(211)
fig = sm.graphics.tsa.plot_acf(data,title="ACF",ax=ax1)
ax2 = fig.add_subplot(212)
fig = sm.graphics.tsa.plot_pacf(data,title="PACF",ax=ax2,method='ywm')
plt.show()

不难看出，根据该序列自相关图，可以视为：自相关图一阶截尾，或偏自相关 2 阶截尾。
也可以视为，自相关图一阶截尾，偏自相关是拖尾的。

4、尝试用多个模型拟合该序列的发展，并考察该序列的拟合模型优化问题

根据自相关和偏自相关，选择对应的模型，这里选择MA(1),AR(2),ARMA(1,1)

#参数检验
model1 = sm.tsa.arima.ARIMA(data, order=(0,0,1))
result1 = model1.fit()
model2 = sm.tsa.arima.ARIMA(data, order=(2,0,0))
result2 = model2.fit()
model3 = sm.tsa.arima.ARIMA(data, order=(1,0,1))
result3 = model3.fit()

print(result1.summary())
print(result2.summary())
print(result3.summary())
print('arma(1,1)参数检验存在问题，参数显著为零，淘汰arma(1,1)')
print('MA(1),AR(2)参数显著非零，残差检验均为白噪声序列，MA(1)模型在AIC和BIC更优，所以选择MA(1)为最优模型')
#AIC BIC HQIC
print('模型      ','AIC               ','BIC              ','HQIC          ')
print('MA(1)    ',arma_mod20.aic,arma_mod20.bic,arma_mod20.hqic)
print('AR(2)    ',arma_mod30.aic,arma_mod30.bic,arma_mod30.hqic)
print('ARMA(1,1)',arma_mod40.aic,arma_mod40.bic,arma_mod40.hqic)

输出结果（各个模型的aic和bic对比）：

各个模型的参数检验（这里只放不合格的arma(1,1)）

ARMA(1,1)参数检验存在问题，参数显著为零，淘汰ARMA(1,1)模型
MA(1),AR(2)参数显著非零，残差检验均为白噪声序列，MA(1)模型在AIC和BIC更优，所以选择MA(1)为最优模型（aic和bic越大越好）

5、利用最优拟合模型预测该城市未来5年的人口死亡率情况。

pred = result1.predict('2022','2027',dynamic=True, typ='levels')
print (pred)

以上代码，如有雷同，纯属搬砖