劉章君，郭生練，何紹坤，巴歡歡，尹家波

（1.武漢大學(xué) 水資源與水電工程科學(xué)國(guó)家重點(diǎn)實(shí)驗(yàn)室，水資源安全保障湖北省協(xié)同創(chuàng)新中心，湖北 武漢　430072；2.江西省水利科學(xué)研究院，江西 南昌　330029）

1　研究背景

水文不確定性處理器（Hydrologic Uncertainty Processor，HUP）作為貝葉斯預(yù)報(bào)系統(tǒng)（Bayesian Forecasting System，BFS）［1］的一個(gè)重要組件，用來(lái)量化除定量降水預(yù)報(bào)不確定性以外的其他所有不確定性，在實(shí)際中應(yīng)用廣泛［2-5］。然而，HUP屬于單變量結(jié)構(gòu)類型，只能獨(dú)立地給出各預(yù)見期實(shí)際流量的貝葉斯后驗(yàn)概率密度，而沒(méi)有考慮它們之間的內(nèi)在相關(guān)性。Krzysztofowicz和Maranzano［6］以單步馬爾科夫轉(zhuǎn)移密度函數(shù)為工具，提出了基于貝葉斯轉(zhuǎn)移預(yù)報(bào)（Bayesian Transition Forecast，BTF）方法的多變量水文不確定性處理器（Multivariate HUP，MHUP），可以在給定確定性預(yù)報(bào)過(guò)程的條件下提供實(shí)際流量過(guò)程的后驗(yàn)聯(lián)合概率密度函數(shù)，考慮了不同預(yù)見期流量之間的相關(guān)性。在此基礎(chǔ)上，針對(duì)傳統(tǒng)的貝葉斯概率洪水預(yù)報(bào)方法均只能提供某一時(shí)刻流量的概率預(yù)報(bào)，無(wú)法滿足決策者對(duì)某一時(shí)段內(nèi)最大流量的預(yù)報(bào)不確定性信息的實(shí)際需求［7-8］，Krzysztofowicz定義了極值概率預(yù)報(bào)的概念，并利用MHUP輸出產(chǎn)品給出了極值變量的超越概率分布函數(shù)，提出了貝葉斯極值預(yù)報(bào)（Bayesian Extremum Forecast，BEF）方法［9］。

BTF方法是MHUP的基本組成構(gòu)件，通常采用基于正態(tài)分位數(shù)轉(zhuǎn)換和線性-正態(tài)假設(shè)的亞高斯模型進(jìn)行計(jì)算，但這種正態(tài)分位數(shù)轉(zhuǎn)換在外推極端事件時(shí)效果不穩(wěn)健，且逆轉(zhuǎn)換時(shí)也可能使結(jié)果偏離最優(yōu)值，影響了該法的適用性［10-11］。鑒于Copula函數(shù)能較好地模擬水文過(guò)程的非線性和非正態(tài)特征，在構(gòu)建多變量聯(lián)合分布和條件分布中具有巨大優(yōu)勢(shì)［12-15］，劉章君等［11］提出了基于Copula函數(shù)的HUP（CHUP），不需要進(jìn)行線性-正態(tài)假設(shè)，適用范圍更廣，應(yīng)用更加靈活。事實(shí)上，BTF方法中的先驗(yàn)轉(zhuǎn)移密度和似然函數(shù)本質(zhì)上也都是條件概率密度函數(shù)，但目前沒(méi)有文獻(xiàn)將Copula函數(shù)引入BTF方法、MHUP及BEF方法的研究中。

因此，本文將采用Copula函數(shù)描述BTF方法中的先驗(yàn)轉(zhuǎn)移密度和似然函數(shù)，推導(dǎo)后驗(yàn)轉(zhuǎn)移密度的解析表達(dá)式，提出基于Copula函數(shù)的貝葉斯轉(zhuǎn)移預(yù)報(bào)（CBTF）方法，得到基于Copula函數(shù)的多變量水文不確定性處理器（CMHUP），進(jìn)而利用CMHUP輸出的各預(yù)見期實(shí)際流量后驗(yàn)聯(lián)合概率密度函數(shù)，導(dǎo)出預(yù)見期時(shí)段內(nèi)洪峰流量的超過(guò)概率分布函數(shù)，發(fā)展基于Copula函數(shù)的貝葉斯極值預(yù)報(bào)（CBEF）方法。以三峽水庫(kù)為例進(jìn)行應(yīng)用，基于所提方法實(shí)現(xiàn)三峽水庫(kù)入庫(kù)流量的轉(zhuǎn)移概率預(yù)報(bào)和極值概率預(yù)報(bào)。

2　多變量水文不確定性處理器

假設(shè)H0表示預(yù)報(bào)時(shí)刻的當(dāng)前實(shí)測(cè)流量，分別表示待預(yù)報(bào)的實(shí)際流量和確定性預(yù)報(bào)流量，K為預(yù)見期長(zhǎng)度；h0、hk和sk分別為隨機(jī)變量H0、Hk和Sk的實(shí)現(xiàn)值。多變量水文不確定性處理器的目的是在給定當(dāng)前實(shí)測(cè)流量H0=h0、確定性預(yù)報(bào)過(guò)程隨機(jī)向量的實(shí)現(xiàn)值的條件下，推求得到實(shí)際流量過(guò)程隨機(jī)向量的后驗(yàn)聯(lián)合概率密度函數(shù)。利用概率乘法定理，后驗(yàn)聯(lián)合概率密度函數(shù)可以因式分解為［6］：

相應(yīng)地，式（1）也可以簡(jiǎn)化為：

3　基于Copula函數(shù)的貝葉斯轉(zhuǎn)移預(yù)報(bào)方法

3.1Copula函數(shù)Copula函數(shù)可以將多個(gè)隨機(jī)變量的邊緣分布連接起來(lái)構(gòu)造聯(lián)合分布。令G(x1，x2，…，xn)為一個(gè)n-維分布函數(shù)，其邊緣分布分別為則存在一個(gè)n-Copula函數(shù)C，使得對(duì)任意x∈Rn（x為n維向量，Rn為n維實(shí)數(shù)空間）［13］：

圖1　CBTF方法概念示意圖

式中：θ為Copula函數(shù)的相關(guān)性參數(shù)。

Copula函數(shù)可以將多個(gè)變量的邊緣分布和它們之間的相關(guān)性結(jié)構(gòu)分開來(lái)研究，且對(duì)邊緣分布類型沒(méi)有任何限制。采用Copula函數(shù)進(jìn)行多變量建模時(shí)，通常包括兩個(gè)步驟：①確定各個(gè)變量的邊緣分布；②選擇合適的Copula函數(shù)，描述變量之間的相關(guān)性結(jié)構(gòu)。

本文選用Archimedean Copula函數(shù)族構(gòu)造聯(lián)合分布，采用K-S檢驗(yàn)法對(duì)Copula函數(shù)進(jìn)行擬合檢驗(yàn)的基礎(chǔ)上，基于RMSE準(zhǔn)則評(píng)價(jià)聯(lián)合分布的理論頻率與經(jīng)驗(yàn)頻率擬合情況，選擇RMSE值最小的Copula函數(shù)作為最優(yōu)的Copula函數(shù)，其數(shù)學(xué)表達(dá)式、參數(shù)估計(jì)方法、K-S檢驗(yàn)法及RMSE值等，詳見文獻(xiàn)［12-15］。

3.2先驗(yàn)分布考慮到變量Hk-1和Hk具有相同的邊緣分布，在表示概率分布函數(shù)和密度函數(shù)時(shí)統(tǒng)一為隨機(jī)變量H。利用Copula函數(shù)，Hk-1、Hk的聯(lián)合分布可寫為：

先驗(yàn)轉(zhuǎn)移概率分布為給定Hk-1=hk-1時(shí)，Hk的條件分布函數(shù)為［18-19］：

先驗(yàn)轉(zhuǎn)移密度函數(shù)：

3.3似然函數(shù)令Sk的邊緣分布函數(shù)為相應(yīng)的概率密度函數(shù)為借助Copula函數(shù)，Hk-1、Hk和Sk的聯(lián)合分布可寫為：

相應(yīng)的密度函數(shù)：

3.4后驗(yàn)分布將式（7）和式（10）代入式（3），可以推導(dǎo)得到實(shí)際流量Hk后驗(yàn)轉(zhuǎn)移密度函數(shù)的解析表達(dá)式：

相應(yīng)的后驗(yàn)轉(zhuǎn)移分布函數(shù)為

根據(jù)數(shù)理統(tǒng)計(jì)原理，可以計(jì)算得到期望值作為確定性預(yù)報(bào)結(jié)果，同時(shí)獲取給定置信水平下的流量預(yù)報(bào)區(qū)間。實(shí)測(cè)流量的期望值hke通過(guò)下式求解：

令實(shí)測(cè)流量Hk取值出現(xiàn)在分布兩端的概率為ξ，就可以定義Hk的置信水平為（1-ξ）的區(qū)間估計(jì)。實(shí)測(cè)流量Hk置信區(qū)間的兩個(gè)端點(diǎn)分別由以下兩式給出［22］：

因此有

即［hkl，hku］為實(shí)測(cè)流量Hk置信水平（1-ξ）的區(qū)間估計(jì)，根據(jù)置信區(qū)間可以對(duì)Hk的轉(zhuǎn)移預(yù)報(bào)不確定性進(jìn)行定量評(píng)價(jià)。

4　基于Copula函數(shù)的貝葉斯極值預(yù)報(bào)方法

4.1極值概率預(yù)報(bào)的定義令隨機(jī)變量Zk代表時(shí)段(0，k］的最大實(shí)際流量，這個(gè)新定義的一維連續(xù)隨機(jī)變量可以看成是隨機(jī)變量的極大值變量，它由一個(gè)時(shí)間離散的流量隨機(jī)過(guò)程{H1，…，Hk}衍生而得。

對(duì)于任意給定的流量值h，設(shè)為最大實(shí)際流量Zk的超過(guò)概率分布函數(shù)，通過(guò)下式定義：

將式（17）代入，并進(jìn)行演算可得：

4.2CBEF方法的數(shù)學(xué)描述CBEF方法的概念流程示意圖如圖2所示。從式（19）可以看出，對(duì)于任意給定的流量值h，求解超過(guò)概率分布函數(shù)的關(guān)鍵在于計(jì)算得到聯(lián)合概率分布值

圖2　CBEF方法概念示意圖

則式（19）可以重新寫為：

5　實(shí)例應(yīng)用

三峽工程是開發(fā)與治理長(zhǎng)江的關(guān)鍵性骨干工程，具有防洪、發(fā)電、航運(yùn)、補(bǔ)水等綜合效益。三峽水庫(kù)集水面積約為100萬(wàn)km2，壩址多年平均流量為14 300 m3/s。本文基于CBTF方法實(shí)現(xiàn)三峽水庫(kù)汛期入庫(kù)流量的轉(zhuǎn)移概率預(yù)報(bào)，進(jìn)而基于CBEF方法實(shí)現(xiàn)三峽水庫(kù)汛期入庫(kù)流量的極值概率預(yù)報(bào)。選用2003—2016年汛期（6月1日—9月30日）發(fā)生的22場(chǎng)實(shí)測(cè)入庫(kù)洪水過(guò)程及長(zhǎng)江水文氣象預(yù)報(bào)中心發(fā)布的預(yù)見期1～3 d相應(yīng)的確定性預(yù)報(bào)入庫(kù)洪水過(guò)程資料。實(shí)測(cè)值和預(yù)報(bào)值均為每日8時(shí)入庫(kù)流量，取樣時(shí)以預(yù)見期對(duì)應(yīng)的洪水發(fā)生的實(shí)際時(shí)間為基礎(chǔ)，反推相應(yīng)的預(yù)報(bào)發(fā)布時(shí)刻，同步樣本數(shù)均為707次。

5.1邊緣分布的確定邊緣分布函數(shù)可以采用任意分布，選用的標(biāo)準(zhǔn)是使得假定理論分布與流量資料擬合較好。本文采用P-Ⅲ分布作為實(shí)際流量和預(yù)報(bào)流量的邊緣分布，其密度函數(shù)表達(dá)式如下［23］：

式中：α、β和γ分別為形狀、尺度和位置參數(shù)。

考慮到變量Hk（k=1，2，3）的邊緣概率分布相同，將其統(tǒng)一表示為隨機(jī)變量H。采用P-Ⅲ分布線型分別作為H、Sk（k=1，2，3）的理論線型，利用L-矩法估計(jì)參數(shù)，結(jié)果見表1。采用K-S檢驗(yàn)法分別對(duì)4個(gè)隨機(jī)變量進(jìn)行擬合檢驗(yàn)，單變量K-S檢驗(yàn)統(tǒng)計(jì)量D值也列于表1。結(jié)果表明，在5%的顯著性水平下，K-S檢驗(yàn)統(tǒng)計(jì)量D均小于相應(yīng)的臨界值（0.0511），各變量均通過(guò)了檢驗(yàn)。因此，估計(jì)的實(shí)際流量和預(yù)報(bào)流量邊緣分布是合理可行的。

表1　邊緣分布參數(shù)估計(jì)結(jié)果和K-S檢驗(yàn)

5.2聯(lián)合分布的建立變量?jī)蓛芍g的Kendall秩相關(guān)系數(shù)矩陣見表2?？梢园l(fā)現(xiàn)，3個(gè)變量之間的相關(guān)性是不對(duì)稱的，Hk和Sk的相關(guān)性最高，且Hk-1和Hk、Hk-1和Sk的相關(guān)性基本相當(dāng)。

考慮到Hk-1、Hk的聯(lián)合分布與H0、H1的聯(lián)合分布相同，因而只需要估計(jì)H0、H1的聯(lián)合分布。分別采用二維Gumbel-Hougaard、Clayton和Frank Copula函數(shù)分別構(gòu)造H0、H1的聯(lián)合分布，基于不同預(yù)見期的同步系列數(shù)據(jù)，分別得到相應(yīng)的Kendall秩相關(guān)系數(shù)τ，根據(jù)τ與參數(shù)θ的關(guān)系分別計(jì)算Copula函數(shù)的參數(shù)值，結(jié)果見表3。二維K-S檢驗(yàn)統(tǒng)計(jì)量D和RMSE值也列于表3?？芍?，在5%的顯著性水平下，3種備選Copula函數(shù)的二維K-S檢驗(yàn)統(tǒng)計(jì)量D均小于臨界值（0.0511），并且二維Frank Copula函數(shù)的RMSE值均為最小，因而選擇其構(gòu)造H0、H1的聯(lián)合分布。

表2　變量Hk-1，Hk和Sk兩兩之間的Kendall秩相關(guān)系數(shù)矩陣

表3　H0和H1聯(lián)合分布參數(shù)估計(jì)、檢驗(yàn)及優(yōu)選結(jié)果

分別采用三維非對(duì)稱Gumbel-Hougaard、Clayton和Frank Copula函數(shù)分別構(gòu)造Hk-1，Hk和Sk的聯(lián)合分布，基于不同預(yù)見期的同步系列數(shù)據(jù)，利用極大似然法估計(jì)的三維非對(duì)稱Copula函數(shù)參數(shù)值、三維K-S檢驗(yàn)統(tǒng)計(jì)量D和RMSE值列于表4?？梢园l(fā)現(xiàn)，在5%的顯著性水平下，對(duì)于預(yù)見期1 d和2 d三維非對(duì)稱Gumbel-Hougaard Copula函數(shù)和Frank Copula函數(shù)通過(guò)了K-S檢驗(yàn)（臨界值0.0511），而預(yù)見期3 d只有三維非對(duì)稱Frank Copula函數(shù)通過(guò)了檢驗(yàn)。而且，3個(gè)預(yù)見期三維非對(duì)稱Frank Copula函數(shù)均給出了最小的RMSE值，因此選擇其構(gòu)造Hk-1，Hk和Sk的聯(lián)合分布。

表4　Hk-1、Hk和Sk三維非對(duì)稱Copula函數(shù)參數(shù)估計(jì)、檢驗(yàn)及優(yōu)選結(jié)果

5.3轉(zhuǎn)移概率預(yù)報(bào)結(jié)果示例任意給定實(shí)際流量Hk-1和確定性預(yù)報(bào)流量Sk（k=1，2，3）的取值，就可以求解實(shí)際流量Hk的后驗(yàn)轉(zhuǎn)移概率密度和后驗(yàn)轉(zhuǎn)移分布函數(shù)，實(shí)現(xiàn)入庫(kù)流量的轉(zhuǎn)移概率預(yù)報(bào)。在實(shí)際分析Hk的后驗(yàn)不確定性時(shí)，實(shí)際流量Hk-1的取值是未知的，其取值狀態(tài)可以是可行域中的任何值。以預(yù)報(bào)時(shí)刻2003年7月5日（20030705）為例，預(yù)報(bào)時(shí)刻的實(shí)測(cè)流量h0=33 700 m3/s。圖3和圖4分別給出了h0（為已知值）、h1和h2取各種不同狀態(tài)時(shí)，預(yù)見期1～3 d實(shí)際入庫(kù)流量Hk的后驗(yàn)轉(zhuǎn)移密度函數(shù)曲線和后驗(yàn)轉(zhuǎn)移分布函數(shù)曲線?？芍?，后驗(yàn)轉(zhuǎn)移分布呈現(xiàn)正偏態(tài)特征，CBTF方法可以準(zhǔn)確地描述相鄰兩個(gè)預(yù)見期流量之間的后驗(yàn)轉(zhuǎn)移概率密度和分布函數(shù)，從而定量評(píng)估轉(zhuǎn)移預(yù)報(bào)的不確定性。

圖3　三峽水庫(kù)20030705入庫(kù)流量后驗(yàn)轉(zhuǎn)移密度曲線

圖4　三峽水庫(kù)20030705入庫(kù)流量后驗(yàn)轉(zhuǎn)移分布曲線

為了展示入庫(kù)流量預(yù)報(bào)不確定性在時(shí)間上的演變特征，我們把實(shí)際發(fā)生的流量過(guò)程當(dāng)作已知值（h0=33 700 m3/s，h1=36 100 m3/s，h2=38 600 m3/s和h3=43 000 m3/s），對(duì)應(yīng)的確定性預(yù)報(bào)結(jié)果過(guò)程為（s1=35 500 m3/s，s2=38 000 m3/s和s3=37 500 m3/s）。三峽水庫(kù)20030705入庫(kù)流量預(yù)報(bào)先驗(yàn)和后驗(yàn)不確定性隨時(shí)間的演變過(guò)程見圖5，其中紅色實(shí)點(diǎn)代表實(shí)際值，虛線為實(shí)際流量過(guò)程?？梢园l(fā)現(xiàn)，隨著預(yù)見期的延長(zhǎng)或者時(shí)間的推移，入庫(kù)流量預(yù)報(bào)的先驗(yàn)和后驗(yàn)轉(zhuǎn)移密度函數(shù)由陡變緩，表明預(yù)報(bào)不確定性增大，這與實(shí)際情況相符。但利用確定性預(yù)報(bào)信息經(jīng)過(guò)貝葉斯修正之后的后驗(yàn)轉(zhuǎn)移密度函數(shù)與相應(yīng)的先驗(yàn)密度函數(shù)相比都更加尖瘦和集中，不確定性更小。CBTF方法和CMHUP為分析水文預(yù)報(bào)不確定性在時(shí)間上的演變特征提供了有效工具。

圖6給出了三峽水庫(kù)“20160701”號(hào)洪水過(guò)程預(yù)見期1～3 d的后驗(yàn)期望值和90%置信區(qū)間及實(shí)際值?？梢钥闯觯樗髁康牧考?jí)越大，相應(yīng)的不確定性一般也越大。后驗(yàn)期望值預(yù)報(bào)與實(shí)際流量序列擬合效果總體較好，但擬合效果隨預(yù)見期的延長(zhǎng)而降低。此外，90%預(yù)報(bào)區(qū)間也隨著預(yù)見期的延長(zhǎng)均變寬，預(yù)報(bào)不確定性增大，但基本上可以包住實(shí)際流量，表明概率區(qū)間預(yù)報(bào)是可靠的，可以為水庫(kù)防洪調(diào)度決策提供更多的風(fēng)險(xiǎn)信息。

5.4轉(zhuǎn)移概率預(yù)報(bào)結(jié)果評(píng)價(jià)基于整個(gè)研究數(shù)據(jù)集（樣本數(shù)為707），從后驗(yàn)期望值預(yù)報(bào)的精度和概率預(yù)報(bào)的整體性能對(duì)CMHUP的轉(zhuǎn)移概率預(yù)報(bào)結(jié)果進(jìn)行評(píng)價(jià)，并與CHUP的結(jié)果進(jìn)行對(duì)比分析。從理論定義上講，CHUP描述的是k步流量狀態(tài)轉(zhuǎn)移H0→Hk（k=1，2，3）的不確定性，而CMHUP描述的則是單步流量狀態(tài)轉(zhuǎn)移Hk-1→Hk的不確定性。需要指出的是，預(yù)見期1 d時(shí)兩者結(jié)果是相同的，描述的都是狀態(tài)轉(zhuǎn)移H0→H1的不確定性，以下重點(diǎn)分析預(yù)見期2 d、3 d的情形。

圖5　三峽水庫(kù)20030705入庫(kù)流量預(yù)報(bào)不確定性隨時(shí)間的演變過(guò)程

圖6　三峽水庫(kù)“20160701”號(hào)入庫(kù)洪水流量實(shí)測(cè)值、期望值與90%置信區(qū)間

（1）后驗(yàn)期望值預(yù)報(bào)。三峽水庫(kù)預(yù)見期1～3 d的確定性預(yù)報(bào)與CHUP、CMHUP后驗(yàn)期望值預(yù)報(bào)結(jié)果的精度評(píng)價(jià)指標(biāo)R2和RE分別見表5?？梢钥闯?，相比確定性預(yù)報(bào)，CHUP、CMHUP后驗(yàn)期望值預(yù)報(bào)結(jié)果R2均有一定程度的提高，RE顯著地減小，而且CMHUP的后驗(yàn)期望值預(yù)報(bào)相比CHUP略有改善。由于三峽水庫(kù)的防洪庫(kù)容較大，水庫(kù)運(yùn)行調(diào)度主要受洪水總量控制，因而準(zhǔn)確地預(yù)報(bào)徑流總量具有重要意義，相比于確定性預(yù)報(bào)而言，CHUP、CMHUP后驗(yàn)期望值預(yù)報(bào)在這方面顯示出了很大優(yōu)勢(shì)。還可以發(fā)現(xiàn)隨著預(yù)見期的延長(zhǎng)，由于確定性預(yù)報(bào)精度的下降，導(dǎo)致CHUP、CMHUP后驗(yàn)均值預(yù)報(bào)的精度都相應(yīng)地降低。

表5　不同預(yù)見期的確定性預(yù)報(bào)與CHUP、CMHUP后驗(yàn)期望值預(yù)報(bào)結(jié)果比較（單位：%）

（2）概率預(yù)報(bào)。采用可靠性（α-index）、分辨率（π-index）和連續(xù)概率排位分?jǐn)?shù)（CRPS，確定性預(yù)報(bào)退化為平均絕對(duì)誤差MAE）等指標(biāo)來(lái)評(píng)定概率預(yù)報(bào)的整體性能［24-25］。CHUP、CMHUP不同預(yù)見期的α-index、π-index和CRPS等指標(biāo)值列于表6?？芍?，不同預(yù)見期的α-index值均大于0.9，表明概率預(yù)報(bào)結(jié)果可靠性高。此外，隨著預(yù)見期的延長(zhǎng)，π-index值逐漸減小，概率預(yù)報(bào)的分辨率和精度降低，說(shuō)明入庫(kù)流量的預(yù)報(bào)不確定性增加。

表6　不同預(yù)見期的CHUP、CMHUP概率預(yù)報(bào)評(píng)價(jià)指標(biāo)值比較

隨著預(yù)見期的延長(zhǎng)，CHUP、CMHUP計(jì)算的綜合指標(biāo)CRPS值也呈現(xiàn)不斷增加的趨勢(shì)，意味著概率預(yù)報(bào)的性能和總體效果降低。然而，基于CHUP、CMHUP得到的概率預(yù)報(bào)CRPS值，始終小于相應(yīng)的確定性預(yù)報(bào)，彰顯了概率預(yù)報(bào)的有效性。而且，CMHUP相比于CHUP也有不同程度的改善，表明入庫(kù)流量單步狀態(tài)轉(zhuǎn)移的不確定性要低于多步轉(zhuǎn)移。CMHUP計(jì)算的預(yù)見期1～3 d的CRPS值相比于確定性預(yù)報(bào)分別降低11.58%、18.52%和23.33%。

5.5最大入庫(kù)流量超過(guò)概率分布函數(shù)根據(jù)樣本空間同時(shí)考慮入庫(kù)流量概率分布的特點(diǎn)，取累積概率超過(guò)0.9999的流量分位數(shù)值作為hk取值的上限hkmax=100 000 m3/s。因此，本文hk（k=1，2，3）的取值區(qū)間為10 100～100 000 m3/s，離散間距為100 m3/s，相應(yīng)的離散點(diǎn)據(jù)數(shù)為900。以預(yù)報(bào)時(shí)刻2003年7月5日（20030705）為例，對(duì)于給定的流量值h（離散取值和hk相同），根據(jù)計(jì)算的轉(zhuǎn)移概率預(yù)報(bào)結(jié)果，可以計(jì)算得到不同預(yù)見期內(nèi)最大入庫(kù)流量的超過(guò)概率分布函數(shù)超過(guò)概率分布曲線見圖7。

圖7　三峽水庫(kù)20030705不同預(yù)見期內(nèi)最大入庫(kù)流量超過(guò)概率分布曲線

從圖中可以讀取不同預(yù)見期和給定流量值h對(duì)應(yīng)的最大入庫(kù)流量超過(guò)概率，例如1 d內(nèi)最大入庫(kù)流量超過(guò)35 000 m3/s的概率，2 d內(nèi)最大入庫(kù)流量超過(guò)40 000 m3/s的概率3 d內(nèi)最大入庫(kù)流量超過(guò)50 000 m3/s的概率另外，還可以計(jì)算得到3 d內(nèi)最大入庫(kù)流量對(duì)應(yīng)的90%置信區(qū)間為［35100，56700］m3/s，提供了預(yù)見期時(shí)段內(nèi)最大入庫(kù)流量預(yù)報(bào)的不確定性和風(fēng)險(xiǎn)信息。

6　結(jié)論

本文提出了基于Copula函數(shù)的貝葉斯轉(zhuǎn)移預(yù)報(bào)（CBTF）方法和基于Copula函數(shù)的多變量水文不確定性處理器（CMHUP），進(jìn)而發(fā)展了基于Copula函數(shù)的貝葉斯極值預(yù)報(bào)（CBEF）方法。以三峽水庫(kù)為例進(jìn)行應(yīng)用，主要研究結(jié)論如下：（1）CBTF方法不需要進(jìn)行線性-正態(tài)假設(shè)，可以在原始流量數(shù)據(jù)空間給出貝葉斯后驗(yàn)轉(zhuǎn)移密度函數(shù)的解析表達(dá)式，能夠很好地捕捉流量過(guò)程的非線性和非正態(tài)特征，適用范圍更加廣泛。（2）CMHUP可以在給定當(dāng)前實(shí)測(cè)流量和確定性預(yù)報(bào)過(guò)程的條件下，更加靈活地給出考慮實(shí)際流量過(guò)程內(nèi)在相關(guān)性結(jié)構(gòu)的后驗(yàn)聯(lián)合概率密度函數(shù)，為分析水文預(yù)報(bào)不確定性在時(shí)間上的演變特征提供了有效工具。（3）CBEF方法給決策者提供了預(yù)見期時(shí)段內(nèi)最大入庫(kù)流量預(yù)報(bào)的不確定性和風(fēng)險(xiǎn)信息，對(duì)防洪預(yù)報(bào)預(yù)警和減災(zāi)具有重要的參考價(jià)值。

本文采用Copula函數(shù)取代原有基于正態(tài)分位數(shù)轉(zhuǎn)換和線性-正態(tài)假設(shè)的亞高斯模型，不僅在理論上存在優(yōu)勢(shì)，而且實(shí)用有效。然而，本文方法在假定實(shí)際流量過(guò)程服從一階馬爾科夫過(guò)程的基礎(chǔ)上推導(dǎo)得到，對(duì)于流量過(guò)程服從兩階及以上馬爾科夫過(guò)程的情形本文方法的思路仍然適用，但需要采用更高維的Copula函數(shù)，計(jì)算更加復(fù)雜，有待進(jìn)一步研究。另外，如何將定量降水預(yù)報(bào)不確定性與本文所提方法進(jìn)行耦合，估計(jì)轉(zhuǎn)移預(yù)報(bào)和極值預(yù)報(bào)的總不確定性，也需要進(jìn)一步探討。

參考文獻(xiàn)：

［1］KRZYSZTOFOWICZ R.Bayesian theory of probabilistic forecasting via deterministic hydrologic model［J］.Water Resources Research，1999，35（9）：2739-2750.

［2］KRZYSZTOFOWICZ R，KELLY K S.Hydrologic uncertainty processor for probabilistic river stage forecasting［J］.Water Resources Research，2000，36（11）：3265-3277.

［3］李向陽(yáng)，程春田，林劍藝.基于 BP神經(jīng)網(wǎng)絡(luò)的貝葉斯概率水文預(yù)報(bào)模型［J］.水利學(xué)報(bào)，2006，37（3）：354-359.

［4］鄢來(lái)標(biāo)，康玲.河道洪水概率預(yù)報(bào)方法研究［J］.水力發(fā)電，2008，34（12）：40-41.

［5］蔣曉蕾，梁忠民，王春青，等 .BFS-HUP模型在潼關(guān)站洪水概率預(yù)報(bào)中的應(yīng)用［J］.人民黃河，2015，37（7）：13-15.

［6］KRZYSZTOFOWICZ R，MARANZANO C J.Hydrologic uncertainty processor for probabilistic stage transition forecasting［J］.Journal of Hydrology，2004，293（1）：57-73.

［7 ］KRZYSZTOFOWICZ R.A theory of flood warning systems［J］.Water Resources Research，1993，29（12）：3981-3994.

［8］程曉陶，李娜，王艷艷，等 .防汛預(yù)警指標(biāo)與等級(jí)劃分的比較研究［J］.中國(guó)防汛抗旱，2010，20（3）：26-31.

［9 ］KRZYSZTOFOWICZ R.Probabilistic flood forecast：Exact and approximate predictive distributions［J］.Journal of Hydrology，2014，517（1）：643-651.

［10］MADADGAR S，MORADKHANI H.Improved Bayesian multimodeling：Integration of copulas and Bayesian model averaging［J］.Water Resources Research，2014，50（12）：9586-9603.

［11］劉章君，郭生練，李天元，等.貝葉斯概率洪水預(yù)報(bào)模型及其比較應(yīng)用研究［J］.水利學(xué)報(bào)，2014，45（9）：1019-1028.

［12］SALVADORI G，de MICHELE C.On the use of copulas in hydrology：theory and practice［J］.Journal of Hydrologic Engineering，2007，12（4）：369-380.

［13］郭生練，閆寶偉，肖義，等.Copula函數(shù)在多變量水文分析計(jì)算中的應(yīng)用及研究進(jìn)展［J］.水文，2008，28（3）：1-7.

［14］孫鵬，張強(qiáng)，陳曉宏.基于Copula函數(shù)的鄱陽(yáng)湖流域極值流量遭遇頻率及災(zāi)害風(fēng)險(xiǎn)［J］.湖泊科學(xué)，2011，23（2）：183-190.

［15］馮平，李新 .基于Copula函數(shù)的非一致性洪水峰量聯(lián)合分析［J］.水利學(xué)報(bào)，2013，44（10）：1137-1147.

［16］MARANZANO C J，KRZYSZTOFOWICZ R.Identification of likelihood and prior dependence structures for hydrologic uncertainty processor［J］.Journal of Hydrology，2004，290（1-2）：1-21.

［17］ENGELAND K，STEINSLAND I.Probabilistic postprocessing models for flow forecasts for a system of catchments and several lead times［J］.Water Resources Research，2014，50（1）：182-197.

［18］劉曾美，覃光華，陳子燊，等.感潮河段水位與上游洪水和河口潮位的關(guān)聯(lián)性研究［J］.水利學(xué)報(bào)，2013，44（11）：1278-1285.

［19］陳璐，盧韋偉，周建中，等.水文預(yù)報(bào)不確定性對(duì)水庫(kù)防洪調(diào)度的影響分析［J］.水利學(xué)報(bào)，2016，47（1）：77-84.

［20］侯蕓蕓，宋松柏，趙麗娜，等.基于Copula函數(shù)的三變量洪水頻率研究［J］.西北農(nóng)林科技大學(xué)學(xué)報(bào)（自然科學(xué)版），2010，38（2）：219-228.

［21］梁忠民，郭彥，胡義明，等.基于copula函數(shù)的三峽水庫(kù)預(yù)泄對(duì)鄱陽(yáng)湖防洪影響分析［J］.水科學(xué)進(jìn)展，2012，23（4）：485-492.

［22］劉章君，郭生練，李天元，等 .設(shè)計(jì)洪水地區(qū)組成的區(qū)間估計(jì)方法研究［J］.水利學(xué)報(bào)，2015，46（5）：543-550.

［23］中華人民共和國(guó)水利部.水利水電工程設(shè)計(jì)洪水計(jì)算規(guī)范（SL44-2006）［S］.北京：中國(guó)水利水電出版社，2006.

［24］THYER M，RENARD B，KAVETSKI D，et al.Critical evaluation of parameter consistency and predictive uncertainty in hydrological modeling：A case study using Bayesian total error analysis［J］.Water Resources Research，2009，45（12）：1211-1236.

［25］李明亮，楊大文，陳勁松 .基于采樣貝葉斯方法的洪水概率預(yù)報(bào)研究［J］.水力發(fā)電學(xué)報(bào)，2011，30（3）：27-33.