淘先锋技术网

首页 1 2 3 4 5 6 7
前段时间整理了一个预测的基本思考框架和常见的方法,其中提到了ARIMA模型,在《大数据预测》那本书里,ARIMA是单独开辟一章讲的,比较复杂和难理解的一个模型,自己最近找了点资料粗浅学习了一下理论,并尝试用Python实现一下。

一、时间序列数据及其预处理

1.时序数据

时序数据顾名思义就是随着时间而变动的数据,是指某个个体在不同时间点上收集到的数据。已经被收集(或者叫观察)到的数据其实是时间序列变量的一个观察值,但由于时间的不可逆性,每一个时间点的变量有且仅能有一个观察值,我们用这些观察值拟合预测模型,用来预测未来时刻时间序列变量的值(此时,未发生的时间序列变量是一个随机变量)。

2.时序数据预处理

时序数据的预处理主要是平稳性和随机性检验,根据检验结果的不同,用不同的模型对数据进行建模。 纯随机序列,又称为白噪声序列,序列的变量之间无任何关系,因此没有分析价值;平稳非白噪声序列,其均值和方差是常数,目前有比较成熟的模型对其进行模拟,主要是AR(Autoregressive model)自回归模型,MA(Movingaverage model)移动平均模型、以及ARMA-一个AR和MA结合的模型

(1)平稳性检验

平稳性检验的基本原理是检验不同时间点时序变量之间的相关关系,由于是一个时间序列数据的不同时刻之间的相互比较,因此取名为自协方差函数和自相关系数,但基本原理和协方差函数、相关系数类似,当一个时间序列有为常数的均值和方差,且其自协方差函数只跟时间间隔有关,即γ(t,s) =γ(k,k+s-t) ,其中γ是自协方差函数。 平稳性检验的方法分为两类,一类是偏主观和模糊的图示检验,一类是比较精确的数学检验。 第一类主要是用时序图和自相关图进行观察,如果一个时间序列的时序图围绕某一个常数上下波动,且波动范围有界,可以认为是平稳的,除此之外,如果自相关图中自相关系数比较快速衰减到0,并且在零附近震荡,可以认为序列是平稳序列,这是因为平稳序列通常具有短期相关性,随着延迟期数的增加,变量之间的相关关系会指数降低。

16ff53a5463cbf38f037c04714587468.png

第二类数据的检验方法是单位根检验。单位根又称为单位圆,指模型的特征方程的根和1的关系,如果模型的特征根小于1,则序列非平稳,如果特征根大于1,则序列平稳。

(2)白噪声检验

白噪声序列的各个序列值之间没有任何关系,即γ(k)=0,在实际中,只要序列变量之间的协方差系数在0值附近波动,就可以认为是白噪声序列了。

二、常用模型

平稳时间序列的建模通常是AR、MA和ARMA模型,它们都是一种多元线性回归模型。

1.AR模型

41e117f86c8c3047188a9c99a0a8321b.png

这个模型的定义是t时刻的序列取值和前p期的序列值有关,且是前p期的线性组合。这个比较容易理解,在通常的印象中,我们大都觉得这次发生的事和最近几次发生的事有关联。 并非所有的AR模型都是平稳的,AR模型满足平稳性的条件是系数|φ|<1; 同时自相关系数ACF呈指数的速度衰减,在0附近波动,表现出拖尾性;其次,由于AR模型里在计算t和t-k期之间的相关性时会受到中间k期的影响,剔除k期影响后的自相关系数称为偏自相关系数PACF,平稳AR模型的PACF具有截尾性(尾巴突然断掉),这是因为yt只跟y(t-1)和y(t-p)之间的序列值有关,当时间间隔k>p的时候,yt和y(t-k)之间就没有直接关系了。

2.MA模型

d6e3667a1b8d06f8fc6167a04cb23509.png

这个模型刚开始不太好理解,某个时刻的序列值跟过去历史q期的白噪声有关?are you kidding me?后来看网上有人说MA模型关注的是AR模型中的随机扰动项的累加,通过扰动项的移动平均能够有效消除预测中的随机波动。 MA模型的自相关系数在q阶后全部等于0,具有截尾性,偏自相关系数在q阶后拖尾,这是因为任何一个可逆的MA模型都是可以转换为无穷阶的AR模型(各种嵌套回去),而无穷阶的AR模型的PACF是拖尾的。

3.ARMA模型

aec0ef6865a93fa6c5312adf9436a261.png

ARMA模型是AR和MA模型的线性组合,特别的,当q=0,它是一个AR(p)模型,如果p=0,它是一个MA(q)模型。 平稳条件下的三类模型的自相关系数ACF和偏自相关系数PACF特性如下:

5f7920740cdd09e35c96ff00d0b4323f.png

利用这些特性可以通过图示的方法判断平稳性。

三、建模步骤

时间序列经过数据预处理被判定为平稳非白噪声序列后,就可以利用ARMA进行建模。 1.计算ACF、PACF; 2.ARMA模型定阶,即确定最优p、q值; 3.估计模型参数、进行参数检验; 4.模型检验; 5.模型预测及评价 在这里提及一下ARIMA模型,ARIMA是指时间序列经过I阶差分之后可以满足平稳非白噪声序列的条件,并利用ARMA为其进行建模。因此,ARIMA和ARMA最大的差异就是前者对原数据序列进行了差分,差分次数一般为1次,2次都不太推荐,因为2次差分之后很可能就成为白噪声序列,别问我怎么知道的(╯︵╰)。

四、Python试一试

会列出一些重点步骤和代码。 1.把数据分成了data_time的模型内数据和用于判断模型效果的data_time_test模型外数据
data_time_origin = pd.read_excel(r'C:\*\arima_data_gmv.xlsx',index_col = '日期')data_time = data_time_origin[:'2017-05-20']  #将5-20号之前的数据作为模型数据data_time_test = data_time_origin['2017-05-21':]

2.看一下时序图和自相关函数图判断是否平稳序列

data_time.plot(figsize=(20,8),fontsize=15)#导入相应的统计模型包from statsmodels.graphics.tsaplots import plot_acf #acf计算from statsmodels.tsa.stattools import adfuller as ADF #单位根检测plot_acf(data_time).show()

1bdb5b02eaf206e8fb431a81f4b9826d.png

时序图上能看到还是有比较明显的趋势的,初步判断不是平稳序列。

3b1f91bf736cac2a438f3986096fd770.png

ACF上看自相关系数长期大于0,且后期又长期小于0,呈现一种对称模式,可以判断不是平稳序列。单位根检验也证实不是平稳序列。
adf_test = ADF(data_time['gmv'])print(adf_test)'''统计值=临界值不能拒绝零假设,存在单位根,序列非平稳。可以看到p-value=0.6377显著大于0.05,另外统计量计算值-1.28显著大于1%、5%、10%的临界值,因此序列存在单位根,序列是非平稳的'''result:(-1.281137399884022, 0.6377571909319921, 5, 75, {'1%': -3.520713130074074, '5%': -2.9009249540740742, '10%': -2.5877813777777776}, 573.5874067738569)

3.进行一阶差分运算看其是否能达到平稳

data_time2 = data_time.diff(periods=1)data_time2 = data_time2.dropna()data_time2.columns = ['gmv差分'] data_time2.plot(figsize=(20,8),fontsize = 15)

1b9f63ef046876086d182ce88d4c105c.png

可以看到基本消除了趋势,数据有稳定的均值。
#自相关图plot_acf(data_time2).show()#偏自相关图from statsmodels.graphics.tsaplots import plot_pacfplot_pacf(data_time2).show() #单位根检验adf_test_d = ADF(data_time2['gmv差分'])print(adf_test_d)result:(-8.458458355369716, 1.591087857743347e-13, 4, 75, {'1%': -3.520713130074074, '5%': -2.9009249540740742, '10%': -2.5877813777777776}, 565.8248717064124)

84293f89d91b0d6192d0f77e06f97cc1.png9ff1adcb85e725bee5881565400f2b23.png

自相关图能看出有一些短期的相关性,单位根检验的结果是统计量-8.45小于临界值,p-value也远远小于0.05,差分之后是平稳序列。
#白噪声检验from statsmodels.stats.diagnostic import acorr_ljungboxprint(acorr_ljungbox(data_time2, lags=[3]))result:(array([13.84640456]), array([0.00312186]))
白噪声检验的p-value 0.003小于0.05,因此拒绝原假设,差分后的时序不是白噪声序列。因此上面的检验过程说明差分后的时间序列是平稳非白噪声序列,可以进行平稳时序建模。

4.模型定阶

上面的PACF图,能看到诡异的拖尾性,前半段截尾,后半段拖尾...看着像是1阶AR模型,具体是不是、行不行我也不知道。 暴力用BIC信息值算一下,找最优的pq值。
from statsmodels.tsa.arima_model import ARIMAdata_time['gmv'] = data_time['gmv'].astype(float)#定阶pmax = int(len(data_time)/10) #一般阶数不超过length/10qmax = int(len(data_time)/10) #一般阶数不超过length/10bic_matrix = [] #bic矩阵for p in range(pmax + 1):    temp = []    for q in range(qmax + 1):        try:            temp.append(ARIMA(data_time,(p,1,q)).fit().bic)        except:            temp.append(None)    bic_matrix.append(temp)bic_matrix = pd.DataFrame(bic_matrix)bic_matrix = bic_matrix.stack()p,q = bic_matrix.idxmin()print('BIC值最小的p和q分别为:p=%s, q=%s'  %(p,q))result:BIC值最小的p和q分别为:p=3, q=3
所以定阶结果竟然是3和3,那就试着用这个模拟一下看看。

5.构建最终模型

arima_model = ARIMA(data_time,(p,1,q)).fit()  arima_model.summary()

83c63fbfdde8d70a36e9a0298e740482.png

强烈怀疑这个模型的有效性... wierd,谁能给指条明路啊。用了80天的数据做训练,然后在p=3,q=3的时候,模型结果只给了系数,没给有效性的信息,摊手。 6.既然模型参数没法检验,只能用预测结果看看效果
forecast_n = pd.DataFrame(arima_model.forecast(11)[0],index = pd.date_range('2017-05-21',periods = 11))forecast_n.columns = ['gmv_fore']forecast_nerror_test = pd.merge(data_time_test,forecast_n,left_index = True, right_index = True)error_test['误差'] = error_test['gmv'] - error_test['gmv_fore']plt.figure(figsize=(15,6))plt.plot(error_test.index,error_test['gmv'],'-*',error_test.index,error_test['gmv_fore'],'-go')plt.legend(loc = 'best', fontsize=20)plt.show()print('预测值与真实值之间的绝对误差为:', np.abs(sum(error_test['gmv']) - sum(error_test['gmv_fore'])))result:预测值与真实值之间的绝对误差为:17.937944428429546

d0e5379f91f576aa0581b45738fccc61.png

其中蓝线是原值,绿线是预测值,看结果还行,把上升和下降的趋势都已经预测出来了,用原值的和同预测值的和做差取绝对值,绝对值17.9,在10天的预测区间内,200多的基础数量级上,这个误差我个人觉得还挺nice啊,比我之前做的那些劳什子玩意强多了,而且还是到天啊,到天!以前都是月度的。

7.做一点模型评价

计量经济学里进行模型评价的指标有三种:一是均方误差MSE、一个是绝对平均误差MAE,另外一个是符号正确百分率。符号正确百分率基本用于经济领域,或者说对增长和下降的趋势预测,只是方向对了就算预测不错了,所以此处用前两种。
'''1.均方误差2.绝对平均误差'''#均方根误差&均方误差rmse = np.sqrt(sum(error_test['误差']**2)/len(error_test['误差']))mse = sum(error_test['误差']**2)/len(error_test['误差'])print('均方根误差为:%.2f' % rmse,'\n','均方误差为:%.2f' % mse)#绝对平均误差mae = sum(abs(error_test['误差']))/len(error_test['误差'])print('绝对平均误差为:%.2f' % mae)result:均方根误差为:18.69 均方误差为:349.30绝对平均误差为:13.10
当然这几个指标是需要在不同模型之间比较才能有优劣关系的,也就是说选用不同的p,q,或者甚至是非ARIMA模型的其它模型预测出来的结果之间的相互比较,但鉴于我今天已经不想搞了,所以就不比较了。但无论怎样,我觉得这个预测的结果还真是挺不错的,对于天的预测竟然趋势全对,而且在200的基数量级上最多是5-30日那天预测误差49,预测多了49万,其它时候模型的效果还是蛮好的,我很满意,嗯,就是这样!

243e3d52cebec4fdbabaea9ba3c7bea7.png

参考资料: 1.中国大学MOOC-中南财经政法大学《时间序列分析》. 2.中国大学MOOC-对外经济贸易大学《计量经济学导论》. 3.张良均 等. Python数据分析与挖掘实战[M].北京:机械工业出版社,2015:119-134. 4.网上众多牛人文章.