余 璐
(北京市十三陵水庫(kù)管理處,北京 102200)
全球氣候變暖和人類活動(dòng)可能引發(fā)極端氣候事件發(fā)生強(qiáng)度和頻率增加,極大地影響了流域的水文特征[1- 2]。強(qiáng)降雨事件作為極端天氣中的一種嚴(yán)重致災(zāi)的因子,容易導(dǎo)致城市洪澇、泥石流和山體滑坡等一系列災(zāi)害。因此,強(qiáng)降雨事件分析對(duì)于城市防災(zāi)減災(zāi)有重要現(xiàn)實(shí)意義[3]。
水文中通常采用降雨的概率評(píng)估水文事件的嚴(yán)重程度。通過(guò)降雨量超過(guò)一定量級(jí)的頻率來(lái)度量。然而,單一降雨要素不能表征降雨的多重屬性。流域內(nèi)一場(chǎng)降雨的降雨量越大,降雨歷時(shí)越小,則控制流域的水庫(kù)越難調(diào)節(jié)此場(chǎng)降雨,容易引起洪水災(zāi)害。現(xiàn)有的水庫(kù)防洪研究往往只考慮不同頻率的降雨量對(duì)水庫(kù)防洪的影響,而沒(méi)有考慮到降雨歷時(shí),這容易忽略降雨量大而降雨歷時(shí)小的高強(qiáng)度降雨。運(yùn)用多變量能更好描述降雨事件,Copula函數(shù)由于可以不受邊緣分布類型拘束,靈活描述多元變量之間復(fù)雜的非線性相關(guān)結(jié)構(gòu),在多變量分析中廣泛運(yùn)用。曹偉華等通過(guò)阿基米德Copula函數(shù)構(gòu)建降水持續(xù)時(shí)間和過(guò)程降水量的二元聯(lián)合分布模型,計(jì)算了北京地區(qū)強(qiáng)降水事件條件重現(xiàn)期[4]。何英等基于阿基米德Copula函數(shù)分析了葉爾羌河流域卡群水文站洪峰流量與洪量的聯(lián)合頻率分布特征[5]。張冬冬等應(yīng)用阿基米德Copula函數(shù)探討了洪峰、洪量和歷時(shí)的聯(lián)合概率分布和條件概率分布[6]。
北京市十三陵水庫(kù)作為溫榆河重要的防洪建筑,開展十三陵水庫(kù)控制流域多變量降雨分析對(duì)北京市防洪減災(zāi)工作具有重要意義。本文以北京市十三陵水庫(kù)流域?yàn)槔?,重點(diǎn)考慮降雨量大而降雨歷時(shí)小的高強(qiáng)度降雨,采用25種常見Copula函數(shù)全面分析降雨量與降雨歷時(shí)倒數(shù)之間的二元聯(lián)合分布規(guī)律,并以此為基礎(chǔ)設(shè)計(jì)出對(duì)應(yīng)頻率下的最可能降雨場(chǎng)景,以期為北京地區(qū)防洪降險(xiǎn)提供科學(xué)依據(jù)和推求高強(qiáng)度降雨的設(shè)計(jì)降雨提供技術(shù)支持。
十三陵水庫(kù)位于溫榆河北支東沙河,流域面積223km2,控制東沙河84%的流域面積,占溫榆河上游山區(qū)流域面積1/4。該流域?yàn)榈湫偷呐瘻貛О霛駶?rùn)大陸性季風(fēng)氣候,夏季炎熱多雨,冬季寒冷干燥,年平均氣溫10~12℃,年平均降雨量550mm,降水季節(jié)分配不均,全年降水的75%集中在夏季,7、8月常有暴雨。水庫(kù)總庫(kù)容7450萬(wàn)m3,防洪設(shè)計(jì)洪水標(biāo)準(zhǔn)為100年一遇,校核洪水標(biāo)準(zhǔn)為2000年一遇。水庫(kù)上游有德勝口、錐石口、上下口及老君堂四條溝。
十三陵水庫(kù)流域設(shè)有老君堂、果莊、上口、碓臼峪、十三陵水庫(kù)5個(gè)雨量站。選取5個(gè)雨量站1980—2020年6—9月汛期每場(chǎng)降雨的雨量和歷時(shí),采用算數(shù)平均法計(jì)算5個(gè)雨量站的平均降雨量作為十三陵水庫(kù)流域的降雨量。
邊緣分布函數(shù)的確定是構(gòu)建降雨歷時(shí)倒數(shù)和降雨量聯(lián)合分布的前提。采用應(yīng)用較為廣泛的16種邊緣分布函數(shù)擬合降雨量和歷時(shí)的倒數(shù)的邊緣分布,并運(yùn)用極大似然法估計(jì)邊緣分布的參數(shù)。16種邊緣分布函數(shù)包括Birnbaum-Saunders分布、指數(shù)分布(exponential)、極值分布(extreme value)、伽瑪分布(Gamma)、廣義極值分布(generalized extreme value)、廣義帕累托分布(generalized Pareto)、逆高斯分布(inverse Gaussian)、邏輯分布(logistic)、對(duì)數(shù)邏輯分布(log-logistic)、對(duì)數(shù)正態(tài)分布(lognormal)、Nakagami 分布、正態(tài)分布(normal)、瑞利分布(Rayleigh)、萊斯分布(Rician)、t位置尺度分布(t location scale)和韋伯分布(Weibull),模型詳細(xì)情況參考文獻(xiàn)[7- 10]。
為了識(shí)別出最佳的邊緣分布函數(shù),需要對(duì)擬合的結(jié)果進(jìn)行統(tǒng)計(jì)檢驗(yàn)。本文運(yùn)用赤池信息準(zhǔn)則(AIC)和貝葉斯信息準(zhǔn)則(BIC)對(duì)16種邊緣分布的擬合結(jié)果進(jìn)行優(yōu)選,識(shí)別出降雨量和歷時(shí)的倒數(shù)最適合的邊緣分布函數(shù)。分別以AIC和BIC值最小為最優(yōu)。AIC和BIC公式如下:
AIC=-2In(L)+2k
(1)
BIC=-2In(L)+In(n)×k
(2)
式中,L—極大似然函數(shù);k—模型的變量的個(gè)數(shù);n—數(shù)據(jù)數(shù)量。
Copula函數(shù)描述的是變量之間的相關(guān)性,是一種將變量聯(lián)合累積分布函數(shù)同變量邊緣累積分布函數(shù)連接起來(lái)的函數(shù)。設(shè)H為聯(lián)合分布函數(shù),F(xiàn)和G為其邊緣分布,那么存在唯一的Copula函數(shù)C,有[11]:
H(x,y)=C[U(X),V(Y)]
(3)
Copula函數(shù)在水文領(lǐng)域應(yīng)用廣泛。在研究中使用最多的Copula函數(shù)主要有阿基米德Copula族和橢圓Copula族兩大類[11]。為使備選的Copula更全面,本文選用25種具有代表性的Copula函數(shù)構(gòu)建歷時(shí)倒數(shù)和降雨量的聯(lián)合分布函數(shù),并采用極大似然法估計(jì)參數(shù)θ。25種Copula函數(shù)包括Gaussian、t 2種橢圓Copula函數(shù),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函數(shù),Copula函數(shù)、參數(shù)θ取值范圍及聯(lián)合分布函數(shù)見表1。
表1 Copula函數(shù)基本特征
為選擇最優(yōu)的Copula函數(shù)來(lái)描述變量間的相關(guān)結(jié)構(gòu),需要對(duì)各函數(shù)進(jìn)行擬合度評(píng)價(jià)和檢驗(yàn)。選用極大似然法(ML)、赤池信息準(zhǔn)則(AIC)和貝葉斯信息準(zhǔn)則(BIC)檢驗(yàn)Copula函數(shù)。以ML值最大和AIC、BIC值最小作為擬合最優(yōu)的檢驗(yàn)標(biāo)準(zhǔn)。
為了確定降雨量和歷時(shí)倒數(shù)的邊緣分布,采用極大似然法估計(jì)降雨量與歷時(shí)倒數(shù)的邊緣分布參數(shù),并利用AIC和BIC對(duì)十三陵流域降雨量與歷時(shí)倒數(shù)的16種邊緣分布進(jìn)行擬合度檢驗(yàn)。
十三陵流域的歷時(shí)倒數(shù)和降雨量16種邊緣分布函數(shù)的AIC和BIC檢驗(yàn)結(jié)果見表2—3。由表2可得,歷時(shí)倒數(shù)邊緣分布函數(shù)檢驗(yàn)中g(shù)eneralized pareto分布的AIC值和BIC值最小,由表3可得,降雨量邊緣分布函數(shù)檢驗(yàn)中l(wèi)ognormal分布的AIC值和BIC值最小,即generalized pareto分布和lognormal分布分別作為十三陵流域歷時(shí)倒數(shù)和降雨量最優(yōu)的邊緣分布。繪制優(yōu)選的擬合分布與經(jīng)驗(yàn)分布的匹配程度的比較圖(圖1a、2a)和變量邊緣分布經(jīng)驗(yàn)頻率與理論頻率關(guān)系圖(圖1b、2b),如圖1a、2a所示,點(diǎn)據(jù)均分布在擬合線附近,同時(shí),由圖1b、2b可得,樣本分位數(shù)點(diǎn)和相應(yīng)的理論上擬合的分位數(shù)接近。兩者都直觀得證明擬合結(jié)果較好。
表2 歷時(shí)倒數(shù)的邊緣分布AIC和BIC檢驗(yàn)結(jié)果
表3 降雨量的邊緣分布AIC和BIC檢驗(yàn)結(jié)果
圖1 歷時(shí)倒數(shù)的邊緣分布擬合圖
在優(yōu)選出降雨量和歷時(shí)倒數(shù)的最優(yōu)邊緣分布后,運(yùn)用Copula函數(shù)構(gòu)建降雨量和歷時(shí)倒數(shù)的聯(lián)合分布。本文運(yùn)用極大似然法對(duì)25種Copula函數(shù)進(jìn)行參數(shù)估算,通過(guò)極大似然估計(jì)方法找出使估計(jì)的聯(lián)合概率與相關(guān)的經(jīng)驗(yàn)聯(lián)合概率最匹配的Copula參數(shù)。然后根據(jù)ML和AIC、BIC比較25種Copula函數(shù)的擬合優(yōu)度,從而確定十三陵流域的Copula函數(shù)類型。
十三陵水庫(kù)流域的降雨量和歷時(shí)倒數(shù)的檢驗(yàn)結(jié)果見表3,Roch-Alegre的ML值最大為1150.2,AIC和BIC值最小分別為-2296.4和-2289.0,即ML、AIC、BIC分析法3種評(píng)價(jià)方法均得出一致的結(jié)果,因此,采用Roch-Alegre函數(shù)描述十三陵流域降雨量和歷時(shí)倒數(shù)的聯(lián)合分布。Roch-Alegre函數(shù)參數(shù)值θ1為0.0001,為1.0973,十三陵流域降雨量與歷時(shí)倒數(shù)聯(lián)合分布函數(shù)表達(dá)式如下: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函數(shù)檢驗(yàn)結(jié)果
本文從二元Copula族中選擇最佳的函數(shù)Roch-Alegre來(lái)建立降雨量與降雨歷時(shí)倒數(shù)的依賴結(jié)構(gòu),將只考慮降雨量的單變量的設(shè)計(jì)頻率概念擴(kuò)展到同時(shí)考慮降雨量與降雨歷時(shí)的多變量。這種方法對(duì)于設(shè)計(jì)降雨量大而降雨歷時(shí)短這種極其危險(xiǎn)的降雨具有重要的現(xiàn)實(shí)意義。圖3為Roch-Alegre函數(shù)的聯(lián)合分布圖,圖中的等高線對(duì)應(yīng)著相應(yīng)設(shè)計(jì)頻率的降雨,同一個(gè)等高線上的點(diǎn)構(gòu)成降雨量大而降雨歷時(shí)短二元風(fēng)險(xiǎn)場(chǎng)景。比如,聯(lián)合分布Copula值為0.9的等高線對(duì)應(yīng)著無(wú)數(shù)個(gè)設(shè)計(jì)頻率為0.1的降雨量與降雨歷時(shí)組合。降雨量與降雨歷時(shí)的無(wú)限組合會(huì)在統(tǒng)計(jì)上是相似的,而它們的影響可能截然不同。如圖3中所示的十三陵水庫(kù)的降雨,同一設(shè)計(jì)頻率降雨下的降雨量與降雨歷時(shí)具有無(wú)限種組合,在曲線的右側(cè)對(duì)應(yīng)著降雨歷時(shí)短,降雨量小的場(chǎng)景,而在曲線的上側(cè)對(duì)應(yīng)著降雨歷時(shí)短,降雨量大的場(chǎng)景。選擇哪一種組合是亟需解決的問(wèn)題。從無(wú)限的降雨量與降雨歷時(shí)組合中進(jìn)行選擇的一種直觀方法是根據(jù)它們相關(guān)的Copula密度值為它們分配權(quán)重。Roch-Alegre函數(shù)的概率密度值可以衡量這無(wú)限種組合的權(quán)重,因此本文采用密度函數(shù)極大值處的組合為設(shè)計(jì)頻率下的降雨量與降雨歷時(shí)組合,即選擇一種最可能的事件。
圖3 十三陵水庫(kù)流域降雨量和降雨歷時(shí)倒數(shù)的Roch-Alegre函數(shù)聯(lián)合分布圖
通過(guò)對(duì)Roch-Alegre密度函數(shù)求導(dǎo)找到極值點(diǎn),即設(shè)計(jì)頻率下最可能降雨量與降雨歷時(shí)組合場(chǎng)景。設(shè)計(jì)頻率0.1條件下的最佳組合如圖4所示,降雨量為57mm,降雨歷時(shí)倒數(shù)為1.3(1/h)。而不考慮降雨量與降雨歷時(shí)依賴關(guān)系各自得到設(shè)計(jì)頻率為0.1條件下的降雨量為29mm,降雨歷時(shí)倒數(shù)為1.2mm。可以看出十三陵水庫(kù)控制流域容易出現(xiàn)短歷時(shí)大降雨量,即降雨強(qiáng)度大的情況。而這種情況會(huì)給防洪帶來(lái)很大的風(fēng)險(xiǎn),所以在設(shè)計(jì)降雨時(shí)就需要提前考慮這種風(fēng)險(xiǎn)疊加。
圖4 設(shè)計(jì)頻率0.1條件下降雨量和歷時(shí)倒數(shù)的最佳組合
綜上所述,運(yùn)用本文提出的考慮降雨量與歷時(shí)依賴關(guān)系的方法可以選擇具有聯(lián)合頻率的危險(xiǎn)組合事件,可以用于危險(xiǎn)降雨的設(shè)計(jì)與評(píng)估,在十三陵的防洪中具有重要的現(xiàn)實(shí)意義。
(1)通過(guò)擬合優(yōu)度評(píng)價(jià),得出北京市十三陵水庫(kù)流域降雨量采用lognormal分布擬合效果較好,而降雨歷時(shí)的倒數(shù)采用generalized pareto分布擬合效果較好,建議選擇Roch-Alegre函數(shù)描述降雨量和歷時(shí)倒數(shù)的聯(lián)合分布。
(2)采用雙變量設(shè)計(jì)降雨模型與分開考慮降雨量與歷時(shí)的傳統(tǒng)方法相比,得到的設(shè)計(jì)降雨量偏大,降雨歷時(shí)偏小。
(3)本文所提的考慮降雨量與歷時(shí)依賴關(guān)系設(shè)計(jì)降雨方法,有助于應(yīng)對(duì)降雨量大而降雨歷時(shí)短的二元風(fēng)險(xiǎn),在水庫(kù)防洪中具有現(xiàn)實(shí)意義,以期為推求高強(qiáng)度降雨的設(shè)計(jì)降雨提供技術(shù)支持。本文未探討短歷時(shí)強(qiáng)降雨所形成的洪水過(guò)程,有待繼續(xù)研究。