數(shù)據(jù)源：arima_data.xls
arima_data.xls
數(shù)據(jù)結(jié)果
Paste_Image.png
ARIMA模型包
# 使用該模型需要一些列判別操作
# 平穩(wěn)性檢測 -》 白噪聲檢驗(yàn) ——》 是否差分 -》 AIC和BIC指標(biāo) --> 模型定階 --》預(yù)測
# acf() 計算自相關(guān)系數(shù)
# autocorr = acf(data,unbiased = False,nlags = 40,qstat = False,fft = False ,alpha = None)
# data_觀測值序列(時間序列,dataFrame or Series)
# autocorr_觀測值序列自相關(guān)函數(shù)
# 其他為可選參數(shù)律适，如 qstat = True,同事返回Q統(tǒng)計量和對應(yīng)的p值
# plot_acf()/plot_pacf() 畫自相關(guān)系數(shù)圖和偏相關(guān)系數(shù)圖
# p = plot_acf(data)
# 返回一個Matplotlib對象谅阿，可以用show()方法顯示圖像
# adfuller() 對觀測值序列進(jìn)行單位根檢驗(yàn)(ADF test)
# h = adffuller(Series,maxlag = None,regression = 'c',autolag = 'AIC',store = False,regresults = False)
# 輸入?yún)?shù)Series為一維觀測值序列，依次返回adf,pvalue,usedlag,nobs,critical,values,icbest,regresults,resstore
# diff() 對觀測值序列進(jìn)行差分計算
# D诡蜓。diff()_D為dataframe 或者 Series
# arima 設(shè)置時序模式的建模參數(shù)蕊肥，創(chuàng)建ARIMA時序模型
# arima = ARIMA(data,(p,d,q)).fit()
# data參數(shù)為輸入時間序列膘茎，d為差分次數(shù)
# summary()/summary2() 生成已有模型報告
# 包含模型系數(shù)漆弄，標(biāo)準(zhǔn)誤差，p值弄砍，AIC和BIC值等
# aic/bic/hqic 計算模型的AIC颅筋、BIC、HQIC值
# arima.aic/arima.bic/arima.hqic
# forecast() 預(yù)測
# a,b,c = arima.forecase(num)
# num為預(yù)測的天數(shù),a_返回的預(yù)測值,b_預(yù)測的誤差,c_預(yù)測值置信區(qū)間
# acorr_ljungbox() 檢測是否為白噪聲序列
# acorr_ljungbox(data,lags=1)
# lags 為滯后數(shù),返回統(tǒng)計量和p值
源代碼
import pandas as pd
import statsmodels.api as sm
import matplotlib.pyplot as plt
#參數(shù)初始化
discfile = 'F:/python 數(shù)據(jù)挖掘分析實(shí)戰(zhàn)/Data/arima_data.xls'
forecastnum = 5
#讀取數(shù)據(jù)输枯，指定日期列為指標(biāo)，Pandas自動將“日期”列識別為Datetime格式
data = pd.read_excel(discfile, index_col = u'日期')
#時序圖
import matplotlib.pyplot as plt
plt.rcParams['font.sans-serif'] = ['SimHei'] #用來正常顯示中文標(biāo)簽
plt.rcParams['axes.unicode_minus'] = False #用來正常顯示負(fù)號
data.plot()
plt.show()
#自相關(guān)圖
from statsmodels.graphics.tsaplots import plot_acf
plot_acf(data).show()
# 從該圖說明自相關(guān)系數(shù)長期大于0占贫，序列間有很強(qiáng)的長期相關(guān)性
#平穩(wěn)性檢測
from statsmodels.tsa.stattools import adfuller as ADF
print(u'原始序列的ADF檢驗(yàn)結(jié)果為：', ADF(data[u'銷量']))
#返回值依次為adf桃熄、pvalue、usedlag、nobs瞳收、critical values碉京、icbest、regresults螟深、resstore
# P值顯著大于0.05谐宙，該序列是非平穩(wěn)的時間序列
#差分后的結(jié)果，對原始數(shù)據(jù)進(jìn)行一階拆分界弧，并進(jìn)行平穩(wěn)性和白噪聲檢驗(yàn)
D_data = data.diff().dropna()
D_data.columns = [u'銷量差分']
D_data.plot() #時序圖
plt.show()
plot_acf(D_data).show() #自相關(guān)圖,在均值附近比較平穩(wěn)的波動凡蜻，自相關(guān)圖有很強(qiáng)的短期性，顯示一階截尾
print(u'原始序列的ADF檢驗(yàn)結(jié)果為：', ADF(D_data[u'銷量差分']))
# 單位根檢驗(yàn)P值小于0.05垢箕，所以一階差之后的序列是平穩(wěn)序列
#白噪聲檢驗(yàn)
from statsmodels.stats.diagnostic import acorr_ljungbox
print(u'差分序列的白噪聲檢驗(yàn)結(jié)果為：', acorr_ljungbox(D_data, lags=1)) #返回統(tǒng)計量和p值
#輸出的P值遠(yuǎn)小于0.05划栓，所以一階差分之后的序列是平穩(wěn)非白噪聲序列
from statsmodels.graphics.tsaplots import plot_pacf
plot_pacf(D_data).show() #偏自相關(guān)圖，顯示出拖尾性条获，本來從自相關(guān)圖和偏自相關(guān)圖可以看出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) #一般階數(shù)不超過length/10
qmax = int(len(D_data)/10) #一般階數(shù)不超過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天的預(yù)測萍悴，返回預(yù)測結(jié)果头遭、標(biāo)準(zhǔn)誤差、置信區(qū)間癣诱。
參考資料：《Python數(shù)據(jù)分析與挖掘?qū)崙?zhàn)》
