時(shí)間序列：ARIMA模型

1、.實(shí)驗(yàn)：建立 ARIMA模型（綜合性實(shí)驗(yàn)）實(shí)驗(yàn)題目： 某城市連續(xù) 14年的月度嬰兒出生率數(shù)據(jù)如下表所示：26.66323.59826.93124.74025.80624.36424.47723.90123.17523.22721.67221.87021.43921.08923.70921.66921.75220.76123.47923.82423.10523.11021.75922.07321.93720.03523.59021.67222.22222.12323.95023.50422.23823.14221.05921.57321.54820.00022.42420.61521.76122

2、.87424.10423.74823.26222.90721.51922.02522.60420.89424.67723.67325.32023.58324.67124.45424.12224.25222.08422.99123.28723.04925.07624.03724.43024.66726.45125.61825.01425.11022.96423.98123.79822.27024.77522.64623.98824.73726.27625.81625.21025.19923.16224.70724.36422.64425.56524.06225.43124.63527.00926

3、.60626.26826.46225.24625.18024.65723.30426.98226.19927.21026.12226.70626.87826.15226.37924.71225.68824.99024.23926.72123.47524.76726.21928.36128.59927.91427.78425.69326.88126.21724.21827.91426.97528.52727.13928.98228.16928.05629.13626.29126.98726.58924.84827.54326.89628.87827.39028.06528.14129.04828

4、.48426.63427.73527.13224.92428.96326.58927.93128.00929.22928.75928.40527.94525.91226.61926.07625.28627.66025.95126.39825.56528.86530.00029.26129.01226.99227.897（ 1）選擇適當(dāng)模型擬和該序列的發(fā)展（ 2）使用擬合模型預(yù)測(cè)下一年度該城市月度嬰兒出生率實(shí)驗(yàn)內(nèi)容：給出實(shí)際問(wèn)題的非平穩(wěn)時(shí)間序列，要求學(xué)生利用R 統(tǒng)計(jì)軟件，對(duì)該序列進(jìn)行分析，通過(guò)平穩(wěn)性檢驗(yàn)、差分運(yùn)算、白噪聲檢驗(yàn)、擬合ARMA模型，建立ARIMA模型，在此基礎(chǔ)上進(jìn)行預(yù)測(cè)。實(shí)驗(yàn)要求：處

5、理數(shù)據(jù)， 掌握非平穩(wěn)時(shí)間序列的ARIMA建模方法， 并根據(jù)具體的實(shí)驗(yàn)題目要求完成實(shí)驗(yàn)報(bào)告，并及時(shí)上傳到給定的FTP 和課程網(wǎng)站。實(shí)驗(yàn)步驟：第一步：編程建立R 數(shù)據(jù)集；第二步：調(diào)用plot.ts 程序?qū)?shù)據(jù)繪制時(shí)序圖。第三步：從時(shí)序圖中利用平穩(wěn)時(shí)間序列的定義判斷是否平穩(wěn)？第四步：若不滿足平穩(wěn)性，則可利用差分運(yùn)算是否能使序列平穩(wěn)？重復(fù)第三步步驟第五步： 根據(jù) Box.test純隨機(jī)檢驗(yàn)結(jié)果， 利用 LB 統(tǒng)計(jì)量和白噪聲特性檢驗(yàn)最后處理的'.時(shí)間序列是否為純隨機(jī)序列？第六步：在序列判斷為平穩(wěn)非白噪聲序列后，求出該觀察值序列的樣本自相關(guān)系數(shù)（ACF）和樣本偏自相關(guān)系數(shù)（PACF）的值，選擇階數(shù)

6、適當(dāng)?shù)腁RIMA（ p,d,q ）模型進(jìn)行擬合，并估計(jì)模型中未知參數(shù)的值。第七步：檢驗(yàn)?zāi)Ｐ偷挠行浴Ｈ绻麛M合模型通不過(guò)檢驗(yàn)，轉(zhuǎn)向步驟6，重新選擇模型再擬合。第八步：模型優(yōu)化。如果擬合模型通過(guò)檢驗(yàn)，仍然轉(zhuǎn)向步驟6，充分考慮各種可能建立多個(gè)擬合模型，從所有通過(guò)檢驗(yàn)的擬合模型中選擇最優(yōu)模型。第九步：利用最優(yōu)擬合模型，預(yù)測(cè)下一年度該城市月度嬰兒出生率。ex5.2=ts(scan("ex5.2.txt"), frequency=4)Read 168 itemsplot.ts(ex5.2)從圖中看出序列一開(kāi)始有下降趨勢(shì)，后面有明顯上升趨勢(shì)，所以序列不平穩(wěn)。d12ex5.2 = diff

7、(ex5.2,lag=12)acf(d12ex5.2,48)plot(d12ex5.2)'.從上面的自相關(guān)圖中可以看出改做滯后12 期差分后為平穩(wěn)。Box.test(d12ex5.2, lag=17, type="Ljung-Box")Box-Ljung testdata:d12ex5.2X-squared = 147.9254, df = 17, p-value < 2.2e-16P 值小于 0.05，可以認(rèn)為是非白噪聲序列。par(mfrow=c(2,1); acf(d12ex5.2, 48); pacf(d12ex5.2, 48)'.ARIMA（

8、 0,0,3 ）、 ARIMA（ 0,0,4 ）、 ARIMA（ 1,0,3 ）、 ARIMA（ 1,0,4 ）四個(gè)模型分別進(jìn)行擬合檢驗(yàn)(rec.ols = arima(d12ex5.2,order=c(0,0,3)Call:arima(x = d12ex5.2, order = c(0, 0, 3)Coefficients:ma1ma2ma3intercept0.79490.44800.11560.2150s.e.0.08390.08320.08850.1744sigma2 estimated as 0.8621:log likelihood = -210.12,aic = 430.25re

9、c.pr = predict(rec.ols, n.ahead=5)U = rec.pr$pred + 1.96*rec.pr$seL = rec.pr$pred - 1.96*rec.pr$seminx = min(d12ex5.2,L)maxx = max(d12ex5.2,U)ts.plot(d12ex5.2, rec.pr$pred, ylim=c(minx,maxx)'.lines(rec.pr$pred, col="red", type="o")lines(U, col="blue", lty="dash

10、ed")lines(L, col="blue", lty="dashed")qqnorm(rec.ols$resid)qqline(rec.ols$resid)'.shapiro.test(rec.ols$resid)Shapiro-Wilk normality testdata:rec.ols$residW = 0.9777, p-value = 0.0125用 shapiro 檢驗(yàn)，發(fā)現(xiàn)p 值為 0.0125，在 5%的顯著性水平下顯著，所以為ARIMA（ 0,0,3 ）模型不合理。(rec.ols = arima(d12ex5.2

11、,order=c(0,0,4)Call:arima(x = d12ex5.2, order = c(0, 0, 4)Coefficients:ma1ma2ma3ma4 intercept0.83060.49430.22540.20700.2041s.e.0.09020.11580.09250.08890.1994'.sigma2 estimated as 0.828:log likelihood = -207.07,aic = 426.15rec.pr = predict(rec.ols, n.ahead=5)U = rec.pr$pred + 1.96*rec.pr$seL = re

12、c.pr$pred - 1.96*rec.pr$seminx = min(d12ex5.2,L)maxx = max(d12ex5.2,U)ts.plot(d12ex5.2, rec.pr$pred, ylim=c(minx,maxx)lines(rec.pr$pred, col="red", type="o")lines(U, col="blue", lty="dashed")lines(L, col="blue", lty="dashed")qqnorm(rec.ols$

13、resid)qqline(rec.ols$resid)'.shapiro.test(rec.ols$resid)Shapiro-Wilk normality testdata:rec.ols$residW = 0.9689, p-value = 0.001363用 shapiro 檢驗(yàn)，發(fā)現(xiàn) p 值為 0.001363，在 5%的顯著性水平下顯著，所以為 ARIMA（ 0,0,4 ）模型不合理。(rec.ols = arima(d12ex5.2,order=c(1,0,3)Call:arima(x = d12ex5.2, order = c(1, 0, 3)Coefficients:a

14、r1ma1ma2ma3intercept0.9288-0.1369-0.2156-0.15860.0240s.e.0.06690.10650.09210.08790.4984sigma2 estimated as 0.7986:log likelihood = -204.35,aic = 420.7'.rec.pr = predict(rec.ols, n.ahead=5)U = rec.pr$pred + 1.96*rec.pr$seL = rec.pr$pred - 1.96*rec.pr$seminx = min(d12ex5.2,L)maxx = max(d12ex5.2,U)

15、ts.plot(d12ex5.2, rec.pr$pred, ylim=c(minx,maxx)lines(rec.pr$pred, col="red", type="o")lines(U, col="blue", lty="dashed")lines(L, col="blue", lty="dashed")qqnorm(rec.ols$resid)qqline(rec.ols$resid)'.shapiro.test(rec.ols$resid)Shapiro-Wi

16、lk normality testdata:rec.ols$residW = 0.9783, p-value = 0.01454(rec.ols = arima(d12ex5.2,order=c(1,0,4)Call:arima(x = d12ex5.2, order = c(1, 0, 4)Coefficients:ar1ma1ma2ma3ma4 intercept0.9084-0.1288-0.2457-0.15110.13090.0493s.e.0.07250.10780.10210.07780.09600.4684'.sigma2 estimated as 0.7891:log likelihood = -203.47,aic = 420.94rec.pr = predict(rec.ols, n.ahead=5)U = rec.pr$pred + 1.96*rec.pr$seL = rec.pr$pred - 1.96*rec.pr$seminx = min(d12ex5.2,L)maxx = max(d12ex5.2,U)ts.plot(d12ex5.2, rec.pr$pred, ylim=c(minx,maxx)lines(rec.pr$pred, col="red", type="o&q