1 基本定义

数据平稳性的图判断:

平稳时间序列的均值和方差都为常数,因此平稳时间序列的时序图应该围绕一条水平线上下波动,而且波动范围有界。

(a)非平稳:有明显的周期性,趋势性

平稳时间序列的序列值之间有短期相关性,则其表现特征是:自相关函数会很快地衰减到零附近

(b)非平稳:自相关函数衰减到零附近的速度比较慢

(c)非平稳:自相关图典型特征,三角对称关系(图1.13)

(d)非平稳 :自相关系数衰减到零的速度非常缓慢,而且呈现明显的周期规律(图1.14)

python 简单预测模型 python建立预测模型_python 简单预测模型

平稳AR(1)模型例图:

        (1)自相关指数按负指数衰减到零附近;

        (2)自相关函数呈现正负相见地衰减;

        (3)自相关函数具有‘伪周期’的衰减特征;

        (4)自相函数衰减速度非常快(短期相关性)。

python 简单预测模型 python建立预测模型_python_02

随机游走序列:

        (b)自回归系数

python 简单预测模型 python建立预测模型_python 简单预测模型_03

1=1,所以噪声

python 简单预测模型 python建立预测模型_python_04

对xt的影响不随时间的推移而减弱;随机游走在一段时间具有上升或下降趋势,也就是运动方向具有一定的持续性,与(a)相比随机游走序列不围绕某个值上下波动。

        (c)方差逐渐增大。

        (d)序列值有明显的减小趋势。

        (z)故 b,c,d都不是平稳序列。

    

python 简单预测模型 python建立预测模型_Data_05

2 序列预处理

2.1判断序列平稳性

(1)画时间序列图和自相关函数图

# 导入相关模块
import pandas as pd
import matplotlib.pyplot as plt
from statsmodels.stats.diagnostic import acorr_ljungbox
%matplotlib inline
plt.rcParams['font.family'] = 'simsun'
plt.rcParams['axes.unicode_minus'] = False
# 绘制时间序列图
Data = pd.read_csv('datas/Guwu.csv',index_col=0).squeeze('columns') 
fig = plt.figure(figsize=(12,4),dpi=150)
ax = fig.add_subplot(111)
ax.plot(Data,marker='o',linestyle='-',color='red')
ax.set_xlabel(xlabel='year',fontsize=17)
ax.set_ylabel(ylabel='output',fontsize=17)
plt.xticks(fontsize=15);plt.yticks(fontsize=15)
plt.title('time_1')
plt.tight_layout()
# 绘制自相关图
def ACF(ts, lag):
    ts_acf = acf(ts, nlags=lag, fft=False)
    fig = plt.figure(figsize=(12,4),dpi=150)
    ax = fig.add_subplot(111)
    ax.vlines(x=list(range(lag+1)),ymin=np.zeros(lag+1),ymax=ts_acf,color='b',linewidth=2.0)
    
    ax.axhline(y=0, linestyle=':', color='b')
    ax.axhline(y=-1.96/np.sqrt(len(ts)),linestyle='--',color='r')
    ax.axhline(y=1.96/np.sqrt(len(ts)),linestyle='--',color='r')
    ax.set_xlabel(xlabel='lag',fontsize=17)
    ax.set_xlabel(xlabel='lag',fontsize=17)
    plt.title('ACF_2',fontsize=20)
    plt.xticks(fontsize=15);plt.yticks(fontsize=15)
    fig.tight_layout()

ACF(guwu_ts,lag=16)

(2) 白噪声检验

# 6,12阶白噪声检验
acorr_ljungbox(Data, return_df=True, lags=[6,12])

判断序列为平稳非白噪声序列,再进行选择模型

2.2选择合适模型

已经画完自相关图,再画偏自相关图(PACF)

# 绘制偏自相关图
def PACF(ts,lag=20,fname='',xlabel=''):
    lag_acf = pacf(ts, nlags=lag)  # 与自相关函数不同处
    plt.vlines(x=[range(lag+1)], ymax=lag_acf, ymin=np.zeros(lag+1), linewidth=2.0, color='b')
    
    plt.axhline(y=0,linestyle=':', color='b')
    plt.axhline(y=-1.96/np.sqrt(len(ts)),linestyle='--', color='r')
    plt.axhline(y = 1.96/np.sqrt(len(ts)), linestyle='--', color='r')
    
    plt.xlabel(xlabel=xlabel, fontsize=17)
    plt.ylabel(ylabel='pacf', fontsize=17)
    plt.xticks(fontsize=15)
    plt.yticks(fontsize=15)
    plt.title('PACF_3',fontsize=20)
    plt.tight_layout()
    if not fname:
        plt.savefig(fname=fname)


PACF(guwu_ts,lag=8)

根据自相关函数图,偏自相关函数图的拖尾性和截尾性选择合适模型

2.3估计模型参数

# 估计未知参数
guwu_est = ARIMA(guwu_ts,order=(1,0,0)).fit()
print(guwu_est.summary().tables[1])

根据估计值写出模型口径

2.4模型检验

# 残差检验
guwu_resid = guwu_est.resid
acorr_ljungbox(guwu_resid,lags=[6,12],return_df=True, boxpierce=True)

残差检验的p值都大于显著性水平,则拟合模型的残差序列为白噪声序列,然后再进行系数检验

2.5参数检验

from scipy.stats import t

# 过度拟合检验
def pt(x,df):
    if x <= 0:
        pv = t.cdf(x,int(df))
    else:
        pv = t.sf(x,int(df))
    return pv


# ar.L1系数检验
t1 = 0.3681/0.109;x = pt(t1,15)
# 常数const检验
t0 = 0.8491 /0.055;z = pt(t0,15)
print(x,z)

2.6模型优化

模型优化,从上面自相关图,偏自相关图看出也可选择ARMA(1,2)当作模型,下面对该模型进行模型检验。

guwu_est2 = ARIMA(guwu_ts,order=(1,0,2)).fit()
guwu_resid2 = guwu_est2.resid
acorr_ljungbox(guwu_resid2,lags=[6,12],return_df=True,boxpierce=True)

模型残差检验通过,进一步提取AIC,BIC,AICc,QHIC的值与原来的进行比较,选择值小的

print(guwu_est.aic<guwu_est2.aic, guwu_est.bic<guwu_est2.bic, guwu_est.aicc<guwu_est2.aicc,guwu_est.hqic<guwu_est2.hqic)

2.7预测

from statsmodels.tsa.arima.model import ARIMA

#进行预测最终模型
guwu_pred = guwu_est.get_prediction(start=75,end=79)
confint1 = guwu_pred.conf_int(alpha=0.20)  # 80%置信区间
confint2 = guwu_pred.conf_int(alpha=0.05)  # 95%置信区间

confint1.columns = ['lower_80', 'upper_80']  # 重命名列
confint2.columns = ['lower_95', 'upper_95']  # 重命名列

confint = pd.concat([confint1, confint2], axis=1)


# 绘制预测图
fig = plt.figure(figsize=(12,4), dpi=150)
ax = fig.add_subplot(111)
ax.plot(guwu_ts[1:],color='y')
guwu_pred.predicted_mean.plot(ax=ax,color='b')
ax.fill_between(confint1.index,confint1.iloc[:,0],confint1.iloc[:,1],color='k',alpha=.1)
ax.fill_between(confint2.index,confint2.iloc[:,0],confint2.iloc[:,1],color='k',alpha=.1)
ax.set_xlabel(xlabel="time", fontsize=17)
plt.legend(loc=2,labels=['real','predicted'],fontsize=13)
plt.xticks(fontsize=15); plt.yticks(fontsize=15)
fig.tight_layout()