基于Copula函數的高強度降雨設計研究

余 璐

(北京市十三陵水庫管理處,北京 102200)

全球氣候變暖和人類活動可能引發極端氣候事件發生強度和頻率增加，極大地影響了流域的水文特征[1- 2]。強降雨事件作為極端天氣中的一種嚴重致災的因子，容易導致城市洪澇、泥石流和山體滑坡等一系列災害。因此，強降雨事件分析對于城市防災減災有重要現實意義[3]。

水文中通常采用降雨的概率評估水文事件的嚴重程度。通過降雨量超過一定量級的頻率來度量。然而，單一降雨要素不能表征降雨的多重屬性。流域內一場降雨的降雨量越大，降雨歷時越小，則控制流域的水庫越難調節此場降雨，容易引起洪水災害。現有的水庫防洪研究往往只考慮不同頻率的降雨量對水庫防洪的影響，而沒有考慮到降雨歷時，這容易忽略降雨量大而降雨歷時小的高強度降雨。運用多變量能更好描述降雨事件，Copula函數由于可以不受邊緣分布類型拘束，靈活描述多元變量之間復雜的非線性相關結構，在多變量分析中廣泛運用。曹偉華等通過阿基米德Copula函數構建降水持續時間和過程降水量的二元聯合分布模型，計算了北京地區強降水事件條件重現期[4]。何英等基于阿基米德Copula函數分析了葉爾羌河流域卡群水文站洪峰流量與洪量的聯合頻率分布特征[5]。張冬冬等應用阿基米德Copula函數探討了洪峰、洪量和歷時的聯合概率分布和條件概率分布[6]。

北京市十三陵水庫作為溫榆河重要的防洪建筑，開展十三陵水庫控制流域多變量降雨分析對北京市防洪減災工作具有重要意義。本文以北京市十三陵水庫流域為例，重點考慮降雨量大而降雨歷時小的高強度降雨，采用25種常見Copula函數全面分析降雨量與降雨歷時倒數之間的二元聯合分布規律，并以此為基礎設計出對應頻率下的最可能降雨場景，以期為北京地區防洪降險提供科學依據和推求高強度降雨的設計降雨提供技術支持。

1 研究區概況

十三陵水庫位于溫榆河北支東沙河，流域面積223km2，控制東沙河84%的流域面積，占溫榆河上游山區流域面積1/4。該流域為典型的暖溫帶半濕潤大陸性季風氣候，夏季炎熱多雨，冬季寒冷干燥，年平均氣溫10～12℃，年平均降雨量550mm，降水季節分配不均，全年降水的75%集中在夏季，7、8月常有暴雨。水庫總庫容7450萬m3，防洪設計洪水標準為100年一遇，校核洪水標準為2000年一遇。水庫上游有德勝口、錐石口、上下口及老君堂四條溝。

十三陵水庫流域設有老君堂、果莊、上口、碓臼峪、十三陵水庫5個雨量站。選取5個雨量站1980—2020年6—9月汛期每場降雨的雨量和歷時，采用算數平均法計算5個雨量站的平均降雨量作為十三陵水庫流域的降雨量。

2 研究方法

2.1 邊緣分布函數

邊緣分布函數的確定是構建降雨歷時倒數和降雨量聯合分布的前提。采用應用較為廣泛的16種邊緣分布函數擬合降雨量和歷時的倒數的邊緣分布，并運用極大似然法估計邊緣分布的參數。16種邊緣分布函數包括Birnbaum-Saunders分布、指數分布(exponential)、極值分布(extreme value)、伽瑪分布(Gamma)、廣義極值分布(generalized extreme value)、廣義帕累托分布(generalized Pareto)、逆高斯分布(inverse Gaussian)、邏輯分布(logistic)、對數邏輯分布(log-logistic)、對數正態分布(lognormal)、Nakagami 分布、正態分布(normal)、瑞利分布(Rayleigh)、萊斯分布(Rician)、t位置尺度分布(t location scale)和韋伯分布(Weibull)，模型詳細情況參考文獻[7- 10]。

為了識別出最佳的邊緣分布函數，需要對擬合的結果進行統計檢驗。本文運用赤池信息準則(AIC)和貝葉斯信息準則(BIC)對16種邊緣分布的擬合結果進行優選，識別出降雨量和歷時的倒數最適合的邊緣分布函數。分別以AIC和BIC值最小為最優。AIC和BIC公式如下：

AIC=-2In(L)+2k

(1)

BIC=-2In(L)+In(n)×k

(2)

式中，L—極大似然函數；k—模型的變量的個數；n—數據數量。

2.2 Copula函數

Copula函數描述的是變量之間的相關性，是一種將變量聯合累積分布函數同變量邊緣累積分布函數連接起來的函數。設H為聯合分布函數，F和G為其邊緣分布，那么存在唯一的Copula函數C，有[11]：

H(x，y)=C[U(X)，V(Y)]

(3)

Copula函數在水文領域應用廣泛。在研究中使用最多的Copula函數主要有阿基米德Copula族和橢圓Copula族兩大類[11]。為使備選的Copula更全面，本文選用25種具有代表性的Copula函數構建歷時倒數和降雨量的聯合分布函數，并采用極大似然法估計參數θ。25種Copula函數包括Gaussian、t 2種橢圓Copula函數，Clayton、Frank、Gumbel、Independence、AMH、Joe、FGM、Plackett、Cuadras-Auge、Raftery、Shih-Louis、Linear-Spearman、Cubic、Burr、Nelsen、Galambos、Marshal-Olkin、Fischer-Hinzmann、Roch-Alegre、Fischer-Kock、BB1、BB5、Tawn 23種Archimedean Copula函數，Copula函數、參數θ取值范圍及聯合分布函數見表1。

表1 Copula函數基本特征

為選擇最優的Copula函數來描述變量間的相關結構，需要對各函數進行擬合度評價和檢驗。選用極大似然法(ML)、赤池信息準則(AIC)和貝葉斯信息準則(BIC)檢驗Copula函數。以ML值最大和AIC、BIC值最小作為擬合最優的檢驗標準。

3 結果與分析

3.1 降雨量與歷時倒數的邊緣分布

為了確定降雨量和歷時倒數的邊緣分布，采用極大似然法估計降雨量與歷時倒數的邊緣分布參數，并利用AIC和BIC對十三陵流域降雨量與歷時倒數的16種邊緣分布進行擬合度檢驗。

十三陵流域的歷時倒數和降雨量16種邊緣分布函數的AIC和BIC檢驗結果見表2—3。由表2可得，歷時倒數邊緣分布函數檢驗中generalized pareto分布的AIC值和BIC值最小，由表3可得，降雨量邊緣分布函數檢驗中lognormal分布的AIC值和BIC值最小，即generalized pareto分布和lognormal分布分別作為十三陵流域歷時倒數和降雨量最優的邊緣分布。繪制優選的擬合分布與經驗分布的匹配程度的比較圖(圖1a、2a)和變量邊緣分布經驗頻率與理論頻率關系圖(圖1b、2b)，如圖1a、2a所示，點據均分布在擬合線附近，同時，由圖1b、2b可得，樣本分位數點和相應的理論上擬合的分位數接近。兩者都直觀得證明擬合結果較好。

表2 歷時倒數的邊緣分布AIC和BIC檢驗結果

表3 降雨量的邊緣分布AIC和BIC檢驗結果

圖1 歷時倒數的邊緣分布擬合圖

3.2 二元Copula分布選擇

在優選出降雨量和歷時倒數的最優邊緣分布后，運用Copula函數構建降雨量和歷時倒數的聯合分布。本文運用極大似然法對25種Copula函數進行參數估算，通過極大似然估計方法找出使估計的聯合概率與相關的經驗聯合概率最匹配的Copula參數。然后根據ML和AIC、BIC比較25種Copula函數的擬合優度，從而確定十三陵流域的Copula函數類型。

十三陵水庫流域的降雨量和歷時倒數的檢驗結果見表3，Roch-Alegre的ML值最大為1150.2，AIC和BIC值最小分別為-2296.4和-2289.0，即ML、AIC、BIC分析法3種評價方法均得出一致的結果，因此，采用Roch-Alegre函數描述十三陵流域降雨量和歷時倒數的聯合分布。Roch-Alegre函數參數值θ1為0.0001，為1.0973，十三陵流域降雨量與歷時倒數聯合分布函數表達式如下：H[u，v]=exp{1-[(((1-ln(u))0.0001-1)1.0973+((1-ln(v))0.0001-1)1.0973)1/1.0973+1]1/0.0001}

圖2 降雨量的邊緣分布擬合圖

表3 Copula函數檢驗結果

3.3 不同設計頻率下的降雨量和降雨歷時分析

本文從二元Copula族中選擇最佳的函數Roch-Alegre來建立降雨量與降雨歷時倒數的依賴結構，將只考慮降雨量的單變量的設計頻率概念擴展到同時考慮降雨量與降雨歷時的多變量。這種方法對于設計降雨量大而降雨歷時短這種極其危險的降雨具有重要的現實意義。圖3為Roch-Alegre函數的聯合分布圖，圖中的等高線對應著相應設計頻率的降雨，同一個等高線上的點構成降雨量大而降雨歷時短二元風險場景。比如，聯合分布Copula值為0.9的等高線對應著無數個設計頻率為0.1的降雨量與降雨歷時組合。降雨量與降雨歷時的無限組合會在統計上是相似的，而它們的影響可能截然不同。如圖3中所示的十三陵水庫的降雨，同一設計頻率降雨下的降雨量與降雨歷時具有無限種組合，在曲線的右側對應著降雨歷時短，降雨量小的場景，而在曲線的上側對應著降雨歷時短，降雨量大的場景。選擇哪一種組合是亟需解決的問題。從無限的降雨量與降雨歷時組合中進行選擇的一種直觀方法是根據它們相關的Copula密度值為它們分配權重。Roch-Alegre函數的概率密度值可以衡量這無限種組合的權重，因此本文采用密度函數極大值處的組合為設計頻率下的降雨量與降雨歷時組合，即選擇一種最可能的事件。

圖3 十三陵水庫流域降雨量和降雨歷時倒數的Roch-Alegre函數聯合分布圖

通過對Roch-Alegre密度函數求導找到極值點，即設計頻率下最可能降雨量與降雨歷時組合場景。設計頻率0.1條件下的最佳組合如圖4所示，降雨量為57mm，降雨歷時倒數為1.3(1/h)。而不考慮降雨量與降雨歷時依賴關系各自得到設計頻率為0.1條件下的降雨量為29mm，降雨歷時倒數為1.2mm。可以看出十三陵水庫控制流域容易出現短歷時大降雨量，即降雨強度大的情況。而這種情況會給防洪帶來很大的風險，所以在設計降雨時就需要提前考慮這種風險疊加。

圖4 設計頻率0.1條件下降雨量和歷時倒數的最佳組合

綜上所述，運用本文提出的考慮降雨量與歷時依賴關系的方法可以選擇具有聯合頻率的危險組合事件，可以用于危險降雨的設計與評估，在十三陵的防洪中具有重要的現實意義。

4 結論

(1)通過擬合優度評價，得出北京市十三陵水庫流域降雨量采用lognormal分布擬合效果較好，而降雨歷時的倒數采用generalized pareto分布擬合效果較好，建議選擇Roch-Alegre函數描述降雨量和歷時倒數的聯合分布。

(2)采用雙變量設計降雨模型與分開考慮降雨量與歷時的傳統方法相比，得到的設計降雨量偏大，降雨歷時偏小。

(3)本文所提的考慮降雨量與歷時依賴關系設計降雨方法，有助于應對降雨量大而降雨歷時短的二元風險，在水庫防洪中具有現實意義，以期為推求高強度降雨的設計降雨提供技術支持。本文未探討短歷時強降雨所形成的洪水過程，有待繼續研究。