復權馬爾可夫鏈及其在怒江水沙預測中的應用

劉新有+彭海英+吳捷+謝飛帆

摘要：[HJ18mm]針對傳統馬爾可夫鏈及其改進的預測方法只能進行狀態預測的局限，根據相依隨機變量的特點，在以傳統馬爾可夫鏈預測方法求得各狀態預測概率的基礎上，進一步以狀態預測概率為權重與狀態平均值加權求和，實現了馬爾可夫鏈預測方法從狀態預測到數值預測的關鍵性改進。利用我國西南國際大河怒江干流道街壩水文站1957-2010年徑流和1964-2010年懸移質輸沙序列為分析期，2011-2015年徑流和懸移質輸沙為驗證期，對所建立的復權馬爾可夫鏈預測方法步驟進行驗證表明，復權馬爾可夫鏈預測方法具有較高的數值預測精度，能夠滿足隨機時間序列短期數值預測的需要。

關鍵詞：復權馬爾可夫鏈；數值預測；徑流；懸移質輸沙；怒江

中圖分類號：P333文獻標識碼：A文章編號：

16721683（2017）06002607

Abstract：In view of the limitations of traditional Markov chain and its improved prediction methods which can only predict the state，in this paper we realized a critical improvement of the Markov chain forecasting method to being able to conduct numerical predictionWe did so by using weighted summation of the average value of each state multiplied by the corresponding predicted probability，on the basis of obtaining the predicted probability of each state with the traditional Markov chain forecasting method according to the characteristics of dependent stochastic variablesThe data of this study were collected from Daojieba hydrological station on the Nujiang river，which is a famous international river in southwest ChinaWe used the runoff series from 1957 to 2010 and the suspended sediment series from 1964 to 2010 for analysis，and used the runoff and suspended sediment series from 2011 to 2015 for validationResults showed that the reweighted Markov chain forecasting had a high accuracy in numerical prediction and could meet the demand of shortterm numerical prediction in stochastic time series

Key words：reweighted Markov chain；numerical prediction；runoff；suspended sediment；Nujiang river

馬爾可夫鏈是俄羅斯數學家馬爾可夫1906-1912年間提出的一種隨機事件預測的重要方法，在教育、經濟、生物、農業、災害、水文氣象、環境預測等眾多領域得到了廣泛應用。尤其在水文氣象預測中，馬爾可夫鏈預測方法應用非常廣泛，并在應用過程中不斷得以改進，加權馬爾可夫鏈[112]、灰色馬爾可夫鏈[1317]、疊加馬爾可夫鏈[18]、時間序列馬爾可夫模型[19]、基于多重轉移概率的馬爾可夫模型[20]均取得了較好的預測精度。夏樂天等[2122]系統研究了各種馬爾可夫鏈預測方法在水文預測中的應用，并對比了三種常用馬爾可夫鏈預測方法的優劣，認為加權馬爾可夫鏈預測方法精度最高。這些研究為馬爾可夫鏈預測方法的應用和發展起到了積極作用，但這些改進方法仍然沒有超出對隨機事件狀態預測的范疇。因此，如何根據馬爾可夫鏈預測狀態概率分布得到預測值仍然有待解決[12]。本文在加權馬爾可夫鏈預測、基于絕對分布的馬爾可夫鏈預測和疊加馬爾可夫鏈預測方法的基礎上，進一步以狀態預測概率為權重，結合狀態平均值進行加權求和，實現了馬爾可夫鏈預測方法從狀態預測到數值預測的關鍵性改進，并通過怒江水沙預測實例對復權馬爾可夫鏈預測方法的數值預測精度進行驗證。

1復權馬爾可夫鏈預測方法

馬爾可夫鏈通過統計隨機事件過去一定時期內的狀態轉移概率來預測將來狀態變化的概率，其中時間參數集T={0，1，2，…}及狀態參數集E={0，1，2，…}稱為馬爾可夫鏈。在實際應用中，一般采用齊次馬爾可夫鏈，即對任意參數u，k∈T，有

Pij（u；k）∈E（1）

式中：Pij（u；k）表示隨機事件u時段所處的狀態i，經過k步狀態轉移后變為狀態j的概率。

傳統齊次馬爾可夫鏈的狀態轉移步長一般取1，即利用初始分布推算未來狀態的絕對分布，沒有考慮各種步長馬爾可夫鏈的絕對分布在預測中所起的作用。為彌補這一缺陷，一些學者將各種步長馬爾可夫鏈求得的狀態絕對分布疊加起來進行狀態預測，但在疊加過程中沒有考慮各種步長在權重上的差異。因此，利用各種步長自相關性的強弱確定不同步長權重的加權馬爾可夫鏈進行狀態預測更符合實際[9]。但由于加權馬爾可夫鏈得到的預測結果仍然是狀態，在實際應用中受到一定的限制。復權馬爾可夫鏈在之前的研究基礎上，進一步以各狀態的預測概率為權重，結合其對應狀態均值進行加權求和，從而實現從狀態預測到數值預測的跨越。endprint

2復權馬爾可夫鏈預測方法步驟

復權馬爾可夫鏈以馬爾科夫鏈求得的各狀態的預測概率為基礎，因此步驟（1）至（9）與加權馬爾可夫鏈基本一致。但為方便對復權馬爾可夫鏈預測的理解，本研究以加權馬爾科夫鏈為基礎，完整介紹復權馬爾可夫鏈預測方法步驟。

（1）初步判斷對象序列是否是隨機變量。若受大型水利工程等人為控制則不適用于馬爾可夫鏈，反之則可能適用于馬爾可夫鏈，最終確定是否適用于馬爾可夫鏈有待馬氏性檢驗結果。

（2）建立序列狀態分級標準，確定資料序列的對應狀態。常用的狀態分級方法有聚類分析法、樣本均值標準差分級法、頻率曲線法等。水文分析中常用PIII型頻率曲線法來確定各年份的豐枯狀態，且為使樣本序列具有代表性，一般要求樣本序列不應少于30年。根據狀態分級標準，即可確定資料序列所對應的狀態。

（3）用fij表示指標值序列x1，x2，…xn中從狀態i經過一步或多步轉移到達狀態j的頻數，i，j∈E。對資料序列所對應的狀態進行統計計算，得到各狀態的轉移規律，進而建立各階（步長）的狀態轉移頻數矩陣。

（4）將狀態轉移頻數矩陣（fij）i，j∈E的第i行第j列元素fij除以各行的總和所得的值稱為“轉移概率”，記為Pij，i，j∈E，即

Pij=[SX（]fij[]∑[DD（]m[]j=1[DD）]fij（2）

式中：m為指標值序列包含的可能的狀態。

（5）對隨機變量進行馬氏檢驗。

將轉移概率矩陣（pij）的第j列之和除以各行各列的總和所得的值稱為“邊際概率”，記為Pj，即

Pj=[SX（]∑[DD（]m[]i=1[DD）]fij[]∑[DD（]m[]i=1[DD）]∑[DD（]m[]j=1[DD）]fij（3）

則當序列長度n充分大時，統計量

X2=2∑[DD（]m[]i=1[DD）]∑[DD（]m[]j=1[DD）]fij[JB（|]lg[SX（]Pij[]Pj[JB）|]（4）

給定顯著性水平α，查表可得分位點X2α（（m-1）2）的值，計算后得統計量X2的值，若X2>X2α·（（m-1）2），則可以認為{Xi}符合馬氏性，否則可以認為該序列不可作為馬氏鏈來處理。

[JP3]（6）計算各階（步長）自相關系數。計算公式如下：

rk=∑[DD（]n-k[]l=1[DD）]（xl-[AKx-]）（xl+k-[AKx-]）∑[DD（]n[]l=1[DD）]（xl-[AKx-]）2（5）

式中：rk為第k步長自相關系數；xl為序列的第l個值；[AKx-]為序列均值；n為序列長度。

（7）各步長自相關系數規范化。計算公式如下：

wk=|rk|∑[DD（]c[]k=1[DD）]|rk|（6）

式中：wk為規范化后的各步長自相關系數，即各步長的馬爾可夫鏈權重；c為按預測需要的最大步長。

（8）以各種步長為初始狀態，結合其對應的轉移概率矩陣，預測其狀態概率Pki。

（9）將同一狀態的各預測概率加權求和，得到該狀態的預測概率，即

Pi=∑[DD（]m[]k=1[DD）]wkPki（7）

（10）以各狀態的預測概率Pi為權重，與其對應狀態的均值[AKx-]i加權求和，得到預測值d，即

d=∑Pi[AKx-]i（8）

將預測值加入原序列，再重復以上步驟，即可進行下一步的數值預測。

基于絕對分布的復權馬爾可夫鏈預測和疊加復權馬爾可夫鏈預測方法與基于加權馬爾可夫鏈的復權馬爾可夫鏈預測方法相似，即在各自的狀態預測概率基礎上[21]加上步驟（10）求得預測值。

3怒江水沙預測應用實例

怒江薩爾溫江是全球最典型的南北向發育國際大河，其上游中國境內稱為怒江。怒江流域屬峽谷地形，南北跨度大，獨特地理環境和氣候條件使其成為全球生物多樣性最突出的地區之一，怒江也蘊藏了極為豐富的水能資源。但由于多種原因，怒江干流水電開發一直未能實施，其水文過程至今沒有受到水利工程等人類活動的控制。本文以怒江干流道街壩水文站1957-2015年徑流和1964-2015年懸移質輸沙序列為數據基礎，并將1957-2010年徑流和1964-2010年懸移質輸沙序列作為預測方法分析期，將2011-2015年徑流和懸移質輸沙作為預測方法的驗證期，以說明復權馬爾可夫鏈預測方法的具體應用并檢驗預測精度。道街壩水文站控制流域面積1102[KG-7]萬[KG-9]km2，占中國境內怒江干流流域面積的883%，該站徑流和懸移質輸沙變化基本能代表怒江干流徑流和懸移質輸沙變化特征。

以怒江道街壩站1957-2010年54年徑流量和1964-2010年47年懸疑質輸沙為例，預測2011年徑流量和懸疑質輸沙量，以基于加權馬爾可夫鏈為基礎的復權馬爾可夫預測為例，詳細介紹其計算過程。

（1）初步判斷道街壩站年徑流和年懸移質輸沙序列是否是隨機變量。怒江干流水電梯級開發尚未實施，徑流和懸移質輸沙沒有受到人為控制。同時，怒江干流流域云南段涉及5個縣區，但2014年末總人口僅15965萬人，社會經濟發展落后，加之山高水低，耕地少且分散，區域內農業以自然耕種為主，產流產沙條件基本保持天然狀態。因此，可初步判斷怒江徑流和懸移質輸沙序列屬隨機變量。

（2）建立道街壩站年徑流和年懸移質輸沙序列分級標準。徑流和年懸移質輸沙序列長度超過30年，樣本具有代表性，宜采用PIII型分布頻率曲線法來確定其所處狀態。分別以保證率0～125%、>125%～375%、>375%～625%、>625%～875%、>875%～100%將年徑流和年懸移質輸沙分為豐、偏豐、平、偏少、少5級，對應狀態E={1，2，3，4，5}。年徑流和年懸移質輸沙PIII型分布各保證率對應的數值見表1。endprint

（3）按照分級標準，確定年徑流和年懸移質輸沙序列對應的狀態（表2）。

（4）據表2進行統計分析，得到1至5階（步長）狀態轉移頻數矩陣（表3、表4）。

注：矩陣a，b，c，d，e分別為步長1，2，3，4，5的馬爾科夫轉移頻數矩陣，下同。

（5）對1至5階（步長）狀態轉移頻數矩陣進行統計分析，得到各個步長的馬爾科夫鏈轉移概率矩陣（表5、表6）。

（6）結合步長為1的轉移概率矩陣和式（3）、式（4），求得怒江道街壩站64年徑流量和47年懸疑質輸沙量序列對應的邊際概率和統計量x2，計算得到x2值分別為3029和3471，大于α=005顯著性

水平下分位點X2α（16）的值26296，因此該徑流量和懸疑質輸沙序列滿足馬氏性。

（7）按照式（5）、式（6）分別計算各步長自相關系數和馬爾可夫鏈權重，結果如表7所示。

（8）以各種滯時為初始狀態，結合相應的轉移概率矩陣預測其狀態概率。依據2010、2009、2008、2007、2006年的年徑流和年懸移質輸沙量及其相應的狀態轉移概率矩陣，結合式（7）將同一狀態的各預

測概率加權求和即可對2011年的年徑流和年懸移質輸沙狀態概率進行預測（表8、表9）。

（9）將各狀態的預測概率作為權重，與其對應狀態的均值依據式（8）進行加權求和，即可得到2011年徑流和年懸移質輸沙的預測值（表10）。

（10）由表10可知，2011年徑流預測值1 766 m3s與實測值對比，相對預測誤差為208%；2011年懸移質輸沙預測值1 234 kgs與實測值對比，相對預測誤差為508%。

將預測值加入原序列，重復以上步驟，即可得到2012-2015年徑流和懸疑質輸沙量的預測值。基于絕對分布的馬爾可夫鏈預測和疊加馬爾可夫鏈預測方法的復權馬爾可夫預測方法與基于加權馬爾科夫鏈的復權馬爾科夫預測方法相似，即先求得各狀態的預測概率，再以預測概率為權重，結合數據序列中各對應狀態的均值加權求和，即求得數值預測結果。

表11是在加權馬爾可夫鏈預測、基于絕對分布的馬爾可夫鏈預測和疊加馬爾可夫鏈預測方法的基礎上，以各狀態預測概率為權重，結合狀態平均值進行加權求和的復權馬爾科夫預測方法的數值預測結果。由表11可知，2011-2015年預測值的徑流和懸疑質輸沙量序列的馬爾可夫檢驗統計量X2均大于26926，說明預測所用的時間序列在005顯著性水平下均滿足馬氏性。預測值與實測值對比表明，徑流數值預測精度總體高于懸移質輸沙數值預測精度，這可能是由于相對于徑流量而言，懸移質輸沙受人類活動影響更大，導致怒江干流懸移質輸沙狀態之間的數值跨度較大，1964-2015年期間怒江輸沙極值比達719。2011-2014年怒江干流徑流數值預測精度相對較高，而2015年徑流數值預測精度較低；2011-2013年怒江干流輸沙數值預測精度相對較高，而2014-2015年輸沙數值預測精度較低，這可能與復權馬爾科夫鏈更適合短期預測有關，隨著預測時間的延長，預測誤差可能被逐步放大。本研究短期數值預測結果精度與馬占青等[23]基于馬爾可夫鏈預測模型的杭州市降水量數值預測結果精度相當，高于馬建琴等[24]的改進型灰色馬爾可夫鏈模型對三門峽入庫年徑流的預測精度，為較少受人類活動控制的河流的徑流和輸沙的數值預測提供了一條值得探索的途徑。

4結論

已有的馬爾可夫鏈預測方法多限于進行狀態預測，而本文建立的復權馬爾可夫鏈預測方法能夠進行數值預測，實現了對馬爾可夫鏈預測方法的關鍵性改進，不僅提高了預測精度，也擴展了該方法的應用范圍。不受人為控制的隨機性序列和足夠的序列長度，是適用于馬爾馬爾可夫鏈的前提條件。復權馬爾可夫鏈在馬爾可夫鏈前期研究的基礎上，進一步以各狀態的預測概率為權重，結合狀態平均值進行二次加權求和，從而實現數值預測。與其他馬爾可夫鏈改進方法相比，復權馬爾可夫鏈能更充分地挖掘隨機序列的信息。怒江干流水沙預測實例表明，所建立的復權馬爾可夫鏈預測方法思路清晰、物理概念明確、計算簡便，為提高隨機變量的數值預測精度提供了一種可行的途徑。

參考文獻（References）：

[1]賀娟，王曉松，王彩云加權馬爾可夫鏈模型在密云水庫入庫流量中的應用[J]南水北調與水利科技，2015，13（4）：618621（HE J，WANG X S，WANG C YApplication of the weighted Markov chain model in the inflow prediction of the Miyun Reservoir[J]SouthtoNorth Water Transfers and Water Science & Technology，2015，13（4）：618621（in Chinese）） DOI：1013476jcnkinsbdqk201504003

[2]王濤，錢會，李培月加權馬爾可夫鏈在銀川地區降雨量預測中的應用[J]南水北調與水利科技，2010，8（1）：7881（WANG T，QIAN H，LI P YPrediction of precipitation based on the weighted Markov chain in Yinchuan area[J]SouthtoNorth Water Transfers and Water Science & Technology，2010，8（1）：7881（in Chinese）） DOI：103969jissn16721683201001021

[3]王亞雄，黃淑嫻，劉祖發，等變化環境下北江下游年徑流量的加權馬爾可夫鏈預測[J]生態環境學報，2011，20（4）：754760（WANG Y X，HUANG S X，LIU Z F，et alForecast of yearly river runoff in lower reaches of Beijiang River by Weighted MarkovChain Method in changing environments[J]Ecology and Environmental Sciences，2011，20（4）：754760（in Chinese）） DOI：103969jissn16745906201104030endprint