沈 婕，梁忠民，胡義明，王 軍，李彬權(quán)

(河海大學(xué)水文水資源學(xué)院，南京 210098)

0 引 言

水文模型都是對自然子流域的一種概化，無論結(jié)構(gòu)是否復(fù)雜、參數(shù)是否準(zhǔn)確，都只是對水文過程的近似，無法真實(shí)地還原自然的水文過程，因此實(shí)時(shí)洪水的預(yù)報(bào)中不可避免地存在不確定性，用于描述這種不確定性的洪水概率預(yù)報(bào)逐漸得到了重視和應(yīng)用[1,2]。

洪水概率預(yù)報(bào)途徑可以概括為兩類[3]：一是總誤差分析途徑，直接在確定性預(yù)報(bào)結(jié)果的基礎(chǔ)上，量化分析最終的預(yù)報(bào)不確定性，進(jìn)而實(shí)現(xiàn)概率預(yù)報(bào)；二是不確定性要素耦合途徑，即分析預(yù)報(bào)過程中各環(huán)節(jié)的主要因子的不確定性，估計(jì)其概率分布，再將這些因子的概率分布耦合到洪水預(yù)報(bào)模型中，從而實(shí)現(xiàn)概率預(yù)報(bào)。從實(shí)時(shí)洪水預(yù)報(bào)角度，第一類的總誤差分析途徑計(jì)算相對簡單，更滿足時(shí)效性要求，因此在洪水作業(yè)預(yù)報(bào)中具有廣泛應(yīng)用。當(dāng)然，不管是哪類途徑，都是在確定性預(yù)報(bào)(定值預(yù)報(bào))的基礎(chǔ)上實(shí)現(xiàn)概率預(yù)報(bào)，因此，如何提高原始定值預(yù)報(bào)結(jié)果的精度，是降低不確定性、提升概率預(yù)報(bào)可靠性或合理性的關(guān)鍵。

在總誤差分析途徑框架下，由Krzysztofowicz等[1]提出的貝葉斯預(yù)報(bào)系統(tǒng)(BFS)是最典型的概率預(yù)報(bào)方法，其中水文不確定性處理器[2](HUP)在貝葉斯理論方法的基礎(chǔ)上，對預(yù)報(bào)變量進(jìn)行了正態(tài)化變換與線性假設(shè)，推求預(yù)報(bào)量后驗(yàn)分布函數(shù)。王善序[4]對貝葉斯概率水文預(yù)報(bào)理論做了簡單的介紹，認(rèn)為其可以不需任何附加假定與要求地定量考慮預(yù)報(bào)的不確定性，同時(shí)也指出了其中的線性-正態(tài)假設(shè)在水文過程中未必適用。邢貞相等[5]采用BP網(wǎng)絡(luò)描述了水文過程的非線性特征，并應(yīng)用AM-MCMC方法進(jìn)行概率預(yù)報(bào)，在此過程中仍采用了正態(tài)性假設(shè)。Todini等[6]提出了模型條件處理器(MCP)，采用截尾正態(tài)聯(lián)合分布揭示了不同流量量級(jí)時(shí)預(yù)報(bào)值與實(shí)測值的關(guān)系。梁忠民等提出了考慮誤差異分布[7](Heterogeneity of Error Distributions, EHDA)方法，分析了預(yù)報(bào)誤差的異分布性，進(jìn)而實(shí)現(xiàn)洪水概率預(yù)報(bào)。

上述這些方法中，均需對原始變量進(jìn)行正態(tài)化變換，并/或?qū)φ龖B(tài)空間中的變量或似然函數(shù)進(jìn)行線性假設(shè)。然而，正態(tài)化變換過程中可能會(huì)造成信息的丟失或偏差[8]，變量或似然函數(shù)的線性假設(shè)也具有較強(qiáng)的主觀性。為避免正態(tài)化變換和線性假設(shè)，劉章君等[9]推導(dǎo)了基于Copula函數(shù)的貝葉斯轉(zhuǎn)移預(yù)報(bào)(BTF)方法的后驗(yàn)密度函數(shù)，提出了CBTF與CMHUP方法，取得了更好的概率預(yù)報(bào)效果。本文采用這種思路，將Copula函數(shù)與MCP方法結(jié)合，構(gòu)建了Copula-MCP洪水概率預(yù)報(bào)模型，并以淮河王家壩斷面為例進(jìn)行了應(yīng)用研究。

1 基于Copula-MCP的洪水概率預(yù)報(bào)模型

1.1 Copula函數(shù)

Copula函數(shù)[10,11]是將隨機(jī)變量的聯(lián)合分布函數(shù)與各自的邊緣分布函數(shù)相連接的連接函數(shù)。設(shè)X、Y為連續(xù)隨機(jī)變量，其邊緣分布函數(shù)分別為Fx(x)與Fy(y)，聯(lián)合分布函數(shù)為F(x,y)，若Fx(x)與Fy(y)連續(xù)，則存在唯一確定的Copula函數(shù)Cθ(u,v)使得：

F(x,y)=Cθ(u,v)

(1)

式中：θ為Copula函數(shù)待定參數(shù)，u=Fx(x)，v=Fy(y)。

Copula函數(shù)將變量的相關(guān)性結(jié)構(gòu)與邊緣分布分開處理，可以考慮變量間的非線性與非正態(tài)關(guān)系，沒有任何限制。

本文采用Archimedean Copula函數(shù)族來構(gòu)造實(shí)測流量與模擬流量的聯(lián)合分布函數(shù)。常見的Archimedean Copula函數(shù)有Gumbel-Hougaard Copula函數(shù)、Clayton Copula函數(shù)和Frank Copula函數(shù)，其二維表達(dá)式如下：

Gumbel-Hougaard Copula函數(shù)：

(2)

Clayton Copula函數(shù)：

(3)

Frank Copula函數(shù)：

(4)

式中：θ為Copula函數(shù)待定參數(shù)；u與v分別為單變量的邊緣分布函數(shù)。

為了對Copula函數(shù)的擬合優(yōu)度進(jìn)行檢驗(yàn)，本文選用K-S檢驗(yàn)法，其統(tǒng)計(jì)量D計(jì)算公式如下：

(5)

聯(lián)合分布的經(jīng)驗(yàn)累積分布與理論分布的擬合情況通過RMSE準(zhǔn)則來評(píng)價(jià)，并根據(jù)最小RMSE值來優(yōu)選Copula函數(shù)。RMSE準(zhǔn)則計(jì)算公式如下：

(6)

式中：n為樣本容量；Fc(xi)為樣本集的累積分布函數(shù)；F(xi)為理論分布。

1.2 模型條件處理器(MCP)

(7)

式中：n為歷史數(shù)據(jù)的數(shù)量；i為變量的排位順序。

(8)

(9)

(10)

預(yù)測不確定性一般定義為以模型預(yù)測值為條件的預(yù)報(bào)量的條件分布，在已知變量聯(lián)合分布和邊緣分布時(shí)，通過計(jì)算聯(lián)合分布與邊緣分布的比值，就可以推得預(yù)測不確定性，即預(yù)報(bào)量的條件概率分布：

(11)

由此可知，在正態(tài)空間中，預(yù)測不確定性相當(dāng)于均值和方差如下的正態(tài)分布：

(12)

1.3 Copula-MCP

(13)

(14)

條件概率分布函數(shù)對應(yīng)的概率密度函數(shù)為：

(15)

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

如圖1所示。王家壩站是淮河上游的總控制站，其洪水由上游干流的息縣、左右岸控制站班臺(tái)和潢川、息-潢-班至王家壩區(qū)間共四部分組成。息潢班-王家壩區(qū)間集水面積為7 110 km2，干流長度為360 km。本文進(jìn)行王家壩洪水預(yù)報(bào)時(shí)，先采用經(jīng)驗(yàn)降雨徑流模型API，對息潢班-王家壩區(qū)間洪水進(jìn)行預(yù)報(bào)，再與經(jīng)馬斯京根法演算至王家壩的息縣、潢川、班臺(tái)實(shí)測流量過程疊加，得到確定性的洪水預(yù)報(bào)結(jié)果。在此基礎(chǔ)上，采用Copula-MCP模型推求以確定性預(yù)報(bào)結(jié)果為條件的預(yù)報(bào)變量的概率分布，實(shí)現(xiàn)王家壩站的洪水概率預(yù)報(bào)。

圖1 淮河息潢班～王家壩區(qū)間流域示意圖

2.1 確定性預(yù)報(bào)

2.1.1 區(qū)間洪水預(yù)報(bào)

采用加權(quán)平均法計(jì)算王家壩、光山等15個(gè)雨量站的區(qū)間平均雨量，并采用流量起漲前30天雨量計(jì)算流域前期影響雨量Pa，徑流深R采用降雨徑流相關(guān)圖法進(jìn)行計(jì)算，降雨徑流相關(guān)圖如圖2所示。

圖2 淮河息潢班～王家壩區(qū)間降雨徑流相關(guān)圖

以降雨中心位置為依據(jù)，采用不同的單位線進(jìn)行區(qū)間匯流計(jì)算。單位線參數(shù)如圖3所示。其中，a線適用于暴雨中心在淮北；b線適用于降雨均勻；c線適用于降雨中心在淮南。

圖3 淮河息潢班～王家壩區(qū)間單位線圖

2.1.2 河道洪水預(yù)報(bào)

采用馬斯京根法，選用先演后合的方式對息縣、潢川、班臺(tái)(大洪河)上游來水進(jìn)行河道洪水演算。相關(guān)參數(shù)見表1。

表1 河道洪水演算參數(shù)表

2.2 API模型的確定性預(yù)報(bào)結(jié)果

對息潢班～王家壩區(qū)間1990-2010年間的25場洪水進(jìn)行模擬，其中20場作為率定，5場作為驗(yàn)證，計(jì)算時(shí)段Δt=6 h，精度統(tǒng)計(jì)見表2?？梢钥闯觯_定性預(yù)報(bào)的結(jié)果良好，洪峰的合格率(相對誤差絕對值在20%以內(nèi))為68%，確定性系數(shù)的優(yōu)秀率(確定性系數(shù)大于0.85)為84%。

表2 API模型精度統(tǒng)計(jì)表

2.3 Copula-MCP方法應(yīng)用

2.3.1 邊緣分布的確定

對實(shí)測及模型計(jì)算的洪水序列(或兩個(gè)隨機(jī)變量)，一般可假設(shè)其邊際分布服從Log-Weibull分布[13]。本次采用線性矩法估計(jì)兩個(gè)變量邊際分布的統(tǒng)計(jì)參數(shù)，為檢驗(yàn)參數(shù)估計(jì)結(jié)果是否合理，采用K-S檢驗(yàn)法分別對兩個(gè)隨機(jī)變量經(jīng)驗(yàn)累積分布與理論分布進(jìn)行擬合優(yōu)度檢驗(yàn)，邊際分布參數(shù)估計(jì)值與K-S檢驗(yàn)統(tǒng)計(jì)量D值結(jié)果見表3。由表3可知，在5%的顯著性水平下，兩個(gè)變量的K-S檢驗(yàn)統(tǒng)計(jì)量D均小于相應(yīng)的臨界值(0.039 5)，均通過擬合優(yōu)度檢驗(yàn)。實(shí)測流量 與預(yù)報(bào)流量 樣本點(diǎn)據(jù)與其邊際分布擬合情況如圖4所示。由圖4可知，理論邊際分布與樣本點(diǎn)據(jù)擬合較好，則使用該分布作為邊際分布是合理的。

表3 y和邊際分布參數(shù)估計(jì)值和K-S檢驗(yàn)

圖4 理論Weibull與樣本點(diǎn)據(jù)擬合圖

2.3.2 聯(lián)合分布的建立

(16)

2.3.3 洪水概率預(yù)報(bào)結(jié)果及分析

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

采用Copula-MCP方法進(jìn)行洪水概率預(yù)報(bào)，其結(jié)果見表5，其中提供了90%置信區(qū)間覆蓋率、洪水概率分布期望值過程的確定性系數(shù)(DC)、洪峰期望值相對誤差等結(jié)果；為對比分析，表中亦列出了API及MCP模型的結(jié)果。圖5為其中一場洪水(19980509)的API、MCP、Copula-MCP預(yù)報(bào)結(jié)果與實(shí)測流量的對比圖。由表5可知，Copula-MCP模型的洪峰相對誤差合格率為92%，高于MCP模型的84%與API模型的68%；Copula-MCP模型的平均確定性系數(shù)為0.93，高于MCP模型的0.92與確定性模型的0.91，Copula-MCP模型確定性系數(shù)的優(yōu)秀率為92%，亦高于MCP模型的88%與API模型的84%；從90%置信區(qū)間的覆蓋率來看，Copula-MCP模型也優(yōu)于MCP模型。

表5 API、MCP與Copula-MCP三個(gè)模型結(jié)果對比表

圖5 19980509場次洪水Copula-MCP概率預(yù)報(bào)結(jié)果

3 結(jié) 論

(1)與目前常用的洪水概率預(yù)報(bào)面向相較，采用Copula函數(shù)直接構(gòu)建實(shí)測流量與預(yù)報(bào)流量的聯(lián)合分布函數(shù)，不需要進(jìn)行變量的正態(tài)化變換及線性假設(shè)，避免變量間相關(guān)信息的丟失和洪水過程非線性特征的改變，由此構(gòu)建的洪水概率預(yù)報(bào)模型，理論上應(yīng)更具優(yōu)勢。

(2)洪水過程的不同發(fā)展階段或不同流量量級(jí)的洪水，其預(yù)報(bào)誤差規(guī)律可能不同，將描述這種異誤差分布的MCP模型與Copula函數(shù)結(jié)合，構(gòu)建的Copula-MCP洪水概率預(yù)報(bào)模型具有更好的適用性。本文研究實(shí)例表明，Copula-MCP模型整體上優(yōu)于MCP模型，可進(jìn)一步提高洪水預(yù)報(bào)精度、增加概率預(yù)報(bào)的可靠性。

□