苗正偉,丁志宏,路 梅,楊學智
(1.河北水利電力學院 水利工程學院,河北 滄州 061001; 2.中水北方勘測設計研究有限責任公司,天津 300222)
洪水頻率是水文統(tǒng)計的重點內容,傳統(tǒng)的單變量頻率分析或者研究洪峰頻率,或者研究洪量頻率,或者假定二者同頻率,這既違背了洪水是多變量隨機事件的科學事實,也忽略了洪水特征變量之間普遍存在的較強相關性。有不少學者基于多變量頻率分析方法來描述洪水、干旱等水文事件,但傳統(tǒng)多變量聯(lián)合分布對邊緣分布的苛刻要求與實際水文現象的高度復雜性格格不入,極大地限制了傳統(tǒng)多變量聯(lián)合分布的應用[1]。近年來,基于Copula函數的聯(lián)合分布在水文領域越來越受關注,Sklar[2]在1959年將Copula概念引入統(tǒng)計領域,繼而由Genest等[3]、Harry等[4]進一步發(fā)展、完善了Copula相關理論,隨后,Nelsen[5]在其專著中對Copula函數理論與成果進行了總結,Marius Hofert等[6]進一步總結了Copula函數的性質、應用及其R語言的實現。眾多研究結果表明,Copula函數很適于描述水文特征變量間的相關性,包括非線性相關關系,且對邊緣分布沒有限制,已被廣泛應用于暴雨、洪水、干旱等水文事件的多變量聯(lián)合分布研究[7-11]。
重現期是水文設計的核心概念之一。在統(tǒng)計意義上,重現期是發(fā)生指定事件的時間間隔的平均值。對單變量事件而言,重現期和設計值是一一對應的。但對于多變量事件,因變量間存在無窮多種組合,其重現期的定義以及設計值的計算都要復雜得多,這已成為近年來水文領域的熱點問題[12-14]。定義雙變量洪水重現期的核心問題是危險域的界定,所謂危險域[15]就是所有危險事件的集合,危險域之外的二維實數域即為安全域,但一般將安全域與危險域的分界線稱為臨界域。給定洪水峰量(峰量即洪峰流量和洪量)組合,若其落在危險域內,則將其判定為超標洪水,屬于危險事件,反之則是安全的或處于臨界狀態(tài)。根據雙變量框架下危險域界定方式的不同,目前對雙變量洪水事件重現期的定義也有很多種,其中,“或”(“OR”)和 “且”(“AND”)重現期應用最廣泛, OR重現期將洪水峰量中至少有一個超過指定閾值的事件集合定義為危險域;而AND重現期將危險域定義為兩變量同時超過指定閾值的事件集合。根據不同的危險域進而計算相應的重現期。但是,這2種重現期對于危險事件的判斷均不盡合理[16-18]。為此,Salvadori 等[17]根據Kendall函數得到Copula等值線,并將峰量聯(lián)合分布概率超過給定閾值的事件集合定義為危險域,進而提出了Kendall 重現期的概念,這使得危險事件的判斷不再有邏輯上的矛盾[19-21]。但該方法給出的安全域卻是無限的,這與事實明顯不符。因此,Salvadori 等[17]于2013年又提出了生存 Kendall重現期的概念,將危險域定義為洪水峰量的聯(lián)合生存概率小于給定閾值的事件集合。生存Kendall重現期不僅對危險域的判定科學合理,而且它給出的安全域也是有界的,與現實相符,是當前雙變量水文事件相對嚴謹的重現期定義方式。
與單變量情況不同,給定的雙變量重現期會對應無窮多種變量組合,由此構成重現期等值線,這意味著,雙變量框架下聯(lián)合洪水設計值并非唯一,而不同的洪水峰量組合將得到不同的洪水過程線,從而對調洪演算的結果產生直接影響,并進一步影響水庫特征水位及水庫規(guī)模,因此,多變量框架下如何科學合理地計算聯(lián)合洪水設計值是至關重要的[22]。目前,計算聯(lián)合洪水設計值的常用方法主要包括同頻率法[23]、極大似然法[18,24]、條件最可能組合法[25-27]等。
當前關于洪水的雙變量設計主要基于OR、AND以及Kendall重現期,較少涉及生存Kendall重現期,而且缺乏洪水設計值計算方法的對比分析。為此,本文以崗南水庫為例,基于Copula函數,采用OR、AND、Kendall及生存Kendall 4種重現期標準以及同頻率、極大似然、條件最可能組合3種洪水設計值計算方法,對洪水的峰量聯(lián)合分布展開研究,以期為水文聯(lián)合設計提供科學參考,為崗南水庫的防洪安全提供理論依據。
位于海河流域子牙河系滹沱河中游的崗南水庫控制流域面積15 900 km2(如圖1所示),總庫容15.71億m3,其中9.71億m3為防洪庫容,主要擔負著下游京九鐵路、華北油田及石家莊等的防洪保護任務,同時承擔向灌區(qū)及石家莊等地區(qū)供水的任務,該水庫對支撐河北省水資源可持續(xù)利用及社會經濟發(fā)展至關重要。本研究所用數據主要是1949—2012年崗南水庫流域年最大洪峰流量和六日洪量,此外還包括一場重現期為3 000 a的古洪水及11場歷史特大洪水,數據來自水文年鑒及海河水利委員會。
圖1 崗南水庫流域Fig.1 Drainage map of Gangnan reservoir
對崗南水庫流域洪峰流量及六日洪量系列,由非連續(xù)樣本的線性矩法初估參數進而進行適線,經驗頻率點據、擬合的統(tǒng)計參數及P-III分布曲線如圖2所示。由此可見,適線效果良好。
圖2 P-III分布擬合結果Fig.2 Fitting results of P-III distribution
表1 Copula函數參數估計和統(tǒng)計檢驗結果Table 1 Estimated parameters of Copula functions and statistical test results
圖3 聯(lián)合觀測值的經驗分布與理論分布比較Fig.3 Comparison between empirical and theoretical distributions of joint observations
給定重現期水平20、50、100、200 a,在同一坐標系中繪制4種重現期標準的等值線,如圖4所示。由此可知,對于不同標準的重現期,OR與Kendall重現期等值線走勢一致;AND與生存Kendall重現期走勢一致。盡管重現期等值線趨勢相似,但它們確定的危險域截然不同,OR、AND重現期等值線上每一點都對應不同的危險域,正因如此,這二者對危險域的判斷存在邏輯上的矛盾,而Kendall、生存Kendall重現期等值線上各點都有一樣的危險域,這使得兩者對危險域的判定更為合理。圖4還顯示: Kendall重現期給出的安全域是無界的,這顯然與實際相悖,與之相反,生存Kendall重現期給出了有界的安全域,更契合工程安全的內在要求,綜上所述,4種重現期標準中,生存Kendall重現期最為科學合理。由圖4還可知,4種給定重現期水平下,不同重現期標準等值線間的關系是一致的,均表現為:OR重現期等值線在最上、AND重現期等值線在最下,Kendall與生存Kendall重現期等值線居中;AND、OR、Kendall 3種重現期等值線無交點;生存Kendall重現期等值線與Kendall重現期等值線有2個對稱的交點,與OR、AND重現期等值線無交點。
圖4 重現期等值線Fig.4 Isolines of return periods
3.4.1 不同計算方法的洪水設計值比較
采用同頻率、極大似然、條件最可能組合3種方法計算各重現期標準的設計值,結果如表2、表3所示。由此可知,對所有給定的重現期水平和重現期標準,無論是洪峰流量還是六日洪量,都有如下規(guī)律成立:同頻率法與極大似然法所計算的設計值差別非常小,條件最可能組合法計算的洪峰流量略大于極大似然法和同頻率法的洪峰流量,而六日洪量略小于極大似然法和同頻率法的六日洪量。不失一般性,以100 a生存Kendall重現期為例,同頻率法與極大似然法求得的洪峰流量設計值僅相差0.334 m3/s,條件最可能組合法比同頻率法的結果大差值為65.62 m3/s,相對差異僅為0.69%。六日洪量設計值的差異與洪峰設計值的差異類似。由此可見,對于洪水設計而言,極大似然法、條件最可能組合法、同頻率法3種不同計算方法對洪水設計值并無明顯影響,其中尤以同頻率法計算最為簡單,因此建議使用同頻率法進行雙變量洪水設計值計算。
3.4.2 不同重現期標準的洪水設計值比較
由表2、表3可知,給定重現期水平,無論何種洪水設計值計算方法,計算所得的洪峰、洪量設計值都遵循以下規(guī)律:OR重現期標準的洪水設計值>生存Kendall重現期洪水設計值>單變量洪水設計值>Kendall重現期洪水設計值>AND重現期洪水設計值;各重現期標準的洪水設計值與單變量洪水設計值的相對差異基本呈現隨重現期水平增大而減小的趨勢,而且AND重現期、OR重現期洪水設計值的相對差異明顯大于Kendall和生存Kendall重現期洪水設計值的相對差異。因此,以OR重現期標準設計工程將偏于保守,工程規(guī)模偏大,投資高;而以AND重現期標準進行設計,則工程規(guī)模偏小,投資少,但安全性降低。由此可見,AND和OR是工程設計的兩種極端情形,它們在工程規(guī)模與投資、安全性方面是失衡的,相較而言,根據生存Kendall或Kendall重現期進行工程設計更為科學合理。又因為生存Kendall重現期的定義更加嚴謹,因此,生存Kendall重現期是兼顧設計洪水計算的合理性以及工程的經濟性、安全性的一個較好選擇,更利于洪水的風險管理。
由表2、表3可知:本文計算的不同重現期標準的洪水設計值與單變量洪水設計值相對差異的絕對值絕大部分在9%以內,差別較小,尤其是當重現期水平超過1 000 a時,各重現期標準、各計算方法所得的聯(lián)合洪水設計值與單變量洪水設計值的相對差異都在3%以內。當將Copula參數減小為3.000 0,聯(lián)合洪水設計值與單變量洪水設計值相對差異見表4和表5(聯(lián)合洪水設計值計算方法為同頻率法)。由此可見,當Copula參數減小時,相對單變量洪水設計值而言,Kendall、AND、OR 3種重現期標準聯(lián)合洪水設計值的差異都明顯增大,增幅在100%左右,而生存Kendall重現期的聯(lián)合洪水設計值差異也有小幅增加。根據二維對稱Gumbel Copula參數與Kendall相關系數的關系進一步分析可知,當Copula參數減小時,Kendall相關系數也是減小的,即洪水的峰量相關性減弱,而本研究采用的Copula參數為5.286 2,相應的Kendall相關系數高達0.810 8,當Copula參數為3.000 0時,Kendall相關系數減小為0.666 7,這說明:洪水峰、量相關性越弱,基于Copula的雙變量聯(lián)合洪水設計值與單變量洪水設計值的差異越明顯。
表2 洪峰流量設計值比較Table 2 Comparison of design peak discharge
表3 六日洪量設計值比較Table 3 Comparison of design 6d flood volume
表4 Copula參數對洪峰流量設計值的影響Table 4 Influence of parameter of Copula on design peak discharge
表5 Copula參數對六日洪量設計值的影響Table 5 Influence of parameter of Copula on design 6 d flood volume
圖4顯示,生存Kendall重現期等值線與Kendall重現期等值線圍成的閉合區(qū)域在Kendall重現期概念下,被判定為危險域,而在生存Kendall概念下被判定為安全域,相較于洪量很大而洪峰很小、或者洪峰很大而洪量很小的情形,該交集區(qū)域內的峰量關系更契合現實洪水峰量的相關性,即這類洪水事件更可能在現實中發(fā)生。因此,該區(qū)域究竟被判定為危險域還是安全域,對防洪工程的規(guī)模具有直接影響。因為生存Kendall重現期將該閉合區(qū)域判定為安全域,所以生存Kendall重現期確定的洪水設計值比Kendall重現期的洪水設計值更大(圖4中的同頻率設計值及表2、表3的結果都證明了該種情況)。因此,基于生存Kendall重現期確定的工程規(guī)模也比Kendall重現期的更大,即以生存Kendall重現期進行工程設計比Kendall重現期進行工程設計更安全。
本文以崗南水庫為例,采用4種重現期標準、3種聯(lián)合洪水設計值計算方法對雙變量洪水設計值進行了研究,主要結論如下。
(1)AND和OR重現期在危險域和安全域的識別上存在局限性,相對而言,Kendall重現期更合理,但其安全域是無界的,這與實際不符。生存Kendall重現期則界定了有界的安全域,使得重現期的概念在邏輯上更科學合理。
(2)同頻率法、極大似然法、條件最可能組合法3種洪水設計值計算方法所得結果相差不大,尤其是同頻率法與極大似然法,二者結果幾乎相等,從簡單實用的角度出發(fā),推薦采用同頻率法計算聯(lián)合洪水設計值。
(3)基于Copula計算給定重現期水平的洪水設計值,既受重現期標準的影響,也受設計值計算方法的影響。但不同設計值計算方法對結果的影響較小,因此,設計結果主要取決于重現期標準,且遵循以下規(guī)律:OR重現期標準的洪水設計值>生存Kendall重現期洪水設計值>單變量洪水設計值>Kendall重現期洪水設計值>AND重現期洪水設計值。從計算過程科學合理、結果安全可靠且兼顧經濟性的角度來說,推薦采用生存Kendall重現期進行雙變量洪水設計。
(4)聯(lián)合洪水設計值與單變量洪水設計值的差異受變量間相關性的影響較大,且變量相關性越弱,差異越大。