freeai logo
<<返回python首页 python

《Python 应用案例》

Python ARIMA时间序列模型预测销量

本节假设目前业务没有重大变化以及外部重大影响，比较稳定，业务部门找到数据部门，希望做未来10天的数据预测，以合理安排资源做好促销活动、广告投放费用预算。因为目前业务稳定，可以不用考虑太多因素对销量的影响，认为只和时间有关系，可以考虑时间序列ARIMA模型预测。

Python ARIMA时间序列模型预测销量

本文转自微信公众号：Python数据科学修炼之路 作者:老杨啊小强

本转载已获作者授权，禁止二次转载

目 录

  • 1、时间序列 & ARIMA模型简介

①时间序列简介

②AR自回归模型

③MA移动平均模型

④ARIMA差分自回归移动平均模型

  • 2、导入库，加载数据

  • 3、数据处理

  • 4、绘制原始序列趋势

  • 5、做d阶差分，平稳性检验

  • 6、网格搜索确定最优的p和q

①赤池信息准则

②贝叶斯准则

  • 7、模型检验

①杜宾-瓦特森检验（W-D检验）

②残差自相关系数ACF和偏相关系数PACF检验

③混成检验/Ljung-Box test(LB检验)

  • 8、模型结果预测

时间序列

时间序列是用来研究数据随时间变化趋势而变化的一类算法，常用于经济预测、股市预测、天气预测等偏宏观或者没有可控自变量的场景下。

时间序列的出发点是：认为事务的发展都具有连续性，都按照本身固有的规律进行，只要规律不变，则认为未来的发展趋势仍然会延续下去。相信“历史不会简单重复，但会重演”。

ARIMA

01、AR模型

AR是自回归模型，自回归即自己的当前值和自身的历史数据做回归。当前值是自己前一天、前p天的值作用的结果。

p阶自回归过程定义公式： 参数说明：

：阶数，当前值用前p期数据预测

：当前值

:是常数项

扰动项/误差项（服从独立同分布），如果是白噪声(期望为0，方差为1)，则AR为纯阶自回归序列。

【注意】：

①自回归模型必须满足平稳性要求。

②自相关要求有相关性，如果相关系数小于0.5，认为相关性不大，不宜使用。

③自回归模型首先要确定一个阶数p，表示用前几期的历史预测当前值。

④平稳性要求序列的均值和方差不发生变化

02、MA移动平均模型

在AR模型中，不是白噪声，通常认为是一个q阶移动平均，即：

【注意】： ①移动平均模型关注的是自回归模型AR误差的累计影响。 ②移动平均能有效消除预测中的随机波动。

03、ARIMA差分自回归移动平均模型

ARIMA(p,d,q)模型全称差分自回归移动平均模型,是AR(p)和MA(q)的结合，得到一个一般的自回归移动平均模型：

参数说明：

是时间序列转化平稳序列是的差分阶数

是自回归，是自回归阶数

是移动回归，是移动回归阶数。

注意：

①一个随机时间序列可以通过一个自回归移动平均模型来表示，即序列可以由自身的过去或者滞后项以及随机扰动项来解释。

②如果该序列平稳，即他的行为并不会随时间的推移而变化，这样我们就可以通过该序列的过去行为预测未来。

ARIMA模型的基本步骤：

  1. 绘制原始序列，进行单位根ADF检验，观察是否平稳。现实中很多都是不平稳的序列，故一般不要进行阶差分，转化为平稳时间序列。
  2. 根据得到的平稳序列，要分别求平稳序列的自相关系数以及偏相关系数，结合自相关图和偏自相关图得到最佳的、。
  3. 上面已经得到了，就可以建立模型，最后对模型进行检验。

02.导入库，加载数据

安装statsmodels库

!pip3 install statsmodels -i https://pypi.tuna.tsinghua.edu.cn/simple

导入相关库及配置

import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import matplotlib as mpl
import datetime

import statsmodels.api as sm
from statsmodels.graphics.api import qqplot
from statsmodels.tsa.arima_model import ARIMA
from statsmodels.tsa.seasonal import seasonal_decompose
from statsmodels.tsa.stattools import adfuller as ADF

## 设置字符集，防止中文乱码
mpl.rcParams['font.sans-serif']=[u'simHei']
mpl.rcParams['axes.unicode_minus']=False

引入库说明：

其他的很熟悉了，就不用说了，介绍不认识的库：

  1. seasonal_decompose：用来将原始序列分解，认为原始序列是趋势序列、残差序列、季节性序列综合影响的结果，这里将各个序列分解出来单独分析。
  2. qqplot：用来绘制Q-Q图，所有点都能很好的拟合在对角线上，认为序列服从正态分布
  3. ADF：迪克-富勒（ADF）检验，也叫单根检验。检验序列中是否存在单位根，如果存在单位根则认为该序列是非平稳序列。
  4. ARIMA就是我们的ARIMA模型库了。
data=pd.read_csv("/share/datasets/arima_time_series.txt")
data.info()

加载数据，查看数据发现，总共两列时间序列Date和销售额Orders，总共149个观测，值得注意的是Date并非时间类型，这里我们需要处理，转换为时间类型。

03.数据处理

1、数据类型转换，将字符时间转化为时间类型

fx=lambda x: pd.datetime.strptime(x,"%Y/%m/%d")
data["Date"]=data["Date"].apply(fx)

data.info()

2、将DataFrame转为Series，重置数据索引为世时间序列

df=data.set_index("Date") #将Date设置为索引
df=pd.Series(df.Orders,index=df.index )  #将销售额转化为Series,索引这设置为Date
df

以上，我们做了两步处理，首先将字符时间用lambda表达式转化为时间类型，第二步，将数据框转化为以时间为索引的Series序列。

04 绘制原始序列趋势

fig=plt.figure(figsize=(10,4))
ax=fig.add_subplot(111)
plt.plot(df)
ax.set(title="Orders of Date",
      ylabel="Orders",
      xlabel="Date")
plt.show()

由上图我们看到，公司2018-01-01至2018-05-29期间的销售额情况，1月回落，2月有所好转，3月回落震荡，4月回升之后又下降，5月回暖之后又下降然后回升。整体在稳定震荡，可能是业绩不佳后紧急采用措施。

ARIMA要求序列平稳，我们需要将原始序列进行差分处理，让他变得平稳，然后再做模型。

05 做差分，平稳性检验

差分，单位根检验

## 做差分，检查平稳性
def diff(timeseries):
    time_diff1=timeseries.diff(1).fillna(0) #1阶差分
    time_diff2=time_diff1.diff(1).fillna(0) #2阶差分

    time_adf=ADF(timeseries)
    time_diff1_adf=ADF(time_diff1)
    time_diff2_adf=ADF(time_diff2)

    return [time_diff1_adf,time_diff2_adf]

diff(df)

【结果解析】： 差分我们做了，并做了单位根检验，输出的检验的结果，但是如何确定序列是否平稳呢？ 首先拿1阶差分的检验结果看：

(-3.683025128820224,
  0.004358356299291195,
  10,
  138,
  {'1%': -3.47864788917503,
   '5%': -2.882721765644168,
   '10%': -2.578065326612056},
  1423.5325819802563),

1、其中-5.871902924306193就是计算的到的ADF值，1%、5%、10%分别表示不同程度的拒绝原假设的统计值，1%表示极为显著地拒绝原假设，5%拒绝原假设，10%勉强拒绝原假设。ADF Test result的值同时小于3个level即说明非常好地拒绝原假设。我们的到的ADF值和level 1相差不大，但是均小于level 2和level 3，可以认为数据平稳。

2、的值为0.004358356299291195，接近0。

1阶差分检验结果认为数据是平稳的。

自回归系数和偏自回归系数

def autocorr(time_series,lags):
    fig=plt.figure(figsize=(12,8))
    ax1=fig.add_subplot(211)
    sm.graphics.tsa.plot_acf(time_series,lags=lags,ax=ax1)

    ax2=fig.add_subplot(212)
    sm.graphics.tsa.plot_pacf(time_series,lags=lags,ax=ax2)

    plt.show()

time_diff1=df.diff(1).fillna(0) 
autocorr(time_diff1,30)

Python ARIMA时间序列模型预测销量

q的确定是ACFq阶后趋于0，落在95%的的置信区间内；p类似，PACFp阶后截尾。由上可知，$p取3，q取5。以上我们p、d（1阶差分）、q都确定了，可以构建模型了。

【知识扩展】相关系数、自相关系数、偏自相关系数

1、相关系数

协方差：两个变量，的值与其各自均值的偏离程度，

相关系数（皮尔逊相关系数），衡量两个随机变量之间线性相关程度的指标。

，表示和之间是正相关

，表示和之前是负相关

，表示和之间不存在任何线性相关性。

2、自相关系数ACF

自相关和两个样本的相关类似，不过只是同一件事物在两个不同时期之间的相关程度，计算公式如下：

对于序列

06 网格搜索确定最优的p和q

上面我们通过观察请确定了p和q的值，通过拖尾或者截尾对模型定阶，有很强的主观性。我们对模型参数的估计方法是通过损失函数和正则项的加权进行评估，平衡预测误差和模型复杂度，我们可以根据信息准则来确定模型阶数。

赤池信息准则：

L表示模型的极大似然函数，K表示模型参数的个数

赤池信息准则认为：

参数越少，AIC的值越小，模型越好

样本越多，AIC的值越小，模型越好。

当样本量很大时，AIC拟合误差的信息随之放大，不能收敛于真实模型。贝叶斯准则弥补了AIC的不足。

贝叶斯准则：

为样本量

网格交叉验证选择最优参数

以下操作非常耗时，需耐心等待

data_eva=sm.tsa.arma_order_select_ic(df,ic=["aic","bic"],trend="nc",max_ar=7,max_ma=7)
print("data_AIC",data_eva.aic_min_order)
print("data_BIC",data_eva.bic_min_order)

data_AIC (7, 7)
data_BIC (1, 1)

如上，我们定义最大和最大的阶数为7，进行网格搜索，得到两个极端值，AIC确定的p和q都为7，BIC确定的都为1，较为尴尬。

不过没关系，至少给我们指定了和的方向，无非就是7和1的组合，加上我们通过拖尾截尾得到的p和q，通过验证，选择最优参数。

arma_77=sm.tsa.SARIMAX(df,order=(7,1,7)).fit()
print("arma_77",arma_77.aic,arma_77.bic,arma_77.hqic)

arma_71=sm.tsa.SARIMAX(df,order=(1,1,7)).fit()
print("arma_71",arma_71.aic,arma_71.bic,arma_71.hqic)

arma_11=sm.tsa.SARIMAX(df,order=(1,1,1)).fit()
print("arma_11",arma_11.aic,arma_11.bic,arma_11.hqic)

arma_35=sm.tsa.SARIMAX(df,order=(3,1,5)).fit()
print("arma_35",arma_35.aic,arma_35.bic,arma_35.hqic)

我们通过组合组合，得到四组p和q的组合，分别验证AIC、BIC和HQIC，综合考虑AIC、BIC、HQIC同时取最小，最后选择。

07 模型检验

折腾这么久，我们确定了、、的值，建立ARIMA模型，但是效果怎么样，我们需要验证：

首先，绘制组合图，观察

arma_77.plot_diagnostics(figsize=(12,8))

我们关心的是，残差（预测值和实际值的差）是不是不相关，是不是均值为0或者正态分布

1、左上角，标准化的残差图也在0轴附近，没有明显趋势性和季节性特征。

2、右上图中概率核密度曲线图，我们看到红色KDE线与N(0,1)是正态分布（均值0 ，标准差为1），残差符合正态分布。

3、左下角Q-Q图，所有点也基本拟合在对角线，也说明符合正态分布。

4、右下角相关系数图，都在95%的置信区间，说明时间序列残差与其本身的滞后项具有相对低的相关性。

杜宾-瓦特森检验（D-W检验）

检验变量的1阶自相关性，一般来说，越接近2越好，说明变量的自相关不明显，模型设计的也越好。

# D-W检验
# DW趋近2,P=0,不存在自相关性
print(sm.stats.durbin_watson(arma_77.resid.values))

上面我们得到的W值为1.8几，趋近2，说明自变量的 1阶相关性不明显，模型效果还可以

D-W检验参考：https://wenku.baidu.com/view/597a3c2e5022aaea988f0fb9.html

残差自相关系数ACF和偏相关系数PACF检验

resid=arma_77.resid
fig=plt.figure(figsize=(16,12))
ax1=fig.add_subplot(211)
sm.graphics.tsa.plot_acf(resid,lags=15,ax=ax1)  #自相关系数

ax2=fig.add_subplot(212)
sm.graphics.tsa.plot_pacf(resid,lags=15,ax=ax2)  #偏相关系数

自相关图和偏自相关图，同我们

plot_diagnostics得到的结果一致，基本在95%的置信区间，

说明时间序列残差与其本身的滞后项具有相对低的相关性或者不相关。

  • 混成检验/Ljung-Box test(LB检验)

  1. 当（一般都用1%, 5%, 10%）, 拒绝原假设，结果显著，序列相关； 当，接受原假设，结果不显著，序列不相关，认为是白噪序列。
acf,q,p=sm.tsa.acf(resid.values.squeeze(),nlags=20,qstat=True)
data=np.c_[range(1,21),acf[1:],q,p]
table=pd.DataFrame(data,columns=["lag","AC","Q","P-value"])
print(table.set_index("lag"))

运行结果，得到的值均大于，接受原假设，拒绝，该序列是独立的，所以残差序列不存在自相关性。

综合以上，我们确定，，的阶数，以及检验残差和本身得滞后项不是自相关的，且满足正态分布，不随趋势和季节波动，是平稳序列。接下来我们进行预测。

08 模型结果预测

pre=arma_77.predict("2018-05-30","2018-06-10",dynamic=True)
pre

#绘制预测曲线图
fig,ax=plt.subplots(figsize=(12,8))
ax=df.ix["2018-01-01":].plot(ax=ax)
fig=arma_77.predict("2018-05-30","2018-06-10",dynamic=True,ax=ax,plot_insample=False).plot(style="r-.")
plt.title("未来10天的销售额预测",fontsize=20)
plt.show()

预测的趋势和原趋势是不是比较吻合？最后，我们将预测的数据和原数据合并，输出为“series_out2.csv”。

df2=df.append(pre)
df2.to_csv("series_out2.csv",header=["order"])

总结

ARIMA模型预测信息量很大，也涉及的知识或者技术，总想串起来，这样有理有据，可深可广，学习效率也高。也请支持作者。记得"转发"、“收藏"、“在看”暖心三连，原创不易，谢谢支持。

移动端设备除iPad Pro外，其它移动设备仅能阅读基础的文本文字。
建议使用PC或笔记本电脑，浏览器使用Chrome或FireFox进行浏览，以开启左侧互动实验区来提升学习效率，推荐使用的分辨率为1920x1080或更高。
我们坚信最好的学习是参与其中这一理念，并致力成为中文互联网上体验更好的学练一体的IT技术学习交流平台。
您可加QQ群：575806994，一起学习交流技术,反馈网站使用中遇到问题。
内容、课程、广告等相关合作请扫描右侧二维码添加好友。

狐狸教程 Copyright 2021

进入全屏