数据源:arima_data.xls
数据结果
ARIMA模型包
# 使用该模型需要一些列判别操作
# 平稳性检测 -》 白噪声检验 ——》 是否差分 -》 AIC和BIC指标 --> 模型定阶 --》预测
# acf() 计算自相关系数
# autocorr = acf(data,unbiased = False,nlags = 40,qstat = False,fft = False ,alpha = None)
# data_观测值序列(时间序列,dataFrame or Series)
# autocorr_观测值序列自相关函数
# 其他为可选参数,如 qstat = True,同事返回Q统计量和对应的p值
# plot_acf()/plot_pacf() 画自相关系数图和偏相关系数图
# p = plot_acf(data)
# 返回一个Matplotlib对象,可以用show()方法显示图像
# adfuller() 对观测值序列进行单位根检验(ADF test)
# h = adffuller(Series,maxlag = None,regression = 'c',autolag = 'AIC',store = False,regresults = False)
# 输入参数Series为一维观测值序列,依次返回adf,pvalue,usedlag,nobs,critical,values,icbest,regresults,resstore
# diff() 对观测值序列进行差分计算
# D。diff()_D为dataframe 或者 Series
# arima 设置时序模式的建模参数,创建ARIMA时序模型
# arima = ARIMA(data,(p,d,q)).fit()
# data参数为输入时间序列,d为差分次数
# summary()/summary2() 生成已有模型报告
# 包含模型系数,标准误差,p值,AIC和BIC值等
# aic/bic/hqic 计算模型的AIC、BIC、HQIC值
# arima.aic/arima.bic/arima.hqic
# forecast() 预测
# a,b,c = arima.forecase(num)
# num为预测的天数,a_返回的预测值,b_预测的误差,c_预测值置信区间
# acorr_ljungbox() 检测是否为白噪声序列
# acorr_ljungbox(data,lags=1)
# lags 为滞后数,返回统计量和p值
源代码
import pandas as pd
import statsmodels.api as sm
import matplotlib.pyplot as plt
#参数初始化
discfile = 'F:/python 数据挖掘分析实战/Data/arima_data.xls'
forecastnum = 5
#读取数据,指定日期列为指标,Pandas自动将“日期”列识别为Datetime格式
data = pd.read_excel(discfile, index_col = u'日期')
#时序图
import matplotlib.pyplot as plt
plt.rcParams['font.sans-serif'] = ['SimHei'] #用来正常显示中文标签
plt.rcParams['axes.unicode_minus'] = False #用来正常显示负号
data.plot()
plt.show()
#自相关图
from statsmodels.graphics.tsaplots import plot_acf
plot_acf(data).show()
# 从该图说明自相关系数长期大于0,序列间有很强的长期相关性
#平稳性检测
from statsmodels.tsa.stattools import adfuller as ADF
print(u'原始序列的ADF检验结果为:', ADF(data[u'销量']))
#返回值依次为adf、pvalue、usedlag、nobs、critical values、icbest、regresults、resstore
# P值显著大于0.05,该序列是非平稳的时间序列
#差分后的结果,对原始数据进行一阶拆分,并进行平稳性和白噪声检验
D_data = data.diff().dropna()
D_data.columns = [u'销量差分']
D_data.plot() #时序图
plt.show()
plot_acf(D_data).show() #自相关图,在均值附近比较平稳的波动,自相关图有很强的短期性,显示一阶截尾
print(u'原始序列的ADF检验结果为:', ADF(D_data[u'销量差分']))
# 单位根检验P值小于0.05,所以一阶差之后的序列是平稳序列
#白噪声检验
from statsmodels.stats.diagnostic import acorr_ljungbox
print(u'差分序列的白噪声检验结果为:', acorr_ljungbox(D_data, lags=1)) #返回统计量和p值
#输出的P值远小于0.05,所以一阶差分之后的序列是平稳非白噪声序列
from statsmodels.graphics.tsaplots import plot_pacf
plot_pacf(D_data).show() #偏自相关图,显示出拖尾性,本来从自相关图和偏自相关图可以看出p,q...看不懂
"""
#第二步、确定p、q,---》 0,1
res = sm.tsa.arma_order_select_ic(
D_data.dropna(),
max_ar=8,
max_ma=8,
ic=['aic', 'bic', 'hqic'],
trend='nc'
)
"""
from statsmodels.tsa.arima_model import ARIMA
data[u'销量'] = data[u'销量'].astype(float)
#定阶
pmax = int(len(D_data)/10) #一般阶数不超过length/10
qmax = int(len(D_data)/10) #一般阶数不超过length/10
bic_matrix = [] #bic矩阵
for p in range(pmax+1):
tmp = []
for q in range(qmax+1):
try: #存在部分报错,所以用try来跳过报错。
tmp.append(ARIMA(data, (p,1,q)).fit().bic)
except:
tmp.append(None)
bic_matrix.append(tmp)
bic_matrix = pd.DataFrame(bic_matrix) #从中可以找出最小值
p,q = bic_matrix.stack().idxmin() #先用stack展平,然后用idxmin找出最小值位置。
print(u'BIC最小的p值和q值为:%s、%s' %(p,q))
model = ARIMA(data, (p,1,q)).fit() #建立ARIMA(0, 1, 1)模型
model.summary2() #给出一份模型报告
model.forecast(5) #作为期5天的预测,返回预测结果、标准误差、置信区间。
参考资料:《Python数据分析与挖掘实战》