R語(yǔ)言時(shí)間序列ARIMA模型方法

2019-11-08 03:17:24

字體：大中小

供稿：網(wǎng)友

原理什么的百度一搜一堆，看不明白，先學(xué)會(huì)用這個(gè)工具吧！　　ARIMA：全稱為自回歸積分滑動(dòng)平均模型(Autoregressive Integrated Moving Average Model,簡(jiǎn)記ARIMA)，是由博克思(Box)和詹金斯(Jenkins)于70年代初提出一著名時(shí)間序列預(yù)測(cè)方法，所以又稱為box-jenkins模型、博克思-詹金斯法。其中ARIMA（p，d，q）稱為差分自回歸移動(dòng)平均模型，AR是自回歸， p為自回歸項(xiàng)； MA為移動(dòng)平均，q為移動(dòng)平均項(xiàng)數(shù)，d為時(shí)間序列成為平穩(wěn)時(shí)所做的差分次數(shù)。所謂ARIMA模型，是指將非平穩(wěn)時(shí)間序列轉(zhuǎn)化為平穩(wěn)時(shí)間序列，然后將因變量?jī)H對(duì)它的滯后值以及隨機(jī)誤差項(xiàng)的現(xiàn)值和滯后值進(jìn)行回歸所建立的模型。ARIMA模型根據(jù)原序列是否平穩(wěn)以及回歸中所含部分的不同，包括移動(dòng)平均過(guò)程（MA）、自回歸過(guò)程（AR）、自回歸移動(dòng)平均過(guò)程（ARMA）以及ARIMA過(guò)程。(這是我從百度百科挖過(guò)來(lái)的-_-) 步驟總結(jié)如下：

#加載時(shí)間序列程序包library(tseries)library(forecast)#使用該包自帶的程序，是指航空乘客的分布air <- AirPassengers#作這個(gè)時(shí)間序列的圖,通過(guò)圖作一個(gè)直觀判斷plot(air)

air

#也可以直接使用tsdisplay來(lái)觀察，它包含了時(shí)序圖，以及acf、pacf兩個(gè)相關(guān)圖tsdisplay(air)

air2

#可以拆掉最后一年來(lái)做樣本的訓(xùn)練集，再將最后一年做樣本的測(cè)試集sair<-ts(as.vector(air[1:132]),frequency=12,start=c(1949,1))#同樣可以看一下拆掉之后的訓(xùn)練集圖形tsdisplay(sair)

sair

#sair的明顯存在一個(gè)向上的趨勢(shì)，用差分方法干掉，先看滯后1次的圖形s1<-diff(sair,1)#從圖上看，基本圍繞了0在振動(dòng)，基本平穩(wěn)，進(jìn)一步使用adf檢驗(yàn)，看一下是否存在單位根（驗(yàn)證平穩(wěn)性，若存在則不平穩(wěn)）adf.test(s1)#單位根檢驗(yàn)通過(guò)（p<0.05,顯著拒絕存在單位根），再看一下差分之后的圖形。tsdisplay(s1)

#圖形顯示acf圖存在拖尾，q=0，pacf圖超出虛線較多，但從整體上看，從16階之后截?cái)啵厦娴木€可以大概看到（1，8，10，12，16）這五個(gè)數(shù)超出虛線甚多，所以，可以每個(gè)都測(cè)試一下。形成(1,1,0)(8,1,0),10,1,0),(12,1,0),(16,1,0),然后判斷一下各自的AIC值，取最小值即可。#且acf圖存在一個(gè)12階的季節(jié)性影響因素，然后通過(guò)作滯后12階的圖看一下是否消除tsdisplay(diff(s1,12))

s12

#圖形可以看出，已經(jīng)消除了季節(jié)性的影響，所以ARIMA季節(jié)參數(shù)中周期應(yīng)該為12，由于季節(jié)性圖的acf圖和pacf上看到最后一個(gè)超出虛線的階數(shù)均為1，季度參數(shù)可以這么設(shè)置(1,1,1)[12]#擬合模型，看一下哪個(gè)模型最好，注意，擬合模型使用的data是原始數(shù)據(jù)，并不是差分之后的數(shù)據(jù)！arima(sair,order=c(1,1,0),seasonal=list(order=c(1,1,1),period=12))#AIC=899.95arima(sair,order=c(8,1,0),seasonal=list(order=c(1,1,1),period=12))#AIC=907.17arima(sair,order=c(10,1,0),seasonal=list(order=c(1,1,1),period=12))#AIC=909arima(sair,order=c(12,1,0),seasonal=list(order=c(1,1,1),period=12))#AIC=905.38arima(sair,order=c(16,1,0),seasonal=list(order=c(1,1,1),period=12))#NaNs#從上面的擬合可以看出，選擇(1,1,0)(1,1,1)[12]的模型AIC值最小，也可以使用auto.arima函數(shù)來(lái)自動(dòng)確定這些參數(shù)：auto.arima(sair)#Series: sair #ARIMA(0,1,1)(0,1,0)[12]#Coefficients:# ma1# -0.2263#s.e. 0.0900#sigma^2 estimated as 110.5: log likelihood=-448.34 AIC=900.69 AICc=900.79 BIC=906.24#auto.arima給出的建議是(0,1,1)(0,1,0)[12]，其AIC是900.69，相差不大，可以使用兩個(gè)方案都看一下預(yù)測(cè)的結(jié)果。#先進(jìn)行擬合fit1<-arima(sair,order=c(1,1,0),seasonal=list(order=c(1,1,1),period=12))fit2<-arima(sair,order=c(0,1,1),seasonal=list(order=c(0,1,0),period=12))#然后使用tsdiag看一下各自的結(jié)果，圖中表明殘差標(biāo)準(zhǔn)差基本都在[-1,1]之間，殘差的自回歸都為0（兩虛線內(nèi)），Ljung-Box檢驗(yàn)的p值都在0.05之上，結(jié)果不錯(cuò)。tsdiag(fit1)

fit1

tsdiag(fit2)

fit2

#預(yù)測(cè)f.p1<-forecast(fit1,h=12,level=c(99.5))f.p2<-forecast(fit2,h=12,level=c(99.5))#先看一下fit1【即(1,1,0)(1,1,1)[12]】的效果plot(f.p1,ylim=c(100,700))lines(f.p1$fitted,col="green")lines(air,col="red")

fit1

#再看一下fit2【即(0,1,1)(0,1,0)[12]】的效果plot(f.p2,ylim=c(100,700))lines(f.p2$fitted,col="green")lines(air,col="red")

fit2

#從上面兩個(gè)圖（紅色為air數(shù)據(jù)，綠色為利用sair計(jì)算出來(lái)的原始數(shù)據(jù)擬合值，藍(lán)色段為預(yù)測(cè)值），可以看出，幾乎沒(méi)有區(qū)別，而且預(yù)測(cè)效果也非常好。

時(shí)間序列預(yù)測(cè)，ARIMA的模型小結(jié)到此為止，到現(xiàn)在為止，關(guān)于AR和MA的系數(shù)還有季節(jié)參數(shù)里面的系數(shù)選擇都沒(méi)有一點(diǎn)頭緒，只能用這種窮舉法來(lái)一個(gè)個(gè)找了，費(fèi)時(shí)費(fèi)力，不過(guò)auto.arima貌似很強(qiáng)大，就是給出的系數(shù)和網(wǎng)上看來(lái)的經(jīng)驗(yàn)得出的系數(shù)有出入，明明是acf拖尾，q應(yīng)該為0，即AR模型，偏偏就使用了MA，不懂，也沒(méi)時(shí)間仔細(xì)研究它的理論，就這樣吧。

上一篇：tree樹(shù)型結(jié)構(gòu)統(tǒng)計(jì)每個(gè)層級(jí)的總數(shù)

下一篇：【C語(yǔ)言簡(jiǎn)單說(shuō)】三：整數(shù)變量擴(kuò)展和輸出擴(kuò)展(3)