R語言進(jìn)行ARIMA分析重點(diǎn)講義

1、R學(xué)習(xí)日記時(shí)間序列分析之ARIMA模型預(yù)測今天學(xué)習(xí)ARIMA預(yù)測時(shí)間序列。 指數(shù)平滑法對(duì)于預(yù)測來說是非常有幫助的，而且它對(duì)時(shí)間序列上面連續(xù)的值之間相關(guān)性沒有要求。但是，如果你想使用指數(shù)平滑法計(jì)算出預(yù)測區(qū)間， 那么預(yù)測誤差必須是不相關(guān)的， 而且必須是服從零均值、 方差不變的正態(tài)分布。即使指數(shù)平滑法對(duì)時(shí)間序列連續(xù)數(shù)值之間相關(guān)性沒有要求，在某種情況下， 我們可以通過考慮數(shù)據(jù)之間的相關(guān)性來創(chuàng)建更好的預(yù)測模型。自回歸移動(dòng)平均模型（ ARIMA） 包含一個(gè)確定（explicit） 的統(tǒng)計(jì)模型用于處理時(shí)間序列的不規(guī)則部分，它也允許不規(guī)則部分可以自相關(guān)。首先，先確定數(shù)據(jù)的差分。ARIMA 模型為平穩(wěn)

2、時(shí)間序列定義的。 因此， 如果你從一個(gè)非平穩(wěn)的時(shí)間序列開始， 首先你就需要做時(shí)間序列差分直到你得到一個(gè)平穩(wěn)時(shí)間序列。如果你必須對(duì)時(shí)間序列做 d 階差分才能得到一個(gè)平穩(wěn)序列，那么你就使用ARIMA(p,d,q)模型，其中 d 是差分的階數(shù)。 我們以每年女人裙子邊緣的直徑做成的時(shí)間序列數(shù)據(jù)為例。從 1866 年到 1911 年在平均值上是不平穩(wěn)的。 隨著時(shí)間增加， 數(shù)值變化很大。  > skirts <- scan("Read 46 items> skirtsts<- ts(skirts,start = c(1866) > pl

3、ot.ts(skirtsts)我們可以通過鍵入下面的代碼來得到時(shí)間序列（數(shù)據(jù)存于“skirtsts”） 的一階差分， 并畫出差分序列的圖:> skirtstsdiff<-diff(skirtsts,differences=1) > plot.ts(skirtstsdiff)從一階差分的圖中可以看出，數(shù)據(jù)仍是不平穩(wěn)的。我們繼續(xù)差分。> skirtstsdiff2<-diff(skirtsts,differences=2)> plot.ts(skirtstsdiff2)二次差分（上面）后的時(shí)間序列在均值和方差上確實(shí)看起來像是平穩(wěn)的， 隨著時(shí)間推移， 時(shí)

4、間序列的水平和方差大致保持不變。因此， 看起來我們需要對(duì)裙子直徑進(jìn)行兩次差分以得到平穩(wěn)序列。第二步，找到合適的ARIMA模型 如果你的時(shí)間序列是平穩(wěn)的，或者你通過做 n 次差分轉(zhuǎn)化為一個(gè)平穩(wěn)時(shí)間序列， 接下來就是要選擇合適的 ARIMA模型，這意味著需要尋找 ARIMA(p,d,q)中合適的 p 值和 q 值。為了得到這些，通常需要檢查平穩(wěn)時(shí)間序列的（自）相關(guān)圖和偏相關(guān)圖。  我們使用 R 中的“acf()”和“pacf” 函數(shù)來分別（ 自） 相關(guān)圖和偏相關(guān)圖?！癮cf()”和“pacf 設(shè)定“plot=FALSE” 來得到自相關(guān)和偏相關(guān)的真實(shí)值。 &g

5、t; acf(skirtstsdiff2,lag.max=20) > acf(skirtstsdiff2,lag.max=20,plot=FALSE)Autocorrelations of series skirtstsdiff2, by lag     0      1      2      3      4  

6、60;   5      6      7      8      9     10  1.000 -0.303  0.096  0.009  0.102 -0.453  0.173 -0.025 -0.039  0.073 -0.094  

7、;   11     12     13     14     15     16     17     18     19     20  0.133 -0.089 -0.027 -0.1

8、02  0.207 -0.260  0.114  0.101  0.011 -0.090 自相關(guān)圖顯示滯后1階自相關(guān)值基本沒有超過邊界值，雖然5階自相關(guān)值超出邊界，那么很可能屬于偶然出現(xiàn)的，而自相關(guān)值在其他上都沒有超出顯著邊界， 而且我們可以期望 1 到 20 之間的會(huì)偶爾超出 95%的置信邊界。  > pacf(skirtstsdiff2,lag.max=20)> pacf(skirtstsdiff2,lag.max=20,plot=FALSE)Partial autocorrelations of series ski

9、rtstsdiff2, by lag     1      2      3      4      5      6      7      8     &#

10、160;9     10     11 -0.303  0.005  0.043  0.128 -0.439 -0.110  0.073  0.028  0.128 -0.355  0.095     12     13     14     15  

11、0;  16     17     18     19     20  0.052 -0.094 -0.103 -0.034 -0.021 -0.002  0.074  0.020 -0.034 偏自相關(guān)值選5階。故我們的ARMIA模型為armia（1,2,5）> skirtsarima<-arima(skirtsts,order=c(1,2,5)>

12、 skirtsarimaSSeries: skirtsts ARIMA(1,2,5)                    Coefficients:          ar1     ma1     ma2  

13、;   ma3     ma4      ma5      -0.4345  0.2762  0.1033  0.1472  0.0267  -0.8384s.e.   0.1837  0.2171  0.2198  0.2716  0.1904   0.2888sigma2 estima

14、ted as 206.1:  log likelihood=-183.8AIC=381.6   AICc=384.71   BIC=394.09預(yù)測后5年裙子的邊緣直徑>  skirtsarimaforecast<-forecast.Arima(skirtsarima,h=5,level=c(99.5)> skirtsarimaforecast     Point Forecast  Lo 99.5  Hi 99.51912   

15、;    548.5762 507.1167 590.03571913       545.1793 459.3292 631.02951914       540.9354 396.3768 685.49401915       531.8838 316.2785 747.48921916       529.129

16、6 233.2625 824.9968> plot.forecast(skirtsarimaforecast$residuals)   #謝謝憶水如煙的指正第三步，檢驗(yàn)在指數(shù)平滑模型下， 觀察 ARIMA 模型的預(yù)測誤差是否是平均值為 0 且方差為常數(shù)的正態(tài)分布（服從零均值、方差不變的正態(tài)分布） 是個(gè)好主意，同時(shí)也要觀察連續(xù)預(yù)測誤差是否（自）相關(guān)。  > acf(skirtsarimaforecast$residuals,lag.max=20)> Box.test(skirtsarimaforecast$residuals, lag=2

17、0, type="Ljung-Box")        Box-Ljung testdata:  skirtsarimaforecast$residuals X-squared = 8.5974, df = 20, p-value = 0.9871既然相 關(guān)圖顯示出在滯后1 - 20階（ l a g s 1 - 20 ）中樣本自相關(guān)值都沒有超出顯著（置信）邊 界，而且Ljung-Box檢驗(yàn)的p值為0.99，所以我們推斷在滯后1-20階（lags1-20）中沒明顯證據(jù)說明預(yù)測 誤差是非零自相關(guān)的。 為了調(diào)查預(yù)測誤差是否是平均值為零且方差為常數(shù)的正態(tài)分布（服從零均值、方差不變的正態(tài)分布），我們可以做預(yù)測誤差的時(shí)間曲線圖和直方圖（具有正態(tài)分布曲線）：> plot.ts(skirtsarimaforecast$residuals)> plotForecastErrors(skirtsarimaforecast$residuals)