沈永梅，王 瓊，何哲飛

(常州大學(xué) 數(shù)理學(xué)院，江蘇 常州 213164)

0 引言

馬爾可夫鏈狀態(tài)預(yù)測方法已經(jīng)在水資源科學(xué)中得到了廣泛的應(yīng)用。張宸、林啟太[1]研究了馬爾可夫鏈理論在礦區(qū)降水災(zāi)害預(yù)測中的應(yīng)用；鄭文瑞等[2]研究了馬爾可夫鏈理論在水污染狀態(tài)風(fēng)險(xiǎn)評價(jià)中的應(yīng)用；張文堅(jiān)[3]用馬氏鏈預(yù)測理論研究了城鎮(zhèn)洪澇災(zāi)害的分析；宋印勝[4]和王渺林[5]研究了馬爾可夫鏈理論在水位預(yù)測中的應(yīng)用；馮耀龍，韓文秀[6]研究了馬爾可夫鏈在河流豐枯狀況預(yù)測中的應(yīng)用。近年來，一些學(xué)者開始研究基于馬氏鏈狀態(tài)預(yù)測方法的點(diǎn)值預(yù)測方法[7]。本文擬提出一種簡單的基于馬氏鏈狀態(tài)預(yù)測方法的點(diǎn)值預(yù)測方法，并將其與普遍接受的時(shí)間序列分析預(yù)測方法進(jìn)行基于統(tǒng)計(jì)試驗(yàn)的比較分析?；隈R氏鏈狀態(tài)預(yù)測方法的點(diǎn)值預(yù)測方法的預(yù)測結(jié)果是通過MATLAB軟件實(shí)現(xiàn)的，而時(shí)間序列分析預(yù)測方法的結(jié)果是通過DecisionTime軟件實(shí)現(xiàn)的。

由于我國有關(guān)水文計(jì)算的規(guī)范規(guī)定，我國主要河流的徑流量和降水量均可假定為服從P-III型分布或三參數(shù)對數(shù)正態(tài)分布。因此本文的主要目的是利用隨機(jī)模擬技術(shù)，生成服從P-III型分布和三參數(shù)對數(shù)正態(tài)分布的相依偽隨機(jī)數(shù)序列，然后再將基于馬氏鏈狀態(tài)預(yù)測方法的點(diǎn)值預(yù)測方法與時(shí)間序列分析預(yù)測方法進(jìn)行比較分析。研究的主要問題如下：（1)馬爾可夫鏈點(diǎn)值預(yù)測方法與時(shí)間序列分析預(yù)測方法的精度哪個(gè)高；（2)對于馬爾可夫鏈的三種點(diǎn)值預(yù)測方法，哪種預(yù)測方法較優(yōu)；（3)馬爾可夫鏈點(diǎn)值預(yù)測方法對P-III型分布和三參數(shù)對數(shù)正態(tài)分布水文分布具有適應(yīng)性；（4）序列的變差系數(shù)Cv和偏態(tài)系數(shù)Cs對預(yù)測精度的影響。

1 馬爾可夫鏈點(diǎn)值預(yù)測方法的思想

馬爾可夫鏈的狀態(tài)預(yù)測方法可見文[6]，這里主要說明如何根據(jù)區(qū)間預(yù)測結(jié)果來確定點(diǎn)預(yù)測值。根據(jù)經(jīng)驗(yàn)可知，下一個(gè)時(shí)段的指標(biāo)值與上一個(gè)時(shí)段的指標(biāo)值有很大的關(guān)系，故本文提出了基于本時(shí)段 預(yù)測中值與上一時(shí)段預(yù)測中值的預(yù)測方法，這里的預(yù)測中值是指預(yù)測狀態(tài)所對應(yīng)的狀態(tài)區(qū)間的中點(diǎn)。假設(shè)狀態(tài)空間E={1,2,…,m}，設(shè)各狀態(tài)區(qū)間的中點(diǎn)值分別為M1,M2,…,Mm，通過馬爾可夫鏈的區(qū)間預(yù)測方法得出本時(shí)段的預(yù)測狀態(tài)為s(t)，上一時(shí)段的狀態(tài)為s(t-1)，則本時(shí)段的點(diǎn)預(yù)測值為βMs(t-1)+(1-β)Ms(t)，其中β為上一時(shí)段預(yù)測中值在整個(gè)點(diǎn)預(yù)測值中所占的權(quán)重。β的選取一般可采用如下方法：假設(shè)通過馬氏鏈狀態(tài)預(yù)測方法所得的上一時(shí)段的預(yù)測狀態(tài)為s(t-1)，而轉(zhuǎn)移到該狀態(tài)的概率為Ps(t-1)；本時(shí)段的預(yù)測狀態(tài)為s(t)，轉(zhuǎn)移到該狀態(tài)的概率為Ps(t)，則

基于絕對分布的馬爾可夫鏈預(yù)測(ADMCP)方法對應(yīng)的點(diǎn)值預(yù)測方法記為ADMCPP，疊加馬爾可夫鏈預(yù)測(SPMCP)方法對應(yīng)的點(diǎn)值預(yù)測方法記為SPMCPP，加權(quán)馬爾可夫鏈預(yù)測(WMCP)方法對應(yīng)的點(diǎn)值預(yù)測方法記為WMCPP。

2 最優(yōu)時(shí)間序列分析預(yù)測方法的確定

時(shí)間序列分析預(yù)測方法比較多，其中ARIMA模型是其中較為常見也非常重要的一類數(shù)據(jù)處理和預(yù)測的方法。從已有數(shù)據(jù)序列選擇一個(gè)好的ARIMA模型需要很多步驟，包括模型的初步識別、定階、參數(shù)估計(jì)、模型的檢驗(yàn)等，這些都給模型的選擇帶來了很大的困難。統(tǒng)計(jì)試驗(yàn)中所需處理的數(shù)據(jù)序列很多，不可能進(jìn)行一一分析，所以需要有一種新的方法來解決這個(gè)問題。

本文選用了DecisionTime軟件進(jìn)行數(shù)據(jù)處理和預(yù)測分析，由于數(shù)據(jù)序列很多，所以選擇了專家預(yù)測模式(Expert Modeler)，DecisionTime軟件會(huì)根據(jù)數(shù)據(jù)自身的特征選擇最適合的模型來擬合和預(yù)測數(shù)據(jù)。所選擇的模型除了前面所提到的模型還有如下7種：Simple exponential smoothing、Holt’s exponential smoothing、Brown’s exponential smoothing、Damped exponential smoothing、Seasonal exponential smoothing、Additive Winters’exponential smoothing、Multiplicative Winters’exponential smoothing。其實(shí)這些模型都是包含于ARIMA模型中的，也就是說它們都是ARIMA模型的特例，之所以用他們來擬合和預(yù)測數(shù)據(jù)就是在擬合效果差不多的情況下盡量的節(jié)省數(shù)據(jù)的運(yùn)算和處理的時(shí)間，這些模型雖然對長期的預(yù)測效果不如ARIMA模型，但是他們對短期預(yù)測還是很有效的，而且它們對有線性趨勢和季節(jié)性的序列能處理得更好。

這里需要說明的是，在專家預(yù)測模式下DecisionTime軟件給出了最適合你的數(shù)據(jù)的模型。如果在預(yù)測向?qū)е袥]有你選定的模型，專家預(yù)測模式會(huì)選擇最好的指數(shù)平滑（修勻）模型和最好的單變量的ARIMA模型：首先，專家預(yù)測模式會(huì)選擇用最小的標(biāo)準(zhǔn)BIC值確定最適合的指數(shù)平滑模型；然后，專家預(yù)測模式會(huì)判斷這組數(shù)據(jù)是否需要進(jìn)行變換或者是通過差分使其平穩(wěn)化，再通過ACF和PACF圖像來確定一個(gè)初始模型，這個(gè)模型是調(diào)整過的或者多次調(diào)整過的，通過t值和Ljung-Box統(tǒng)計(jì)量，以及殘差A(yù)CF和殘差PACF圖像；最后專家預(yù)測模式會(huì)在指數(shù)平滑模型和ARIMA模型中選擇標(biāo)準(zhǔn)BIC值更小的模型作為預(yù)測用模型。本文稱通過這種方法得到的最優(yōu)預(yù)測方法為時(shí)間序列預(yù)測方法，記其為“TSAP”法。

3 試驗(yàn)設(shè)計(jì)與試驗(yàn)結(jié)果的綜合分析

3.1 試驗(yàn)設(shè)計(jì)

3.1.1 偽隨機(jī)數(shù)序列的長度和分布參數(shù)的設(shè)計(jì)

由于天然河流的實(shí)測水文資料比較短，一般少于60年，因此只有在小樣本下進(jìn)行預(yù)測方法的比較才有意義。出于這種考慮，本文對于生成的序列，僅考慮長度為100的服從P-III型分布和三參數(shù)對數(shù)正態(tài)分布的相依偽隨機(jī)數(shù)序列，分別用馬爾可夫鏈點(diǎn)值預(yù)測方法和較優(yōu)的時(shí)間序列預(yù)測方法對相同的數(shù)據(jù)進(jìn)行預(yù)測，并比較兩類預(yù)測方法的效果。我們分別對兩種分布采用了5組總體參數(shù)（見表1）共生成了10組長度為80000的數(shù)據(jù)序列，數(shù)據(jù)的生成方法見文[8]。

經(jīng)驗(yàn)表明，大多數(shù)實(shí)際徑流量和降水量資料的參數(shù)都在表1的參數(shù)值的范圍內(nèi)，由于一般的水文序列均為弱相關(guān)序列，所以為了方便起見序列的一階相關(guān)系數(shù)統(tǒng)ρx(1)一規(guī)定為0.2。

3.1.2 計(jì)算方案設(shè)計(jì)及預(yù)測偏差的衡量

按表1中給定的各組參數(shù)分別生成5組長度為80000的服從P-III型分布和對數(shù)正態(tài)分布的相依偽隨機(jī)數(shù)序列，依次取100個(gè)數(shù)據(jù)(分別考慮取50次、200次和800次)。對于每取得的長度為100的相依偽隨機(jī)數(shù)序列，我們可假定前C個(gè)數(shù)據(jù)為實(shí)測資料，利用這C年的資料序列，分別運(yùn)用馬爾可夫鏈點(diǎn)值預(yù)測方法和較優(yōu)的時(shí)間序列預(yù)測方法預(yù)測后100-C年的數(shù)據(jù)值，并比較兩種預(yù)測方法的優(yōu)劣。比較的指標(biāo)設(shè)定如下：可以用生成的數(shù)據(jù)值與預(yù)測值求得絕對預(yù)測誤差向量，這是一個(gè)維數(shù)為100-C的向量。假如重復(fù)試驗(yàn)N（N=50、200、800）次，則馬爾可夫鏈點(diǎn)值預(yù)測方法和較優(yōu)時(shí)間序列分析預(yù)測方法都可以得到N個(gè)絕對預(yù)測誤差向量，這些向量的維數(shù)均為100-C，如果把這些向量作為行向量排成矩陣，可構(gòu)成一個(gè)N×（100-C）的絕對預(yù)測誤差矩陣，記為ΕN×(100-C)，如果對其列求平均值，則可以得到多次試驗(yàn)下的絕對預(yù)測誤差的平均向量，它能反映預(yù)測方法的準(zhǔn)確性，從而反映預(yù)測方法的優(yōu)劣；若對其列求方差，則可以得到在多次試驗(yàn)情況下的絕對預(yù)測誤差的方差向量，其能反映預(yù)測方法的穩(wěn)健性。如果要對兩類方法進(jìn)行比較，還需找到了一個(gè)能反映向量大小的量，本文選用的這個(gè)量為向量的各個(gè)分量的平均值。

為了更清楚地說明上面一段文字，可用數(shù)學(xué)符號把相關(guān)量表示出來。假設(shè)重復(fù)試驗(yàn)的次數(shù)為N，其中第k(k=1,2,…,N)次所取的100個(gè)數(shù)據(jù)值為，假定前C個(gè)實(shí)測值記為，后100-C個(gè)數(shù)據(jù)值記為，用兩類預(yù)測方法得到的預(yù)測結(jié)果統(tǒng)一記為，則絕對預(yù)測誤差矩陣為ΕN×(100-C)=(eki)，其中eki=|x(k,i)[-x[(k,p;i)|,i=1,2,…,100-C。

絕對預(yù)測誤差均值向量mvape（mean vector of absolute predicting error）的第i個(gè)分量定義為：

用其平均值來衡量預(yù)測方法優(yōu)劣，簡記為mmvape(mean of mvape)

絕對預(yù)測誤差方差向量vvape（variance vector of absolute predicting error）的第i個(gè)分量定義為：

其中

用其均值來衡量預(yù)測方法的穩(wěn)健性，簡記為mvvape(mean of vvape)

表1 總體參數(shù)值表

預(yù)測誤差分析表中每一種預(yù)測方法下面有兩列，分別記錄了每組數(shù)據(jù)用該方法時(shí)的mmvape和mvvape。每張表格的倒數(shù)第二行是關(guān)于五組不同參數(shù)的數(shù)據(jù)的mmvape和mvvape平均值，最后一行是以較優(yōu)的時(shí)間序列預(yù)測方法所得的預(yù)測值為基準(zhǔn)，基于馬爾可夫鏈狀態(tài)預(yù)測方法的點(diǎn)值預(yù)測方法的相對誤差下降率。

3.2 預(yù)測結(jié)果匯編

3.2.1 對數(shù)正態(tài)分布數(shù)據(jù)的兩類預(yù)測方法的統(tǒng)計(jì)試驗(yàn)結(jié)果

選定實(shí)測資料數(shù)據(jù)長度C為60，利用這60年的資料序列，分別運(yùn)用馬氏鏈點(diǎn)值預(yù)測方法和較優(yōu)的時(shí)間序列預(yù)測方法預(yù)測后40年的數(shù)據(jù)值。馬爾可夫鏈預(yù)測方法中選取的狀態(tài)分級參數(shù)α=(α1,α2,α3,α4)如表2。

表2 不同參數(shù)值所對應(yīng)的α表

表3 對數(shù)正態(tài)分布數(shù)據(jù)的馬氏鏈點(diǎn)值預(yù)測方法及較優(yōu)時(shí)間序列預(yù)測方法的預(yù)測誤差分析N=50

N=200

N=800

表4 不同參數(shù)值所對應(yīng)的α表

表5 P-Ⅲ型分布數(shù)據(jù)的馬氏鏈點(diǎn)值預(yù)測方法及較優(yōu)時(shí)間序列預(yù)測方法的預(yù)測誤差分析N=50

基于統(tǒng)計(jì)試驗(yàn)的預(yù)測結(jié)果如表3，表3中分別列出了試驗(yàn)次數(shù)N取50、200、800次的條件下的預(yù)測結(jié)果。

3.2.2 P-Ⅲ型分布數(shù)據(jù)的兩類預(yù)測方法的統(tǒng)計(jì)試驗(yàn)結(jié)果

同樣選定實(shí)測資料數(shù)據(jù)長度C為60，利用這60年的資料序列，分別運(yùn)用馬氏鏈點(diǎn)值預(yù)測方法和較優(yōu)的時(shí)間序列預(yù)測方法預(yù)測后40年的數(shù)據(jù)值。馬爾可夫鏈預(yù)測方法中選取的狀態(tài)分級參數(shù)α=(α1,α2,α3,α4)，如表4。

基于統(tǒng)計(jì)試驗(yàn)的預(yù)測結(jié)果如表5。表5中分別列出了試驗(yàn)次數(shù)N取50、200、800次的條件下的預(yù)測結(jié)果。

從服從三參數(shù)對數(shù)正態(tài)分布和P-III型分布的數(shù)據(jù)的預(yù)測誤差分析表3和表5可以得出下面四個(gè)結(jié)論：（1)馬爾可夫鏈點(diǎn)值預(yù)測方法的精度比時(shí)間序列分析預(yù)測方法要高；（2)馬爾可夫鏈的三種點(diǎn)值預(yù)測方法中基于絕對分布的馬爾可夫鏈點(diǎn)值預(yù)測方法的精度最低，其它兩種點(diǎn)值預(yù)測方法預(yù)測精度相差不大；（3)馬爾可夫鏈點(diǎn)值預(yù)測方法對不同的水文分布（三參數(shù)對數(shù)正態(tài)分布和P-III型分布）具有適應(yīng)性；（4)序列的變差系數(shù)Cv和偏態(tài)系數(shù)Cs對兩類預(yù)測方法的精度的影響顯著。精度是隨著Cv的增大而迅速下降，隨著Cs的增大而有所提高。其他某些因素對預(yù)測精度的影響可參見文[9]。

4 結(jié)語

本文的結(jié)論僅限于服從P-III型分布或三參數(shù)對數(shù)正態(tài)分布的水文時(shí)間序列。文中所提出的基于馬氏鏈狀態(tài)預(yù)測的點(diǎn)值預(yù)測方法是一種比較簡單的預(yù)測方法，本文的主旨在于基于統(tǒng)計(jì)試驗(yàn)的方法將其與時(shí)間序列分析預(yù)測方法進(jìn)行比較，并未對點(diǎn)值預(yù)測做深入的研究。結(jié)論說明在水文時(shí)間序列預(yù)測中基于馬氏鏈狀態(tài)預(yù)測方法的點(diǎn)值預(yù)測方法更為有效，故可以更進(jìn)一步地研究如何通過區(qū)間預(yù)測結(jié)果獲得點(diǎn)預(yù)測值的方法，從而更有效地獲得點(diǎn)預(yù)測值。

[1]張宸，林啟太.模糊馬爾可夫鏈狀模型在礦區(qū)降水災(zāi)害預(yù)測中的應(yīng)用[J].國外建材科技，2004，(1).

[2]鄭文瑞，王新代，紀(jì)昆，王漢林.非確定數(shù)學(xué)方法在水污染狀態(tài)風(fēng)險(xiǎn)評價(jià)中的應(yīng)用[J].吉林大學(xué)學(xué)報(bào)，2003，(1).

[3]張文堅(jiān).90年代浙江城鎮(zhèn)洪澇災(zāi)害分析及其展望[J].科技通報(bào)，1996，(1).

[4]宋印勝.馬爾可夫鏈模型在地下水水位預(yù)測中的應(yīng)用[J].山東地質(zhì)，1998，(1).

[5]王渺林.灰色馬爾可夫模型在寸灘站年最高水位預(yù)測中的應(yīng)用[J].2004，(2).

[6]馮耀龍，韓文秀.加權(quán)馬爾可夫鏈在河流豐枯狀況預(yù)測中的應(yīng)用[J].系統(tǒng)工程理論與實(shí)踐,1999，(10).

[7]孫才志，林學(xué)鈺.降水預(yù)測的模糊權(quán)馬爾可夫模型及應(yīng)用[J].系統(tǒng)工程學(xué)報(bào),2003，(4).

[9]沈永梅.基于統(tǒng)計(jì)試驗(yàn)的馬氏鏈點(diǎn)值預(yù)測方法和時(shí)間序列分析預(yù)測方法的比較分析[D].河海大學(xué),2006.