Copula-MCP洪水概率預報方法研究與應用

沈 婕，梁忠民，胡義明，王 軍，李彬權

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

0 引 言

水文模型都是對自然子流域的一種概化，無論結構是否復雜、參數是否準確，都只是對水文過程的近似，無法真實地還原自然的水文過程，因此實時洪水的預報中不可避免地存在不確定性，用于描述這種不確定性的洪水概率預報逐漸得到了重視和應用[1,2]。

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

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

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

1 基于Copula-MCP的洪水概率預報模型

1.1 Copula函數

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

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

(1)

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

Copula函數將變量的相關性結構與邊緣分布分開處理，可以考慮變量間的非線性與非正態關系，沒有任何限制。

本文采用Archimedean Copula函數族來構造實測流量與模擬流量的聯合分布函數。常見的Archimedean Copula函數有Gumbel-Hougaard Copula函數、Clayton Copula函數和Frank Copula函數，其二維表達式如下：

Gumbel-Hougaard Copula函數：

(2)

Clayton Copula函數：

(3)

Frank Copula函數：

(4)

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

為了對Copula函數的擬合優度進行檢驗，本文選用K-S檢驗法，其統計量D計算公式如下：

(5)

聯合分布的經驗累積分布與理論分布的擬合情況通過RMSE準則來評價，并根據最小RMSE值來優選Copula函數。RMSE準則計算公式如下：

(6)

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

1.2 模型條件處理器(MCP)

(7)

式中：n為歷史數據的數量；i為變量的排位順序。

(8)

(9)

(10)

預測不確定性一般定義為以模型預測值為條件的預報量的條件分布，在已知變量聯合分布和邊緣分布時，通過計算聯合分布與邊緣分布的比值，就可以推得預測不確定性，即預報量的條件概率分布：

(11)

由此可知，在正態空間中，預測不確定性相當于均值和方差如下的正態分布：

(12)

1.3 Copula-MCP

(13)

(14)

條件概率分布函數對應的概率密度函數為：

(15)

2 應用實例

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

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

2.1 確定性預報

2.1.1 區間洪水預報

采用加權平均法計算王家壩、光山等15個雨量站的區間平均雨量，并采用流量起漲前30天雨量計算流域前期影響雨量Pa，徑流深R采用降雨徑流相關圖法進行計算，降雨徑流相關圖如圖2所示。

圖2 淮河息潢班～王家壩區間降雨徑流相關圖

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

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

2.1.2 河道洪水預報

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

表1 河道洪水演算參數表

2.2 API模型的確定性預報結果

對息潢班～王家壩區間1990-2010年間的25場洪水進行模擬，其中20場作為率定，5場作為驗證，計算時段Δt=6 h，精度統計見表2。可以看出，確定性預報的結果良好，洪峰的合格率(相對誤差絕對值在20%以內)為68%，確定性系數的優秀率(確定性系數大于0.85)為84%。

表2 API模型精度統計表

2.3 Copula-MCP方法應用

2.3.1 邊緣分布的確定

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

表3 y和邊際分布參數估計值和K-S檢驗

圖4 理論Weibull與樣本點據擬合圖

2.3.2 聯合分布的建立

(16)

2.3.3 洪水概率預報結果及分析

表4 聯合分布參數估計、檢驗及優選結果

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

表5 API、MCP與Copula-MCP三個模型結果對比表

圖5 19980509場次洪水Copula-MCP概率預報結果

3 結 論

(1)與目前常用的洪水概率預報面向相較，采用Copula函數直接構建實測流量與預報流量的聯合分布函數，不需要進行變量的正態化變換及線性假設，避免變量間相關信息的丟失和洪水過程非線性特征的改變，由此構建的洪水概率預報模型，理論上應更具優勢。

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

□