本文關鍵詞：時間序列分析

  更多相關文章： python 時間序列 分析 ARIMA

1 時間序列與時間序列分析

在生產和科學研究中，對某一個或者一組變量 所得到的離散數(shù)字組成的序列集合，稱之為時間序列。
時間序列分析是根據系統(tǒng)觀察得到的時間序列數(shù)據，通過曲線擬合和參數(shù)估計來建立數(shù)學模型的理論和方法。時間序列分析常用于國民宏觀經濟控制、市場潛力預測、氣象預測、農作物害蟲災害預報等各個方面。

2 時間序列建�；静襟E 3 ARIMA實戰(zhàn)解剖

原理大概清楚，實踐卻還是會有諸多問題。相比較R語言，Python在做時間序列分析的資料相對少很多。下面就通過Python語言詳細解析后三個步驟的實現(xiàn)過程。
文中使用到這些基礎庫： 。 對其調用如下

from __future__ import print_function import pandas as pd import numpy as np from scipy import stats import matplotlib.pyplot as plt import statsmodels.api as sm from statsmodels.graphics.api import qqplot 3.1 獲取數(shù)據

這里我們使用一個具有周期性的測試數(shù)據，進行分析。
數(shù)據如下：

dta=[10930,10318,10595,10972,7706,6756,9092,10551,9722,10913,11151,8186,6422,
6337,11649,11652,10310,12043,7937,6476,9662,9570,9981,9331,9449,6773,6304,9355,
10477,10148,10395,11261,8713,7299,10424,10795,11069,11602,11427,9095,7707,10767,
12136,12812,12006,12528,10329,7818,11719,11683,12603,11495,13670,11337,10232,
13261,13230,15535,16837,19598,14823,11622,19391,18177,19994,14723,15694,13248,
9543,12872,13101,15053,12619,13749,10228,9725,14729,12518,14564,15085,14722,
11999,9390,13481,14795,15845,15271,14686,11054,10395]

dta=pd.Series(dta) dta.index = pd.Index(sm.tsa.datetools.dates_from_range('2001','2100')) dta.plot(figsize=(12,8))

3.2 時間序列的差分

ARIMA 模型對時間序列的要求是平穩(wěn)型。因此，當你得到一個非平穩(wěn)的時間序列時，首先要做的即是做時間序列的差分，直到得到一個平穩(wěn)時間序列。如果你對時間序列做d次差分才能得到一個平穩(wěn)序列，那么可以使用ARIMA(p,d,q)模型，其中d是差分次數(shù)。

fig = plt.figure(figsize=(12,8)) ax1= fig.add_subplot(111) diff1 = dta.diff(1) diff1.plot(ax=ax1)

一階差分的時間序列的均值和方差已經基本平穩(wěn)，不過我們還是可以比較一下二階差分的效果

fig = plt.figure(figsize=(12,8)) ax2= fig.add_subplot(111) diff2 = dta.diff(2) diff2.plot(ax=ax2)

可以看出二階差分后的時間序列與一階差分相差不大，并且二者隨著時間推移，時間序列的均值和方差保持不變。因此可以將差分次數(shù)d設置為1。
其實還有針對平穩(wěn)的檢驗，叫“ADF單位根平穩(wěn)型檢驗”，以后再更。

3.3 合適的

現(xiàn)在我們已經得到一個平穩(wěn)的時間序列，接來下就是選擇合適的ARIMA模型，即ARIMA模型中合適的。
第一步我們要先檢查平穩(wěn)時間序列的自相關圖和偏自相關圖。

dta= dta.diff(1)#我們已經知道要使用一階差分的時間序列，之前判斷差分的程序可以注釋掉 fig = plt.figure(figsize=(12,8)) ax1=fig.add_subplot(211) fig = sm.graphics.tsa.plot_acf(dta,lags=40,ax=ax1) ax2 = fig.add_subplot(212) fig = sm.graphics.tsa.plot_pacf(dta,lags=40,ax=ax2)

其中l(wèi)ags 表示滯后的階數(shù)，以上分別得到acf 圖和pacf 圖

通過兩圖觀察得到：
* 自相關圖顯示滯后有三個階超出了置信邊界；
* 偏相關圖顯示在滯后1至7階（lags 1,2,…，7）時的偏自相關系數(shù)超出了置信邊界，從lag 7之后偏自相關系數(shù)值縮小至0
則有以下模型可以供選擇：
1. ARMA(0,1)模型：即自相關圖在滯后1階之后縮小為0，且偏自相關縮小至0，則是一個階數(shù)q=1的移動平均模型；
2. ARMA(7,0)模型：即偏自相關圖在滯后7階之后縮小為0，且自相關縮小至0，則是一個階層p=3的自回歸模型；
3. ARMA(7,1)模型：即使得自相關和偏自相關都縮小至零。則是一個混合模型。
4. …還可以有其他供選擇的模型
現(xiàn)在有以上這么多可供選擇的模型，我們通常采用ARMA模型的AIC法則。我們知道：增加自由參數(shù)的數(shù)目提高了擬合的優(yōu)良性，AIC鼓勵數(shù)據擬合的優(yōu)良性但是盡量避免出現(xiàn)過度擬合(Overfitting)的情況。所以優(yōu)先考慮的模型應是AIC值最小的那一個。赤池信息準則的方法是尋找可以最好地解釋數(shù)據但包含最少自由參數(shù)的模型。不僅僅包括AIC準則，目前選擇模型常用如下準則：
* AIC=-2 ln(L) + 2 k 中文名字：赤池信息量 akaike information criterion
* BIC=-2 ln(L) + ln(n)*k 中文名字：貝葉斯信息量 bayesian information criterion
* HQ=-2 ln(L) + ln(ln(n))*k hannan-quinn criterion
構造這些統(tǒng)計量所遵循的統(tǒng)計思想是一致的，就是在考慮擬合殘差的同時，依自變量個數(shù)施加“懲罰”。但要注意的是，這些準則不能說明某一個模型的精確度，也即是說，對于三個模型Ａ，Ｂ，Ｃ，我們能夠判斷出Ｃ模型是最好的，但不能保證Ｃ模型能夠很好地刻畫數(shù)據，因為有可能三個模型都是糟糕的。

arma_mod20 = sm.tsa.ARMA(dta,(7,0)).fit() print(arma_mod20.aic,arma_mod20.bic,arma_mod20.hqic) arma_mod30 = sm.tsa.ARMA(dta,(0,1)).fit() print(arma_mod30.aic,arma_mod30.bic,arma_mod30.hqic) arma_mod40 = sm.tsa.ARMA(dta,(7,1)).fit() print(arma_mod40.aic,arma_mod40.bic,arma_mod40.hqic) arma_mod50 = sm.tsa.ARMA(dta,(8,0)).fit() print(arma_mod50.aic,arma_mod50.bic,arma_mod50.hqic)

可以看到ARMA(7,0)的aic，bic，hqic均最小，因此是最佳模型。

3.4 模型檢驗

在指數(shù)平滑模型下，觀察ARIMA模型的殘差是否是平均值為0且方差為常數(shù)的正態(tài)分布（服從零均值、方差不變的正態(tài)分布），同時也要觀察連續(xù)殘差是否（自）相關。

3.4.1 我們對ARMA(7,0)模型所產生的殘差做自相關圖 fig = plt.figure(figsize=(12,8)) ax1 = fig.add_subplot(211) fig = sm.graphics.tsa.plot_acf(resid.values.squeeze(), lags=40, ax=ax1) ax2 = fig.add_subplot(212) fig = sm.graphics.tsa.plot_pacf(resid, lags=40, ax=ax2)

3.4.2 做D-W檢驗

德賓-沃森（Durbin-Watson）檢驗。德賓-沃森檢驗,簡稱D-W檢驗，是目前檢驗自相關性最常用的方法，但它只使用于檢驗一階自相關性。因為自相關系數(shù)ρ的值介于-1和1之間，所以 0≤DW≤４。并且DW＝O＝＞ρ＝１　　 即存在正自相關性
DW＝４＜＝＞ρ＝－１　即存在負自相關性
DW＝２＜＝＞ρ＝０　　即不存在（一階）自相關性
因此，當DW值顯著的接近于O或４時，則存在自相關性，而接近于２時，則不存在（一階）自相關性。這樣只要知道ＤＷ統(tǒng)計量的概率分布，在給定的顯著水平下，根據臨界值的位置就可以對原假設Ｈ０進行檢驗。

print(sm.stats.durbin_watson(arma_mod20.resid.values))

檢驗結果是2.02424743723，說明不存在自相關性。

3.4.3 觀察是否符合正態(tài)分布

這里使用QQ圖，它用于直觀驗證一組數(shù)據是否來自某個分布，或者驗證某兩組數(shù)據是否來自同一（族）分布。在教學和軟件中常用的是檢驗數(shù)據是否來自于正態(tài)分布。QQ圖細節(jié)，下次再更。

resid = arma_mod20.resid#殘差 fig = plt.figure(figsize=(12,8)) ax = fig.add_subplot(111) fig = qqplot(resid, line='q', ax=ax, fit=True)

3.4.4 Ljung-Box檢驗

Ljung-Box test是對randomness的檢驗,或者說是對時間序列是否存在滯后相關的一種統(tǒng)計檢驗。對于滯后相關的檢驗，我們常常采用的方法還包括計算ACF和PCAF并觀察其圖像，但是無論是ACF還是PACF都僅僅考慮是否存在某一特定滯后階數(shù)的相關。LB檢驗則是基于一系列滯后階數(shù)，判斷序列總體的相關性或者說隨機性是否存在。
時間序列中一個最基本的模型就是高斯白噪聲序列。而對于ARIMA模型，其殘差被假定為高斯白噪聲序列，所以當我們用ARIMA模型去擬合數(shù)據時，擬合后我們要對殘差的估計序列進行LB檢驗，，判斷其是否是高斯白噪聲，如果不是，那么就說明ARIMA模型也許并不是一個適合樣本的模型。

r,q,p = sm.tsa.acf(resid.values.squeeze(), qstat=True) data = np.c_[range(1,41), r[1:], q, p] table = pd.DataFrame(data, columns=['lag', "AC", "Q", "Prob(>Q)"]) print(table.set_index('lag'))

檢驗的結果就是看最后一列前十二行的檢驗概率（一般觀察滯后1~12階），如果檢驗概率小于給定的顯著性水平，比如0.05、0.10等就拒絕原假設，其原假設是相關系數(shù)為零。就結果來看，如果取顯著性水平為0.05，那么相關系數(shù)與零沒有顯著差異，即為白噪聲序列。

3.5 模型預測

模型確定之后，就可以開始進行預測了，我們對未來十年的數(shù)據進行預測。

predict_sunspots = arma_mod20.predict('2090', '2100', dynamic=True) print(predict_sunspots) fig, ax = plt.subplots(figsize=(12, 8)) ax = dta.ix['2001':].plot(ax=ax) predict_sunspots.plot(ax=ax)

前面90個數(shù)據為測試數(shù)據，最后10個為預測數(shù)據；從圖形來，預測結果較為合理。至此，本案例的時間序列分析也就結束了。

參考文獻與推薦閱讀

本文編號：871998