ARIMA模型——本質上是error和t-?時刻數據差分的線性模型！！！如果數據序列是非平穩的，並存在一定的增長或下降趨勢，則需要對數據進行差分處理!ARIMA（p，d，q）稱為差分自回歸移動平均模型，AR是自回歸， p為自回歸項； MA為移動平均，q為移動平均項數，d為時間序列成為平穩時所做的差分次數

https://www.cnblogs.com/bradleon/p/6827109.html 文章里寫得非常好，需詳細看。尤其是arima的舉例！

可以看到：ARIMA本質上是error和t-?時刻數據差分的線性模型！！！

ARIMA模型全稱為自回歸積分滑動平均模型(Autoregressive Integrated Moving Average Model,簡記ARIMA)，是由博克思(Box)和詹金斯(Jenkins)於70年代初提出一著名時間序列（Time-series Approach）預測方法 ^[1]  ，所以又稱為Box-Jenkins模型、博克思-詹金斯法。其中ARIMA（p，d，q）稱為差分自回歸移動平均模型，AR是自回歸， p為自回歸項； MA為移動平均，q為移動平均項數，d為時間序列成為平穩時所做的差分次數。所謂ARIMA模型，是指將非平穩時間序列轉化為平穩時間序列，然后將因變量僅對它的滯后值以及隨機誤差項的現值和滯后值進行回歸所建立的模型。ARIMA模型根據原序列是否平穩以及回歸中所含部分的不同，包括移動平均過程（MA）、自回歸過程（AR）、自回歸移動平均過程（ARMA）以及ARIMA過程。

 

優點： 模型十分簡單，只需要內生變量而不需要借助其他外生變量。

缺點：
1.要求時序數據是穩定的（stationary），或者是通過差分化(differencing)后是穩定的。
2.本質上只能捕捉線性關系，而不能捕捉非線性關系。
注意，采用ARIMA模型預測時序數據，必須是穩定的，如果不穩定的數據，是無法捕捉到規律的。比如股票數據用ARIMA無法預測的原因就是股票數據是非穩定的，常常受政策和新聞的影響而波動。

 

嚴謹的定義： 一個時間序列的隨機變量是穩定的，當且僅當它的所有統計特征都是獨立於時間的（是關於時間的常量）。
判斷的方法：

1. 穩定的數據是沒有趨勢(trend)，沒有周期性(seasonality)的; 即它的均值，在時間軸上擁有常量的振幅，並且它的方差，在時間軸上是趨於同一個穩定的值的。
2. 可以使用Dickey-Fuller Test進行假設檢驗。
    

ARIMA模型有三個參數:p,d,q。

• p--代表預測模型中采用的時序數據本身的滯后數(lags) ,也叫做AR/Auto-Regressive項
• d--代表時序數據需要進行幾階差分化，才是穩定的，也叫Integrated項。
• q--代表預測模型中采用的預測誤差的滯后數(lags)，也叫做MA/Moving Average項

ARIMA建模基本步驟

1. 獲取被觀測系統時間序列數據；
2. 對數據繪圖，觀測是否為平穩時間序列；對於非平穩時間序列要先進行d階差分運算，化為平穩時間序列；
3. 經過第二步處理，已經得到平穩時間序列。要對平穩時間序列分別求得其自相關系數ACF 和偏自相關系數PACF，通過對自相關圖和偏自相關圖的分析，得到最佳的階層 p 和階數 q
4. 由以上得到的d、q、p，得到ARIMA模型。然后開始對得到的模型進行模型檢驗。
   具體例子會在另一篇文章中給出。

from:https://www.cnblogs.com/bradleon/p/6827109.html

 

預測程序

ARIMA模型預測的基本程序

（一）根據時間序列的散點圖、自相關函數和偏自相關函數圖以ADF 單位根檢驗其方差、趨勢及其季節性變化規律，對序列的平穩性進行識別。一般來講，經濟運行的時間序列都不是平穩序列。

（二）對非平穩序列進行平穩化處理。 如果數據序列是非平穩的，並存在一定的增長或下降趨勢，則需要對數據進行差分處理，如果數據存在異方差，則需對數據進行技術處理，直到處理后的數據的自相關函數值和偏相關函數值無顯著地異於零。

（三）根據 時間序列模型的識別規則，建立相應的模型。若平穩序列的偏相關函數是截尾的，而 自相關函數是拖尾的，可斷定序列適合 AR模型；若平穩序列的偏相關函數是拖尾的，而自相關函數是截尾的，則可斷定序列適合MA模型； 若平穩序列的偏相關函數和自相關函數均是拖尾的，則序列適合ARMA模型。（截尾是指時間序列的自相關函數（ACF）或偏自相關函數（PACF）在某階后均為0的性質（比如AR的PACF）；拖尾是ACF或PACF並不在某階后均為0的性質（比如AR的ACF）。）

（四）進行參數估計，檢驗是否具有統計意義。

（五）進行假設檢驗，診斷殘差序列是否為白噪聲。

（六）利用已通過檢驗的模型進行預測分析。

二、應用ARIMA模型進行預測

每月稅收數據.可以看作是隨着時間的推移而形成的一個隨機時間序列，通過對該時間序列上稅款值的 隨機性、平穩性以及季節性等因素的分析，將這些單月稅收值之間所具有的相關性或依存關系用數學模型描述出來，從而達到利用過去及現在的稅收值信息來預測未來稅收情況的目的。

(一)對序列取對數和作差分處理，形成穩定隨機序列

ARIMA模型建模的 基本條件是要求待預測的數列滿足平穩的條件，即個體值要圍繞序列均值上下波動，不能有明顯的上升或下降趨勢，如果出現上升或下降趨勢，需要對原始序列進行差分平穩化處理。

(二)模型參數的估計

時間序列預測模塊的自相關分析包括對自 相關系數和偏相關系數的分析，通過對比分析從而實現對時間序列特性的識別。從計算結果可知，自相關函數1步截尾，偏自相關函數2步截尾，自相關函數通過白噪聲檢驗。根據變換數列的自相關函數和偏自相關函數的特點，並經過反復測試，對ARIMA模型的參數進行估計.三個參數定為d=l，p=2和q=l。

對參數進行檢驗。從檢驗結果可知，參數估計全部通過 顯著性檢驗.擬合優度 統計量表中給出了殘差序列的方差(0.063367)和 標准誤差(0.251729)，以及按AIC和SBC標准計算的統計量(9.496798)和(18.54752)，這兩個值都較小，表明對預測模型擬合得較好。從殘差的自相關檢驗結果數據中.可以得知殘差通過白噪聲顯著性檢驗。預測模型最終形式為：(14-0.98284B)(1+0.56103B-2)Z=(1-0.34111B)(1+B)u其中，Z=logX。B為后移算子，u為隨機干擾項。

 

使用ARIMA模型對裙子長度預測

from：https://www.cnblogs.com/ECJTUACM-873284962/p/7379717.html

1、加載數據

skirts <- scan("http://robjhyndman.com/tsdldata/roberts/skirts.dat", skip=5)

str(skirts) head(skirts) boxplot(skirts) length(skirts)

2、把數據轉化為是時間序列

skirts_ts <- ts(skirts, start=c(1886), frequency=1)

1)查看時間序列對應的時間

skirts_ts

2)畫出時間序列圖

plot.ts(skirts_ts)

從圖可知：女人裙子邊緣的直徑做成的時間序列數據，從 1866 年到 1911 年在平均值上是不平穩的

3、做差分得到平穩序列

1)做時間序列的一階差分

skirts_diff <- diff(skirts_ts, differences = 1) plot.ts(skirts_diff)

從一階差分的圖中可以看出，數據仍是不平穩的，繼續差分

2)做時間序列的二階差分

skirts_diff2 <- diff(skirts_ts, differences = 2) plot.ts(skirts_diff2)

二次差分后的時間序列在均值和方差上看起來是平穩了

4、找到合適的ARIMA模型

尋找 ARIMA(p,d,q)中合適的 p 值和 q

1)自相關圖ACF

acf(skirts_diff2, lag.max = 20) acf(skirts_diff2, lag.max = 20, plot = F)

自相關圖顯示滯后1階自相關值基本沒有超過邊界值，雖然5階自相關值超出邊界，那么很可能屬於偶然出現的，而自相關值在其他上都沒有超出顯著邊界， 而且我們可以期望 1 到 20 之間的會偶爾超出 95%的置信邊界。 自相關圖5階后結尾

2)偏相關圖PACF

pacf(skirts_diff2, lag.max = 20) pacf(skirts_diff2, lag.max = 20, plot = F)

偏自相關值選1階后結尾
故我們的ARMIA模型為armia（1,2,5

3)使用auto.arima()函數，自動獲取最佳的ARIMA模型

library(forecast) auto.arima(skirts_ts, ic=c("aicc", "aic", "bic"), trace = T)

Best model: ARIMA(1,2,0)

5、建立ARIMA模型：並對比arima(1, 2, 0)與arima(1, 2, 5)模型

1）arima(1, 2, 0)模型

(skirts_arima <- arima(skirts_ts, order = c(1, 2, 0)))

aic = 391.33

2）arima(1, 2, 5)模型

(skirts_arima <- arima(skirts_ts, order = c(1, 2, 5)))

aic = 381.6

AIC是赤池消息准則SC是施瓦茨准則，當兩個數值最小時，則是最優滯后分布的長度。我們進行模型選擇時，AIC值越小越好。所以arima(1, 2, 5)模型較好

6、預測：預測5年后裙子的邊緣直徑

(skirts_forecast <- forecast.Arima(skirts_arima, h=5, level = c(99.5))) plot.forecast(skirts_forecast)

 

ARIMA 模型實踐

from：https://www.jianshu.com/p/4130bac8ebec

模型具體的理論知識就不再做過多說明了，來個實際的例子吧。

ARIMA 模型對湖北省 GDP 的實證分析及預測

這里的例子是采用了一篇論文的數據，【ARIMA模型在湖北省GDP預測中的應用】，可以去中國知網搜索篇名進行下載。

年份	GDP
1978	151.00
1979	188.46
1980	199.38
...	...
2013	24668.49

數據的平穩性處理及檢驗

這里我們用 Python 對數據進行分析處理建模。

畫出原始數據的時間路徑圖

#-*- coding:utf-8 -*- import pandas as pd import numpy as np import matplotlib.pyplot as plt import statsmodels.api as sm from statsmodels.graphics.tsaplots import acf,pacf,plot_acf,plot_pacf from statsmodels.tsa.arima_model import ARMA time_series = pd.Series([151.0, 188.46, 199.38, 219.75, 241.55, 262.58, 328.22, 396.26, 442.04, 517.77, 626.52, 717.08, 824.38, 913.38, 1088.39, 1325.83, 1700.92, 2109.38, 2499.77, 2856.47, 3114.02, 3229.29, 3545.39, 3880.53, 4212.82, 4757.45, 5633.24, 6590.19, 7617.47, 9333.4, 11328.92, 12961.1, 15967.61]) time_series.index = pd.Index(sm.tsa.datetools.dates_from_range('1978','2010')) time_series.plot(figsize=(12,8)) plt.show() 

原始數據.png

由上圖我們可以看出，這個時間序列是呈指數形式的，波動性比較大，不是穩定的時間序列，一般對於這種指數形式的數據，可以對其取對數，將其轉化為線性趨勢。

time_series = np.log(time_series) time_series.plot(figsize=(8,6)) plt.show() 

對原始數據取對數.png

由上圖可以看出，去了對數之后的時間路徑圖明顯具有線性趨勢，為了確定其穩定性，對取對數后的數據進行 adf 檢驗

t=sm.tsa.stattools.adfuller(time_series, )
output=pd.DataFrame(index=['Test Statistic Value', "p-value", "Lags Used", "Number of Observations Used","Critical Value(1%)","Critical Value(5%)","Critical Value(10%)"],columns=['value']) output['value']['Test Statistic Value'] = t[0] output['value']['p-value'] = t[1] output['value']['Lags Used'] = t[2] output['value']['Number of Observations Used'] = t[3] output['value']['Critical Value(1%)'] = t[4]['1%'] output['value']['Critical Value(5%)'] = t[4]['5%'] output['value']['Critical Value(10%)'] = t[4]['10%'] print(output) 

檢驗結果如下

檢驗項	檢驗結果
Test Statistic Value	0.807369
p-value	0.991754
Lags Used	1
Number of Observations Used	31
Critical Value(1%)	-3.66143
Critical Value(5%)	-2.96053
Critical Value(10%)	-2.61932

由上表可知，t 統計量要大於任何置信度的臨界值，因此認為該序列是非平穩的，所以再對序列進行差分處理，發現差分之后的序列基本達到穩定，如下圖所示，並且通過了 ADF 檢驗，檢驗結果見下表。

time_series = time_series.diff(1)
time_series = time_series.dropna(how=any)
time_series.plot(figsize=(8,6))
plt.show()
t=sm.tsa.stattools.adfuller(time_series)
output=pd.DataFrame(index=['Test Statistic Value', "p-value", "Lags Used", "Number of Observations Used","Critical Value(1%)","Critical Value(5%)","Critical Value(10%)"],columns=['value']) output['value']['Test Statistic Value'] = t[0] output['value']['p-value'] = t[1] output['value']['Lags Used'] = t[2] output['value']['Number of Observations Used'] = t[3] output['value']['Critical Value(1%)'] = t[4]['1%'] output['value']['Critical Value(5%)'] = t[4]['5%'] output['value']['Critical Value(10%)'] = t[4]['10%'] print(output) 

差分.png

檢驗項	檢驗結果
Test Statistic Value	-3.52276
p-value	0.00742139
Lags Used	0
Number of Observations Used	31
Critical Value(1%)	-3.66143
Critical Value(5%)	-2.96053
Critical Value(10%)	-2.61932

確定自相關系數和平均移動系數（p,q）

根據時間序列的識別規則，采用 ACF 圖、PAC 圖，AIC 准則（赤道信息量准則）和 BIC 准則（貝葉斯准則）相結合的方式來確定 ARMA 模型的階數, 應當選取 AIC 和 BIC 值達到最小的那一組為理想階數。

plot_acf(time_series)
plot_pacf(time_series)
plt.show()

r,rac,Q = sm.tsa.acf(time_series, qstat=True) prac = pacf(time_series,method='ywmle') table_data = np.c_[range(1,len(r)), r[1:],rac,prac[1:len(rac)+1],Q] table = pd.DataFrame(table_data, columns=['lag', "AC","Q", "PAC", "Prob(>Q)"]) print(table） 

偏自相關圖.png

自相關圖.png

ARIMA.png

根據上面的幾個圖，我們可以先取 p=1, q=2。進行模型估計,結果見下圖。

p,d,q = (1,1,2)
arma_mod = ARMA(time_series,(p,d,q)).fit(disp=-1,method='mle') summary = (arma_mod.summary2(alpha=.05, float_format="%.8f")) print(summary) 

模型結果.png

這里的 p和q 參數可以調整，然后找出最佳的（AIC最小，BIC最小），經過比較， p=0,q=1 為理想階數。
這里有一個自動取 p和q 的函數，如果要自動定階的話，可以采用

(p, q) =(sm.tsa.arma_order_select_ic(dta,max_ar=3,max_ma=3,ic='aic')['aic_min_order']) #這里需要設定自動取階的 p和q 的最大值，即函數里面的max_ar,和max_ma。ic 參數表示選用的選取標准，這里設置的為aic,當然也可以用bic。然后函數會算出每個 p和q 組合(這里是(0,0)~(3,3)的AIC的值，取其中最小的,這里的結果是(p=0,q=1)。 

殘差和白噪聲檢驗

個人感覺這個就是對模型 ARIMA(0,1,1) 的殘差序列 arma_mod.resid 進行 ADF 檢驗。

arma_mod = ARMA(time_series,(0,1,1)).fit(disp=-1,method='mle') resid = arma_mod.resid t=sm.tsa.stattools.adfuller(resid) output=pd.DataFrame(index=['Test Statistic Value', "p-value", "Lags Used", "Number of Observations Used","Critical Value(1%)","Critical Value(5%)","Critical Value(10%)"],columns=['value']) output['value']['Test Statistic Value'] = t[0] output['value']['p-value'] = t[1] output['value']['Lags Used'] = t[2] output['value']['Number of Observations Used'] = t[3] output['value']['Critical Value(1%)'] = t[4]['1%'] output['value']['Critical Value(5%)'] = t[4]['5%'] output['value']['Critical Value(10%)'] = t[4]['10%'] print(output) #結果如下 Test Statistic Value -3.114 p-value 0.025534 Lags Used 1 Number of Observations Used 30 Critical Value(1%) -3.66992 Critical Value(5%) -2.96407 Critical Value(10%) -2.62117 

當然這里也可以畫出 acf 圖和 pacf 圖。

模型預測

arma_model = sm.tsa.ARMA(time_series,(0,1)).fit(disp=-1,maxiter=100) predict_data = arma_model.predict(start=str(1979), end=str(2010+3), dynamic = False) 

預測結果.png

預測結果還原

對預測出來的數據，進行逆差分操作（由原始數據取對數后的數據加上預測出來的數據），然后再取指數即可還原。

預測結果還原.png

年份	2011年	2012年	2013年
實際值	19632.26	22250.45	24668.49
預測值	19314.03	22415.10	26014.08

上圖最后3個為預測值，然后查詢2011年到2013年湖北GDP的實際值，可以進行對照

年份	2011年	2012年	2013年
實際值	19632.26	22250.45	24668.49
預測值	19314.03	22415.10	26014.08

原始數據與預測結果.png

小結

從預測對結果看，2011年到2013年的預測結果和實際的差別不大。這個模型在短期預測結果比較好。模型處理主要還是應用了Python 第三方庫 statsmodels 中的模型算法，其中還有很多細節，可以查閱相關文檔，這里只是簡單的應用了一下，由於代碼都是一小段一小段寫的，很亂，只提供了一些片段供參考。

參考資料

python時間序列分析
時間序列實戰（一）
用ARIMA模型做需求預測
python 時間序列分析之ARIMA
ARIMA模型文檔

作者：熙淺
鏈接：https://www.jianshu.com/p/4130bac8ebec
來源：簡書
簡書著作權歸作者所有，任何形式的轉載都請聯系作者獲得授權並注明出處。