基于非對(duì)稱GH-Copula函數(shù)推求珠海市設(shè)計(jì)暴雨過(guò)程線*

陳子燊，李鴻皓

（1. 中山大學(xué)新華學(xué)院，廣東廣州510520；2. 中山大學(xué)地理科學(xué)與規(guī)劃學(xué)院，廣東廣州510275；3. 珠海市公共氣象服務(wù)中心，廣東珠海51900）

設(shè)計(jì)典型暴雨過(guò)程線是制定排水防洪系統(tǒng)設(shè)計(jì)時(shí)需要考慮的關(guān)鍵問(wèn)題，可為市政排水系統(tǒng)的規(guī)劃和管理、與城市雨洪分析計(jì)算提供科學(xué)依據(jù)。一些研究人員從城市應(yīng)對(duì)極端天氣事件與防災(zāi)減災(zāi)的風(fēng)險(xiǎn)管理角度對(duì)雨型作了探索。蔣明［1］指出，雨型是描述降雨過(guò)程和降雨強(qiáng)度在時(shí)間尺度上的分配過(guò)程，是徑流過(guò)程計(jì)算的基礎(chǔ)。成丹等［2］把設(shè)計(jì)雨型作為制定排水防澇系統(tǒng)設(shè)計(jì)時(shí)的重要因素，應(yīng)用于城市市政排水系統(tǒng)的規(guī)劃和管理及排水分析，為城市流域雨洪調(diào)度計(jì)算提供科學(xué)依據(jù)。葉姍姍等［3］選取宿遷市實(shí)測(cè)的主副型雨峰偏后的暴雨雨型，對(duì)其降雨過(guò)程進(jìn)行同頻率分時(shí)段縮放，采用Copula 函數(shù)的風(fēng)險(xiǎn)聯(lián)合概率模型分析了不同的兩時(shí)段之間出現(xiàn)的暴雨風(fēng)險(xiǎn)。楊星等［4］利用深圳雨量站34 年實(shí)測(cè)逐時(shí)降雨資料，對(duì)比了不同典型暴雨設(shè)計(jì)雨型研究方面的差異，按構(gòu)建的兩變量Copula 推求了深圳市不同重現(xiàn)期雨型的風(fēng)險(xiǎn)率和典型暴雨的特征。陳子燊等［5］使用珠海市32 年不同歷時(shí)暴雨組合推算排水排澇兩級(jí)標(biāo)準(zhǔn)銜接的設(shè)計(jì)暴雨水平，分析各歷時(shí)暴雨組合的遭遇概率，推算各組合的設(shè)計(jì)暴雨值。

Grimaldi & Serincesco［6］、Ganguli & Reddy［7］、Zhang & Singh［8］、陳子燊等［9］分別采用非對(duì)稱的阿基米德Copulas 構(gòu)建了不對(duì)稱三變量洪水要素聯(lián)合分布模型推算設(shè)計(jì)洪水，嘗試應(yīng)用于洪水風(fēng)險(xiǎn)規(guī)劃管理。肖義等［10］、李天元等［11］則分別采用兩變量和三變量的Copula 函數(shù)建立了聯(lián)合分布的設(shè)計(jì)洪水過(guò)程線的推求方法，為城市設(shè)計(jì)暴雨過(guò)程線提供了一種新思路。

華南雨水非常豐沛，汛期降雨可分為以鋒面雨為主的前汛期與臺(tái)風(fēng)雨為主的后汛期。由于城市化進(jìn)程加速，城市建設(shè)管網(wǎng)等排水系統(tǒng)滯后于城市發(fā)展，設(shè)計(jì)標(biāo)準(zhǔn)偏低，導(dǎo)致城市排水系統(tǒng)的排洪能力低下。近年來(lái)一些城市內(nèi)澇問(wèn)題日顯突出，衍生的災(zāi)害對(duì)城市經(jīng)濟(jì)社會(huì)發(fā)展帶來(lái)了重大的負(fù)面影響。本研究以珠海市近數(shù)十年來(lái)的逐時(shí)雨量數(shù)據(jù)，采用非對(duì)稱阿基米德Copulas(Asym?metric Archimedean Copulas)構(gòu)建典型設(shè)計(jì)暴雨過(guò)程線，可為城市防洪排澇工程設(shè)計(jì)和防災(zāi)減災(zāi)的風(fēng)險(xiǎn)管理提出新思路。

1 三變量Copula函數(shù)

1.1 Copula函數(shù)的定義

設(shè)隨機(jī)變量Xi(i= 1，2，…，n)的邊緣分布函數(shù)分別為FXi(xi)=P(Xi≤xi)，其中n為隨機(jī)變量的個(gè)數(shù)，xi為隨機(jī)變量Xi的值。依Sklar 理論，對(duì)于連續(xù)分布函數(shù)FXi(xi)，存在唯一的聯(lián)合分布函數(shù)［12］：

利用Copula 函數(shù)構(gòu)造聯(lián)合概率分布，使得變量的所有信息都存在于邊緣分布函數(shù)里，不會(huì)在轉(zhuǎn)換過(guò)程中產(chǎn)生信息失真。因此，Copula函數(shù)理論是構(gòu)建多變量水文聯(lián)合概率分布的理想工具。

1.2 三變量Archimedean Copula

三變量對(duì)稱的Archimedean Copula 單參數(shù)形式［13］為：

式 中，u1，u2，u3為 邊 緣 分 布，u1=FX1(x1)，u2=FX2(x2)，u3=FX3(x3)，u1，u2，u3∈［0，1］；φθ為Ar?chimedean copula 的生成元；θ為參數(shù)；φ?1θ為φθ的反函數(shù)。對(duì)稱的GH?Copula為：

三變量非對(duì)稱Archimedean Copula兩參數(shù)形式為［13］：

式中，符號(hào)“°”表示函數(shù)組合。非對(duì)稱的GH?Copula為：

2.3 三變量聯(lián)合重現(xiàn)期

3 實(shí)例研究

3.1 基本數(shù)據(jù)

本文根據(jù)珠海市氣象站1961～2018 年逐時(shí)降水記錄數(shù)據(jù)，在提取歷年最大24 h 雨量基礎(chǔ)上，再分別提取雨峰和連續(xù)最大6 h 雨量數(shù)據(jù)，構(gòu)成最大1 h 雨量(R1)、最大6 h 雨量(R6)和最大24 h雨量(R24)作為分析樣本。其中，珠海最大24 h出現(xiàn)在2000 年4 月13 日19 時(shí)～4 月14 日18 時(shí)，累計(jì)降雨630.3 mm，出現(xiàn)在第13 h 的雨峰，雨量為89.3 mm。次最大24 h雨量出現(xiàn)在1994年7月22日7 時(shí)～7 月22 日6 時(shí)，累計(jì)降雨572.6 mm，雨峰靠前，出現(xiàn)在第4 h，雨量102 mm。二者分別作為前汛期和后汛期的代表性極端降水過(guò)程，以此構(gòu)建2000 年和1994 年3 個(gè)時(shí)段雨量聯(lián)合分布的兩場(chǎng)典型設(shè)計(jì)暴雨過(guò)程線。

3.2 邊緣分布與聯(lián)合分布

三個(gè)樣本的邊緣分布均采用廣義極值(GEV)分布函數(shù)：

式中，ξ，β，μ分別為形態(tài)參數(shù)、尺度參數(shù)和位置參數(shù)。

采用廣義極值分布(GEV)對(duì)R1、R6和R24樣本加以擬合。參數(shù)估計(jì)使用線性矩(L?矩)方法。經(jīng)驗(yàn)頻率分布使用Gringorten 公式計(jì)算。擬合結(jié)果采用均方根誤差(RMSE)和概率點(diǎn)據(jù)相關(guān)系數(shù)(PPCC)檢驗(yàn)其擬合優(yōu)度。擬合結(jié)果(表1)表明，GEV 分布適用于3 個(gè)極端降雨尺度R1、R6和R24的邊緣分布。

計(jì)算的R1、R6和R24兩兩間的Kendall 相關(guān)系數(shù)表明不同歷時(shí)暴雨間存在較高的相關(guān)性，其中R1和R6最大，τ=0.576；R1和R24之間的相關(guān)系數(shù)為0.303；R6和R24相關(guān)系數(shù)為0.543。采用Kendall 相關(guān)系數(shù)與Copula 參數(shù)的關(guān)系式，構(gòu)建非對(duì)稱三變量GH?Copula［10］。為了對(duì)比，通過(guò)MLM 法計(jì)算了三變量對(duì)稱GH?Copula 的參數(shù)θ。采用Akaike 信息準(zhǔn)則(AIC)和最小OLS 準(zhǔn)則驗(yàn)證理論聯(lián)合分布函數(shù)與經(jīng)驗(yàn)聯(lián)合分布函數(shù)的擬合程度，結(jié)果見(jiàn)表2?？梢?jiàn)以二維GH?Copula 的三維非對(duì)稱形式的GH?Copula 的OLS 和AIC 值最小，表明非對(duì)稱Copula 略優(yōu)。Nelson［12］和Salvadori 等［13］證明當(dāng)且僅當(dāng)邊緣分布和Copula 函數(shù)均為極值分布時(shí)，構(gòu)造的聯(lián)合分布才是極值分布，而GH?Copula 是Archimedean copula 函數(shù)族中的唯一多變量極值Copula 函數(shù)，適用于極端事件的頻率分析?？紤]到R1、R6和R24之間的相關(guān)性存在較大差別，因此選用非對(duì)稱GH?Copula構(gòu)建珠海雨量站歷年最大24 h暴雨量三個(gè)歷時(shí)R1、R6和R24之間的聯(lián)合分布：

表1 三變量暴雨樣本的GEV分布參數(shù)與優(yōu)度檢驗(yàn)值Table 1 The parameters of the marginal distribution and the values of goodness of fit test

圖1 降雨尺度R1、R6和R24的GEV分布Fig.1 GEV distribution of R1,R6 and R24rainfall scales

表2 非對(duì)稱和對(duì)稱三維GH?Copula參數(shù)估計(jì)及擬合優(yōu)度評(píng)價(jià)Table 2 Results of parameter estimation and goodness of fit tests for three dimensional asymmetric and symmetric GH?Copulas

3.4 典型設(shè)計(jì)暴雨過(guò)程線

典型暴雨過(guò)程線的選擇采用以下原則：1)分別選擇峰高量大具有一定代表性的實(shí)測(cè)暴雨過(guò)程線暴雨過(guò)程；2)設(shè)計(jì)暴雨過(guò)程線采用同頻率放大法，以降水主峰對(duì)城市內(nèi)澇形成為首要影響因子，選定時(shí)段為1 h 的設(shè)計(jì)雨峰作為設(shè)計(jì)標(biāo)準(zhǔn)，放大出來(lái)的過(guò)程線形狀能與原來(lái)的典型過(guò)程一致。

珠海前、后汛期典型年的24 h 暴雨過(guò)程線如圖2 顯示，1994 年最大24 h 暴雨過(guò)程存在主副峰，主峰靠前。2000 年最大24 h 暴雨過(guò)程同樣為主副峰雨型，主峰大致位于最大24 h 暴雨的中間時(shí)刻。兩個(gè)典型暴雨R1、R6、R243 個(gè)時(shí)段最大降水量及重現(xiàn)期見(jiàn)表3。

圖2 兩個(gè)典型年暴雨過(guò)程Fig.2 Rainstormhydrograph of two typical year

從表3 可見(jiàn)，1994 年和2000 年的R1?R6?R24組合雨量的“或”聯(lián)合重現(xiàn)期小于單一時(shí)段雨量重現(xiàn)期，此說(shuō)明考慮多時(shí)段組合條件下單一時(shí)段雨量即可致災(zāi)的可能性最高，同時(shí)出現(xiàn)三時(shí)段組合雨量的“且”聯(lián)合重現(xiàn)期明顯超過(guò)“或”聯(lián)合重現(xiàn)期，屬于可能性小的事件。表4為不同時(shí)段雨量組合的“或”聯(lián)合重現(xiàn)期，同樣可見(jiàn)同頻率下R1?R6?R24三時(shí)段組合雨量的“或”聯(lián)合重現(xiàn)期小于兩時(shí)段組合雨量的“或”聯(lián)合重現(xiàn)期。因此，如果以三時(shí)段雨量組合的“或”聯(lián)合重現(xiàn)期更適用于應(yīng)對(duì)珠海市的內(nèi)澇風(fēng)險(xiǎn)。

3.5 同頻率法推求設(shè)計(jì)暴雨過(guò)程線

按照同頻率放大法的思路，假定3個(gè)時(shí)段降雨變量R1、R6、R24同頻率，令u1=u2=u3，可得到基于某一聯(lián)合重現(xiàn)期Tu1，u2，u3的頻率組合(u1，u2，u3)。根據(jù)此組合，按照各變量的邊緣分布函數(shù)反推可得到三變量不同時(shí)段雨量的聯(lián)合設(shè)計(jì)值組合(r1，r2，r3)，再以此設(shè)計(jì)值組合放大典型暴雨過(guò)程，即得到基于三變量聯(lián)合分布的設(shè)計(jì)暴雨過(guò)程線。非對(duì)稱GH?Copula 推算的R1、R6、R24三變量同頻率分布聯(lián)合設(shè)計(jì)值為［10］：

表3 兩個(gè)典型年R1、R6、R24最大降水量及重現(xiàn)期Table 3 Maximum precipitation and return period of R1，R6，R24 of two typical years

表4 不同時(shí)段雨量組合的“或”聯(lián)合重現(xiàn)期Table 4 ‘OR’joint return periods of rainstorm in different time interval

同 理 可 分 別 推 算 兩 變 量(u1，u2) 的Tu1，u2和(u1，u3)的Tu1，u3的同頻率聯(lián)合設(shè)計(jì)值。

選取當(dāng)?shù)貙?shí)測(cè)的典型暴雨，對(duì)其降雨過(guò)程進(jìn)行同頻率分時(shí)段縮放。放大系數(shù)公式K=X設(shè)計(jì)/X典型。式中X設(shè)計(jì)為不同重現(xiàn)期的設(shè)計(jì)降雨量；X典型為典型暴雨降雨量。以雨峰同頻率放大方法求重現(xiàn)期200 年(P=0.05%)R1?R6?R24三時(shí)段聯(lián)合分布的設(shè)計(jì)暴雨過(guò)程線。為了比較，另推求了R1?R6和R1?R24兩變量聯(lián)合分布以及以雨峰R1單變量同頻率放大的設(shè)計(jì)暴雨過(guò)程線(圖3)。多變量同頻率設(shè)計(jì)值計(jì)算結(jié)果見(jiàn)表5。

多變量方法與單變量方法所推求的200年一遇設(shè)計(jì)暴雨過(guò)程線的比較顯示，考慮了變量之間的相關(guān)性的多變量聯(lián)合重現(xiàn)期設(shè)計(jì)值大于單變量同頻率設(shè)計(jì)值，即多變量方法放大的過(guò)程線對(duì)于城市防雨洪設(shè)計(jì)偏安全或保險(xiǎn)。以多變量聯(lián)合重現(xiàn)期作為城市排水防澇標(biāo)準(zhǔn)更安全，傳統(tǒng)的單變量同頻率方法則達(dá)不到標(biāo)準(zhǔn)。因此，得到的典型設(shè)計(jì)暴雨過(guò)程線更符合城市水文的內(nèi)在規(guī)律和排水防澇工程設(shè)計(jì)的要求，由此設(shè)計(jì)的暴雨過(guò)程線更適用于應(yīng)對(duì)城市內(nèi)澇風(fēng)險(xiǎn)。

圖3 典型暴雨三時(shí)段組合、兩時(shí)段組合和單一時(shí)段方法推求的200年一遇設(shè)計(jì)暴雨過(guò)程線的對(duì)比Fig.3 Comparison among the design rainstorm hydrograph with return period 200 years of a typical rainstorm to three and two durations

表5 單變量暴雨設(shè)計(jì)值和多變量同頻率聯(lián)合設(shè)計(jì)值Table 5 Univariate and multivariate storm design value combined with same frequency

4 結(jié) 論

本文將不同歷時(shí)雨量之間具有相關(guān)關(guān)系的暴雨過(guò)程簡(jiǎn)化為雨峰、6 h 雨量和24 h 雨量三變量聯(lián)合分布，采用非對(duì)稱GH?copula 構(gòu)建珠海市典型暴雨過(guò)程線，并與由兩個(gè)時(shí)段和由單一雨峰的同頻率設(shè)計(jì)暴雨過(guò)程線方法進(jìn)行了比較。研究結(jié)果有以下結(jié)論：

1)多變量聯(lián)合分布模型合理地拓展了Copula函數(shù)在城市水文工程中的實(shí)際應(yīng)用。1994 年和2000 年兩個(gè)典型年的R1?R6?R24三個(gè)時(shí)段雨量聯(lián)合分布的“或”聯(lián)合重現(xiàn)期都小于單一時(shí)段雨量重現(xiàn)期，危險(xiǎn)率最大。三個(gè)時(shí)段雨量聯(lián)合分布的“或”聯(lián)合重現(xiàn)期更適用于應(yīng)對(duì)城市雨洪風(fēng)險(xiǎn)。

2)采用三個(gè)歷時(shí)雨量聯(lián)合分布推求的珠海市設(shè)計(jì)暴雨值大于兩個(gè)時(shí)段和單一時(shí)段設(shè)計(jì)暴雨值，由此放大的設(shè)計(jì)暴雨過(guò)程線，整體效果相對(duì)最優(yōu)，對(duì)設(shè)計(jì)暴雨過(guò)程線的研究方法提供了新思路。按同頻率放大的典型設(shè)計(jì)暴雨過(guò)程線可作為前汛期和后汛期排水防澇工程設(shè)計(jì)的參考依據(jù)。

3)與1994年主峰靠前的暴雨過(guò)程比較，主峰居中的2000 年的暴雨過(guò)程危險(xiǎn)率更大，對(duì)珠海市內(nèi)澇具有更大的威脅，構(gòu)建的典型設(shè)計(jì)暴雨過(guò)程線更具代表性。