一、準(zhǔn)備工作

1、數(shù)據(jù)準(zhǔn)備

所使用的數(shù)據(jù)是TSA包中的co2數(shù)據(jù)，如果沒(méi)有這個(gè)包的話，可以先裝一下

install.packages("TSA")	# 安裝包 TSA

會(huì)有讓你選鏡像的過(guò)程，隨便選就行了。下載好之后，導(dǎo)入并查看數(shù)據(jù)

library(TSA)
data(co2)
win.graph(width = 4.875,height = 3,pointsize = 8)
plot(co2,ylab='CO2')        #繪制原始數(shù)據(jù)

可以看到，原始數(shù)據(jù)明顯有一個(gè)向上的趨勢(shì)和一個(gè)周期趨勢(shì)。

2、基本概念

赤池信息準(zhǔn)則（Akaike’s(1973) Information Criterion, AIC）是建立在熵的概念基礎(chǔ)上，可以權(quán)衡所估計(jì)模型的復(fù)雜度和此模型擬合數(shù)據(jù)的優(yōu)良性。
$AIC=?2log(極大似然估計(jì)值)+2k$
其中，如果模型包含截距或常數(shù)項(xiàng)，那么k=p+q+1；否則k=p+q。AIC越小越好。

Ljung-Box檢驗(yàn)即LB檢驗(yàn)、隨機(jī)性檢驗(yàn)，用來(lái)檢驗(yàn)m階滯后范圍內(nèi)序列的自相關(guān)性是否顯著，或序列是否為白噪聲（或者統(tǒng)計(jì)量服從自由度為m的卡方分布）。若是白噪聲數(shù)據(jù)，則該數(shù)據(jù)沒(méi)有價(jià)值提取，即不用繼續(xù)分析了。

二、數(shù)據(jù)處理

拿到一個(gè)序列之后，首先判斷它是不是平穩(wěn)時(shí)間序列，如果是就進(jìn)行模式識(shí)別；如果不是就扣除趨勢(shì)項(xiàng)將其變成一個(gè)平穩(wěn)時(shí)間序列。接著做模式識(shí)別、參數(shù)估計(jì)、模型診斷和預(yù)測(cè)。

ps: 這是從老師課件上找的流程圖，個(gè)人感覺(jué)模式識(shí)別部分，不應(yīng)包含參數(shù)d，因?yàn)楹?code>d的一般是ARIMA(p,d,q)模型，它是非平穩(wěn)模型，而上一步已將非平穩(wěn)時(shí)間序列轉(zhuǎn)成平穩(wěn)時(shí)間序列了，d應(yīng)該是在上一步確認(rèn)的。也有可能是，扣除趨勢(shì)項(xiàng)和模式識(shí)別的界限根本不可能分那么細(xì)，但是又要用流程圖表示出來(lái)，所以才這么寫的。

1、模式識(shí)別

一般來(lái)講，模式識(shí)別就是判別出ARIMA(p,d,q)中的各階數(shù)p,d,q。模式識(shí)別常用的方法有：acf, pacf, eacf

首先來(lái)看它的自相關(guān)函數(shù)

acf(as.vector(co2),lag.max = 36)                #自相關(guān)函數(shù)

季節(jié)自相關(guān)關(guān)系十分顯著：在滯后12，24，36，……上表現(xiàn)出很強(qiáng)的相關(guān)性。

plot(diff(co2),ylab='1st Diff. of CO2',xlab='Year') #一次差分，消除整體上升趨勢(shì)

可以看到，經(jīng)過(guò)一次差分后，序列中的整體上升趨勢(shì)已經(jīng)消除。再來(lái)看其樣本自相關(guān)函數(shù)

acf(as.vector(diff(co2)),lag.max = 36)          #一次差分后的自相關(guān)函數(shù)

一次差分后，序列中仍存在強(qiáng)烈的季節(jié)性；應(yīng)用季節(jié)差分法應(yīng)該可以得到更為簡(jiǎn)約的模型。

plot(diff(diff(co2),lag=12),xlab='Year',ylab='1st & seasonal Diff.')

繪制其自相關(guān)函數(shù)

acf(as.vector(diff(diff(co2),lag=12)),lag.max=36,ci.type='ma')

可以看到，經(jīng)過(guò)一次差分和季節(jié)差分后的時(shí)間序列已經(jīng)消除了季節(jié)性的大部分影響。根據(jù)樣本自相關(guān)函數(shù)可以看到，除了在滯后1和12上具有自相關(guān)性外，經(jīng)一次和季節(jié)差分后的序列幾乎不再具有自相關(guān)性，所建模型只需要在滯后1和12上具有自相關(guān)性即可。

綜上，考慮構(gòu)建乘法季節(jié) $ARIMA(0,1,1)\times (0,1,1{)}_{12}$模型。

2、參數(shù)估計(jì)

模型建立后，需要估計(jì)模型的參數(shù)。乘法季節(jié)ARIMA模型只是一般ARIMA模型的特例。

m1.co2=arima(co2,order=c(0,1,1),seasonal=list(order=c(0,1,1),period=12))
print(m1.co2)
-------------------------------------------
Coefficients:
          ma1     sma1
      -0.5792  -0.8206
s.e.   0.0791   0.1137

sigma^2 estimated as 0.5446:  log likelihood = -139.54,  aic = 283.08

上面第一行代碼便得到了參數(shù)的極大似然估計(jì)值，參數(shù)估值的標(biāo)準(zhǔn)差為0.5446，對(duì)數(shù)似然值為-139.54，AIC=283.08。模型的參數(shù)估值均為高度顯著，進(jìn)而將對(duì)該模型加以檢驗(yàn)。

ps:為什么能根據(jù)這些指標(biāo)值說(shuō)明參數(shù)估值高度顯著，標(biāo)準(zhǔn)是多少？

3、診斷性檢驗(yàn)

1 殘差序列

首先觀察殘差的時(shí)間序列圖

plot(window(rstandard(m1.co2),start=c(1995,2)),ylab='Standardized Resi.',type='o');
abline(h=0)

除了序列中間存在某些異常行為外，殘差圖中并沒(méi)有表明模型有任何主要的不規(guī)則性。然后對(duì)殘差的樣本自相關(guān)函數(shù)進(jìn)一步觀察

acf(as.vector(window(rstandard(m1.co2),start=c(1995,2))),lag.max=36)

統(tǒng)計(jì)上顯著的相關(guān)系數(shù)位于滯后22，其值僅為-0.194，相關(guān)性非常小，且滯后22上的依賴關(guān)系難以給出合理的解釋。除了滯后22處的邊緣顯著以外，該模型似乎已捕捉到了序列中依賴關(guān)系的本質(zhì)。

注：在acf前打個(gè)print即可輸出滯后各階的自相關(guān)函數(shù)的值。

2 Ljung-Box 檢驗(yàn)

下面進(jìn)行Ljung-Box檢驗(yàn)：

win.graph(width=3,height =3,pointsize = 8)
hist(window(rstandard(m1.co2),start=c(1995,2)),xlab='Standardized Resi.',ylab='Frequency')

直方圖的形狀像鐘形，但并不標(biāo)準(zhǔn)。對(duì)模型進(jìn)行Ljung-Box檢驗(yàn)，給出自由度為22的x=25.59, p=0.27，表明該模型已捕獲時(shí)間序列中的依賴關(guān)系。

why?

R 語(yǔ)言進(jìn)行 Ljung-Box 檢驗(yàn)的函數(shù)如下：

Box.test(x, lag = 1, type = c("Box-Pierce", "Ljung-Box"), fitdf = 0)

• x: 一個(gè)時(shí)間序列，殘差檢驗(yàn)時(shí)，一般是殘差
• lag: 基于自相關(guān)因子得出的lag值
• type: Ljung-Box 檢驗(yàn)就設(shè)置為 Ljung-Box
• fitdf: 如果x是一系列殘差，則需要減去自由度。

調(diào)用函數(shù)后，我們關(guān)心的就是p值，如果p > 0.05，則說(shuō)明是白噪聲序列，是純隨機(jī)性序列。否則數(shù)據(jù)不是白噪聲，具有研究?jī)r(jià)值。

示例如下：

x <- rnorm (100)
Box.test(x, lag = 5)
Box.test(x, lag = 10, type = "Ljung")
a=Box.test(resid(m1.xpole),type="Ljung",lag=20,fitdf=11)

接著繪制分位數(shù)-分位數(shù)圖（qq圖）

win.graph(width=5,height =5,pointsize = 8)
qqnorm(window(rstandard(m1.co2),start=c(1995,2)))
abline(c(0,0),c(1,1),col='red')

QQ 圖的上尾部，再次出現(xiàn)了一個(gè)異常值。但是，Shapiro-Wilk正態(tài)性檢驗(yàn)給出的統(tǒng)計(jì)量W=0.982，進(jìn)而得到p=0.11，且在任何顯著水平上正態(tài)性都未被拒絕。

作為對(duì)模型的進(jìn)一步檢驗(yàn)，考慮用$ARIMA(0,1,2)\times (0,1,1{)}_{12}$模型進(jìn)行過(guò)度擬合。

m2.co2=arima(co2,order=c(0,1,2),seasonal=list(order=c(0,1,1),period=12))
print(m1.co2)
print(m2.co2)
--------------------------------
arima(x = co2, order = c(0, 1, 1), seasonal = list(order = c(0, 1, 1), period = 12))
Coefficients:
          ma1     sma1
      -0.5792  -0.8206
s.e.   0.0791   0.1137
sigma^2 estimated as 0.5446:  log likelihood = -139.54,  aic = 283.08
--------------------------------
arima(x = co2, order = c(0, 1, 2), seasonal = list(order = c(0, 1, 1), period = 12))

Coefficients:
          ma1      ma2     sma1
      -0.5714  -0.0165  -0.8274
s.e.   0.0897   0.0948   0.1224

sigma^2 estimated as 0.5427:  log likelihood = -139.52,  aic = 285.05

可以看到，${\theta }_1$和$\theta$的估計(jì)變化很小（考慮標(biāo)準(zhǔn)差的大小時(shí)）。新參數(shù)${\theta }_2$的估值在統(tǒng)計(jì)上與零無(wú)異。AIC 已經(jīng)增加的情況下，sigma^2和對(duì)數(shù)似然值均無(wú)顯著變化。所以使用$ARIMA(0,1,2)\times (0,1,1{)}_{12}$模型是過(guò)度擬合，根據(jù)“奧卡姆剃刀原理”，如無(wú)必要，勿增實(shí)體。所以使用$ARIMA(0,1,1)\times (0,1,1{)}_{12}$模型即可。

4、預(yù)測(cè)

前置時(shí)間設(shè)為2年，進(jìn)行預(yù)測(cè)2年并繪制預(yù)測(cè)值與預(yù)測(cè)極限。

win.graph(width = 4.875,height = 3,pointsize = 8)
plot(m1.co2,n1=c(2003,1),nahead=24,xlab='Year',type='o',ylab='CO2 Levels')

前置時(shí)間設(shè)為1年，進(jìn)行預(yù)測(cè)4年并繪制預(yù)測(cè)值與預(yù)測(cè)極限。

win.graph(width = 4.875,height = 3,pointsize = 8)
plot(m1.co2,n1=c(2004,1),n.ahead=48,xlab='Year',type='b',ylab='CO2 Levels')

ps: 上面參數(shù)的含義，n1=c(2004,1)代表從2004年1月開始（實(shí)際數(shù)據(jù)到2005年12月結(jié)束），n.ahead=48代表預(yù)測(cè)48個(gè)值（一年12個(gè)值，所以是4年）文章來(lái)源地址http://www.zghlxwxcb.cn/news/detail-469818.html

到了這里，關(guān)于R 語(yǔ)言做時(shí)間序列分析的實(shí)例（模式識(shí)別、擬合、檢驗(yàn)、預(yù)測(cè)）的文章就介紹完了。如果您還想了解更多內(nèi)容，請(qǐng)?jiān)谟疑辖撬阉鱐OY模板網(wǎng)以前的文章或繼續(xù)瀏覽下面的相關(guān)文章，希望大家以后多多支持TOY模板網(wǎng)！