arima预测建模步骤python_ARIMA模型构建、预测——基于Python

1 #-*- coding: utf-8 -*-

2 """

3 Created on Thu Aug 23 20:00:31 20184

5 @author: Luove6 """

7

8 importpandas as pd9 importnumpy as np10 from statsmodels.tsa.arima_model importARIMA11 from statsmodels.tsa.stattools importadfuller as ADF12 from statsmodels.stats.diagnostic importacorr_ljungbox13

14 file1='../data/discdata.xls'

15 #file2='../data/discdata_processed.xls'

16 data=pd.read_excel(file1)17

18 data=data[data['TARGET_ID']==184]19 data_group=data.groupby('COLLECTTIME')20 defattribute_turn(x):21 result=pd.Series(index=['SYS_NAME','CWXT_DB:184:c\\','CWXT_DB:184:d\\','COLLECTTIME'])22 result['SYS_NAME']=x['SYS_NAME'].iloc[0]23 result['CWXT_DB:184:c\\']=x['VALUE'].iloc[0]24 result['CWXT_DB:184:d\\']=x['VALUE'].iloc[1]25 result['COLLECTTIME']=x['COLLECTTIME'].iloc[0]26 returnresult27 data_processed=data_group.apply(attribute_turn)28

29 data_1=data_processed[:-5]30 diff=031 adf=ADF(data_processed['CWXT_DB:184:d\\'])32 while adf[1]>=0.05: #h0:存在单位根，序列不平稳

33 diff+=1

34 adf=ADF(data_processed['CWXT_DB:184:d\\'].diff(diff).dropna())35

36 print('原序列经过%d阶差分平稳，此时p值为%f'%(diff,adf[1]))37 [[lb],p]=acorr_ljungbox(data_1['CWXT_DB:184:d\\'],lags=1) #LB时间序列是否存在滞后相关的一种统计检验，选一阶滞后

38 if p<0.05: #h0序列是白噪声序列

39 print('原序列是非白噪声序列，此时p值为%f'%p)40 else:41 print('原序列是白噪声序列，此时p值为%f'%p)42 [[lb],p]=acorr_ljungbox(data_1['CWXT_DB:184:d\\'].diff().dropna(),lags=1)43 if p<0.05: #h0序列是白噪声序列

44 print('一阶差分序列是非白噪声序列，此时p值为%f'%p)45 else:46 print('一阶差分序列是白噪声序列，此时p值为%f'%p)47

48 xdata=data_1['CWXT_DB:184:d\\']49 pmax=int(len(xdata)/10)50 qmax=int(len(xdata)/10)51 bic_matrix=[]52

53 for p in range(pmax+1):54 tmp=[]55 for q in range(qmax+1): #生成一个在p=p时，q=0,1,2,3,4时的ARIMA模型的BIC，q为之后的列索引

56 try:57 tmp.append(ARIMA(xdata,(p,1,q)).fit().bic)58 except:59 tmp.append(None)60 bic_matrix.append(tmp) #此时的行列索引分别对应p,q

61 bic_matrix=pd.DataFrame(bic_matrix)62 p,q=bic_matrix.stack().idxmin() #将行拆散成列，形成双行索引结构，外层主索引为原行索引内层为列索引转化而来，转化过程默认去除NaN和None

63 print('BIC最小的p和q分别是%d和%d'%(p,q)) #确定了p,i,q，模型设定完毕可以开始回归和预测

64

65 arima=ARIMA(xdata,(p,1,q)).fit()66 xdata_pred=arima.predict(typ='levels')67 pred_error=(xdata_pred-xdata).dropna()68 #arima.resid # 残差

69

70 [lb],p=acorr_ljungbox(pred_error,lags=1) #LB时间序列是否存在滞后相关的一种统计检验，选一阶滞后

71 if p<0.05: #即拒绝原假设，h0是白噪声序列

72 print('ARIMA(0,1,2)残差序列不是白噪声序列，p值为%f'%p)73 else:74 print('ARIMA(0,1,2)残差序列是白噪声序列，p值为%f'%p)75

76 #arima.predict('2014-11-12','2014-11-16',dynamic=True,typ='levels') #预测之后5天的值

77 #ARIMAResults.predict(start=None, end=None, exog=None, typ='linear', dynamic=False)#dynamic，逻辑值，True表样本外预测，默认False样本预测，typ，取值'linear', 'levels'表示根据内生变量的差分做线性预测，预测原数据的水平(源数据的模型预测值)

78 compare_data=pd.concat([data_processed[-5:]['CWXT_DB:184:d\\'],arima.predict('2014-11-12','2014-11-16',dynamic=True,typ='levels')],axis=1)#合并列

79 compare_data.columns=['real','pred']80

81 abs_=(compare_data['real']-compare_data['pred']).abs()82 mae=abs_/mean()#Mean Absolute Error ，平均绝对误差

83 rmse=(abs_**2.mean())**0.5#Root Mean Square Error,均方根误差

84 mape=(abs_/compare_data['real']).mean()#mean absolute percentage error，平均绝对百分比误差

85

86 print('平均绝对误差MAE=%f；\n均方根误差RMSE=%f；\n平均绝对百分比误差MAPE=%f。'%(mae,rmse,mape))87

88 #至此，模型构建、计算、评价已完成；评价即误差在误差设定阈值之内，说明模型可以使用

89 #由此模型来构建预警系统，主要系由历史数据(每日定时采集)运行该模型(主要是P,I,Q的确定)，进而预测未来5天磁盘使用量

90 #由预测的磁盘使用量占总容量比率来划定预警级别，发布预警级别信息

91 #实际中考虑使用容量一般不会突变，每日变化不大，遂该模型调整(主要是P,I,Q的确定)可根据业务情况半月/一个月来修正