時間序列分析作業(yè)

1、時間序列分析與應(yīng)用課程作業(yè)地震數(shù)據(jù)（COP.BHZ-24）時間序列分析一前言本次作業(yè)選取了第24號文件，共1440個數(shù)據(jù)。截取前1200個數(shù)據(jù)進行理分析，然后建立模型。之后再對數(shù)據(jù)進行預(yù)測，然后對1200之后的30個數(shù)據(jù)進行更新，將更新結(jié)果與原觀測值進行比對分析，最后得出結(jié)論。二數(shù)據(jù)處理1. 數(shù)據(jù)讀取與畫圖首先將文件“COP.BHZ.txt”保存到E盤根目錄下，以便于讀取。用scan()函數(shù)將數(shù)據(jù)讀入，并保存到sugar2文件中。如圖1所示。圖1 數(shù)據(jù)讀取然后，畫出該時間序列圖。橫軸表示時間，單位是*10ms，縱軸表示高程，單位是um。代碼及圖示如圖2、圖3所示。圖2 時序圖代碼圖3 前120

2、0個數(shù)據(jù)散點圖2. 平穩(wěn)性檢驗從圖中看出，該組數(shù)據(jù)隨時間變化基本平穩(wěn)，僅有小幅波動。最高點與最低點相差也僅在250um之內(nèi)。通過adf.test()函數(shù)可以驗證該假設(shè)，可以看出該序列是平穩(wěn)的（stationary）。如圖4所示。然后用求平均函數(shù)mean()求出這1200個數(shù)據(jù)的平均值a，可以從圖5看到結(jié)果。圖4 平穩(wěn)性檢驗結(jié)果圖5 求平均值然后，將原始數(shù)據(jù)減去平均值，得到一組零均值的新數(shù)據(jù)，命名為sugar3。3. 數(shù)據(jù)建模分析接下來繪制震前數(shù)據(jù)的自相關(guān)函數(shù)和偏自相關(guān)函數(shù)圖像，初步判斷其大概符合什么模型。圖6為畫出圖像的代碼，新序列sugar3的ACF、PACF圖像如下所示。圖6 ACF、PA

3、CF、EACF圖像代碼 圖7 ACF圖圖8 PACF圖從ACF、PACF圖可以看出，序列一階之后相關(guān)性較強，雖然在第19階滯后處有超限的情況，但從總體來看，兩個圖都是拖尾的情況。因此要借助于EACF圖來做進一步判斷。擴展自相關(guān)函數(shù)EACF圖如下。圖9 EACF圖3 模型識別由EACF圖可以看出此時間序列符合ARMA(0,1)或ARMA(2,2)，根據(jù)以上信息尚不能明確判斷出具體的模型，要建立確定的模型，就需要排除上述模型中的一種，用模型診斷的方法可以實現(xiàn)。模型診斷，或模型評價，涉及檢驗?zāi)Ｐ偷臄M合優(yōu)度，并且如果擬合程度很差，要給出適當(dāng)?shù)恼{(diào)整建議。模型診斷的方法有兩種：分析擬合模型的殘差和分析過度

4、參數(shù)化的模型。下面先使用殘差法。3.1 ARMA(0,1)模型診斷圖10 ARMA(0,1)模型3.1.1 殘差法圖11 畫殘差圖的代碼圖12 ARMA(0,1)模型殘差圖3.1.2 分位數(shù)圖法圖13 分位數(shù)法代碼圖14 ARMA(0,1)模型的分位數(shù)-分位數(shù)圖3.2 ARMA(2,2)模型診斷圖15 ARMA(2,2)模型3.2.1 殘差法圖16 畫殘差圖代碼圖17 ARMA(2,2)模型殘差圖3.2.2 分位數(shù)圖法圖18 分位數(shù)法代碼圖19 ARMA(0,1)模型的分位數(shù)-分位數(shù)圖從兩種待定模型的殘差圖和分位數(shù)-分位數(shù)圖看出：兩種模型都較好符合，難以比較優(yōu)劣，因此較難取舍。下面使用第二種模

5、型診斷的方法，過度擬合法。3.3 利用過度擬合法進行模型診斷3.3.1 ARMA(0,1)的過度擬合圖20 ARMA(0,1)情況過度擬合根據(jù)過度擬合時應(yīng)遵守的原則：第一，在擬合時不能同時增加AR和MA部分的階數(shù)；第二應(yīng)按殘差分析建議的方向來擴展模型。本文中，擬合了MA(1)模型，殘差在2階滯后處仍存在明顯相關(guān)性，則嘗試MA(2)、MA(3)模型，而不是ARMA(1,1)模型。由此得到以上三種模型的參數(shù)估計。比較前兩種模型，可以看出，arima(0,0,2)中的新增的參數(shù)為0.0437，不顯著地不為零。而且兩種模型的共同的參數(shù)相差很小。再結(jié)合第三種模型，同樣可以看出新增參數(shù)為-0.0096，不

6、顯著地不為零，且此三種模型的共同參數(shù)相差無幾，因此可以判斷出ARMA(0,1)模型是較合適的。 另外，結(jié)合以上三種模型的對數(shù)似然值和AIC值，可以看出，ARMA(0,2)模型AIC值較小，且對數(shù)似然值較大，因此選取ARMA(0,2)模型。下面來診斷ARMA(2,2)模型是否合適。3.3.2 ARMA(2,2)的過度擬合圖21 ARMA(2,2)情況過度擬合將ARIMA(2,0,2)模型同下面兩種模型比較可以看出：雖然額外的參數(shù)不顯著地不為零，但是它們公同的參數(shù)發(fā)生了顯著的變化，因此可以得出結(jié)論：通過過度擬合的診斷，ARMA(2,2)模型不合適。通過上面的敘述，得出了二階滑動平均MA(2)模型，

7、即ARIMA(0,0,2)模型，接下來，對該模型的參數(shù)進行估計，使用的估計方法是條件平方和估計法、極大似然估計法，通過對不同的方法得出的參數(shù)估值進行比較，得到最優(yōu)估值。圖22 參數(shù)估計圖可以看出，兩種方法得出的參數(shù)估值、方差、對數(shù)似然值等相差無幾，因此選用哪種方法所得的參數(shù)估值都可以，在這里選用對數(shù)似然值較大的CSS方法。可得模型表達(dá)式：Yt + 0.0004 = et 0.0982et-1 0.0436et-2 1.87834 模型預(yù)測時間序列建模的主要目標(biāo)之一，是預(yù)測該序列未來的取值。圖23 時間序列預(yù)測圖從圖中可知，預(yù)測結(jié)果并不好，好像哪里出了問題。（可能哪里代碼有誤，暫時還沒有發(fā)現(xiàn)。）

8、下面利用新建模型更新原觀測數(shù)值30步，即1201至1230。圖24 預(yù)測更新代碼圖25 預(yù)測更新結(jié)果從上圖更新的30步觀測值可以看出，雖然擬合模型的趨勢（黃色）大致與實際觀測結(jié)果的走向（黑色）大致相同，但是有些點上還是有不小的擬合誤差，例如2號、10號、27號點等。究其原因，可能是因為所選本組數(shù)據(jù)相對平緩數(shù)據(jù)變化相對穩(wěn)定，因此較難擬合出趨勢較明顯的模型。5 心得體會通過一學(xué)期的課程學(xué)習(xí)，我感覺收貨頗豐。首先，在學(xué)習(xí)上我不僅學(xué)會了如何對一組時間序列進行基本的模型識別、參數(shù)估計、模型診斷預(yù)測更新等理論方法，還學(xué)會了如何用R語言去具體實現(xiàn)。同時也鍛煉了自己的自學(xué)及鉆研能力。其次，在課堂上，張老師、李

9、老師除了交給我們課本上的理論知識外，還教導(dǎo)我們?yōu)槿颂幨赖牡览?，特別是張老師，從老師的身上學(xué)到了如何去尊重別人等基本的卻是容易被人忽視的基本禮節(jié)。學(xué)術(shù)固然重要，禮節(jié)卻更可貴。參考書目：1 Jonathan D.Cryer, Kung-Sik Chan.時間序列分析與應(yīng)用M.機械工業(yè)出版社.2 王福林,王吉權(quán),吳昌友,吳秋峰.  實數(shù)遺傳算法的改進研究J. 生物數(shù)學(xué)學(xué)報. 2006(01) 3 孫曉云,高鑫,王鵬.  新型并行遺傳算法及其在參數(shù)估計中的應(yīng)用J. 計算機工程與應(yīng)用. 2005(19) 4 田小梅,龔靜.  實數(shù)編碼遺傳算法的評述J. 湖南環(huán)境生物職業(yè)技術(shù)學(xué)院學(xué)報. 2005(01) 5 陳曉梅,楊成祥.  遺傳進化算法在時間序列建模中的應(yīng)用J. 計算機工程與應(yīng)用. 2005(05) 6 來鵬,陳平.  基于遺傳算法的時序混合模型的參數(shù)估計J. 山西師范大學(xué)學(xué)報(自然科學(xué)版). 2004(03) 7 孫靖,程大章.  基于季節(jié)性時間序列模型的空調(diào)負(fù)荷預(yù)測J. 電工技術(shù)學(xué)報. 2004(03) 8 安德洪,柳湘月,劉嘉焜,許樹荊.  基