氣候變化對巖溶區流域極端徑流頻率分析的影響研究

王大洋，王大剛

（中山大學地理科學與規劃學院，廣州510275）

當今環境下，氣候變化正時刻擾動著地表流域的水文循環過程，使許多流域的水文系統發生了明顯的改變。極端氣候的頻繁發生導致極端水文事件頻發，量級也明顯增大，這其中最為突出的當屬降水和洪水［1，2］。傳統的洪水頻率分析多是基于徑流序列滿足一致性的前提，而氣候變化和人類活動的影響，很可能會使原始洪水序列的一致性遭到破壞。因此，研究氣候變化條件下的洪水頻率分析是非常重要且有現實意義的［3-5］。此外，巖溶區作為中國西南地區常見的地貌類型，其水文循環過程本身就較普通流域復雜，因此，有必要針對巖溶區流域開展氣候變化對其極端徑流頻率分析的影響研究。

1 研究方法

1.1 趨勢分析和突變檢驗

研究分別采用線性趨勢分析和Mann-Kendall 檢驗方法對研究流域的徑流序列進行趨勢分析和突變點確定，從而將徑流序列劃分為兩個時段。突變前時段被認為是未受到氣候變化影響的時段，稱為天然時段；突變之后時段被認為是受到氣候變化影響的時段，稱為氣候變化影響時段。由于線性趨勢分析和Mann-Kendall 檢驗方法較為常用，故本文不再贅述相關理論，具體可參考相關文獻［6-7］。

1.2 徑流頻率分析模型優選

針對極端徑流情況，研究采用廣義極值分布（GEV）和廣義帕累托（GPD）分布。選擇兩者的原因有兩點：①與皮爾遜III（P-III）型曲線在中國的廣泛應用相比，GEV分布和GPD分布在中國的研究起步相對較晚，但近些年的關注度在持續上升；②多數分布都是以兩者的統一形式為基礎進行衍生，如Gumbel分布、Fréchet 分布和Weibull 分布是GEV 分布的形狀參數選取不同符號（正、零、負）的具體表現。Pareto 分布、指數分布和Beta分布分別是GPD 形狀參數為（正、零、負）的具體表現［5］。GEV分布和GPD分布模型的表達式如下。

GEV的概率分布表達式為：

GPD的概率分布表達式為：

式中：μ、σ、ξ分別為模型的位置參數，尺度參數和形狀參數，滿足μ∈R，σ＞0，ξ∈R和1+ξ(x-μ)/σ＞0的條件。

1.3 樣本抽取

樣本抽選是頻率分析的基礎，針對GEV和GPD兩種不同的分布模型，其抽樣方法也是有所差異的。對于GEV，通常選擇序列中的“年最大值（AM）”作為抽樣依據；而對于GPD，通常選擇超過某個定量值（POT）的所有數據來建模，適用于超門限值等系列的頻率分析。本研究中，參考我國學者進行POT 模型分析的經驗，以樣本年限為基數，以年均超過定量值發生次數2倍為原則進行樣本抽取，此外，對于POT 抽樣方法，還應當注意抽取樣本的獨立性，從而避免因樣本依賴導致的抽樣誤差［5］。

1.4 參數估計

由于極大似然估計法（MLE）結果能夠在總體上較好地反映樣本的統計信息，具有良好的統計性質［8，9］。因此本研究采用極大似然估計法進行分布模型的參數估計。用MLE 對GEV 和GPD進行參數估計的似然函數如下。

對于GEV，計算極大似然估計參數需要求解方程組：

對于GPD，其極大似然方程為：

式中：xi表示樣本值；ε表示門限閾值。

為了便于求解，通常將其轉化為對數形式的方程：

根據牛頓迭代法，可以計算得到相關參數。

1.5 模型優選檢驗

模型建立后，通過對觀測樣本進行擬合和參數估計，可得到各模型的擬合結果。結合模型檢驗方法可以優選出最合適的線型。本研究采用Kolmogorov-Smirnov 檢驗［10］（簡稱K-S 檢驗），PP概率圖檢驗，QQ 圖檢驗3種檢驗方法，根據GEV 和GPD模型的擬合結果進行檢驗，從而確定最優分布模型。

（1）Kolmogorov-Smirnov 檢驗。該方法以經驗分布函數和理論分布函數之差為基礎，進行模型擬合效果的檢驗，其檢驗統計量為：

式中：yi是經驗頻率去線上的點；F(yi)是累積分布函數，N為樣本大小。

（2）PP 圖和QQ 圖。PP 累積分布概率（CDF）圖是根據變量的累積概率對應所指定的理論分布概率來繪制散點圖，它可以較為直觀地反映樣本數據與擬合模型分布的一致性。當樣本數據符合某一分布時，樣本點會較為接近45°對角線。QQ 圖則和PP圖比較相似，只是采用變量數據分布的分位數與所指定分布的分位數的關系進行比較。兩種方法的基本原理如下。

若樣本為{x1，x2，…，xN}，若經驗分布函數為(x) ={x1，x2，…，xN}，模型估計的擬合分布函數為(x)，則可以通過(x)和(x)的差異反映模型擬合效果。若用概率進行表示（PP圖），則為：

若用樣本點的分位數對比，即QQ圖：

2 研究區域概況

澄碧河發源于廣西壯族自治區西部的百色市，屬于珠江流域的西江水系，河流長度為127 km，流域集雨面積2 087 km2。高程由西北向東南逐漸走低，流域內分布著較多的典型石灰巖巖溶峰林，伏流河，落水洞較為常見，是世界上最大的連續巖溶地區之一。流域多年平均降水量為1 560 mm，因地形原因和巖溶發育較廣，汛期降水多集中且易引發洪澇災害，造成較為嚴重的社會經濟損失和負擔。

流域內部設有4個水文站，其中部分站點由于遷移、改測和撤銷等原因使水文序列時間較短，且不能滿足一致性要求。研究選取平塘站的1963-2017年的實測徑流數據，該站點以上集雨面積占據整個流域的70%左右，且觀測資料較為完整，因此可以作為研究數據。

3 研究結果及分析

3.1 趨勢分析及突變檢測

氣候變化影響流域水文循環，水文循環中水分經過降水產流匯流過程，最終在河道的徑流環節得以體現。研究選取澄碧河水庫流域的年最大1日徑流量、年最大3日徑流量、年最大5日徑流量和年最大洪峰流量作為極端徑流的研究對象。首先采用線性趨勢分析，對4個極端徑流指標進行趨勢分析，然后采用Mann-Kendall 方法對其進行突變檢驗，從而確定極端徑流指標的突變時間點。

采用線性趨勢分析時，研究發現年最大1日徑流量、年最大3日徑流量、年最大5日徑流量和年最大洪峰流量均表現為下降的趨勢，將其趨勢線方程匯總在表1。從結果可知，最大5日徑流的下降率最為明顯，年最大3日次之，年最大洪峰流量相比之下降率最小。但4 個極端徑流指標總體均表現為下降趨勢，較為一致。

表1 各個水文特征值趨勢線方程

采用Mann-Kendall 檢驗方法對澄碧河水庫流域的上述4個極端徑流指標進行趨勢分析和突變檢驗，結果見圖1。從結果來看，年最大1日、年最大3日、年最大5日徑流的趨勢上具有較強的一致性。UF曲線的正負顯示了特征值的變化趨勢，圖1表明，年最大1日、年最大3日、年最大5日徑流在1990年之后，UF一直小于0，說明三者均表現出明顯的下降趨勢，且這一下降趨勢在1998年左右超出了α= 0.05 的顯著水平，表明下降趨勢顯著。相比之下，年最大洪峰流量，雖然在1983年之后呈現下降趨勢（UF＜0），但未超出置信區間，因此其表現為不顯著下降。此外，由圖可知UF和UB的曲線出現交點，且該交點位于置信區間范圍內部。從交點位置來看，年最大1日徑流量、年最大3日徑流量、年最大5日徑流量的交點位置均出現在1991年，而年最大洪峰流量的交點較多，為1973、1976、2002 和2013年。綜合各檢測結果，確定澄碧河流域的氣候突變發生在1991年。

3.2 序列分段與樣本抽取

根據檢測出的突變點，將整個樣本分為兩個時段：天然時段和氣候變化影響時段。研究中天然時段為1963-1991年，氣候變化影響時段為1992-2017年。因此，分別以兩段時間為基礎進行樣本抽選。

樣本抽選是頻率分析的基礎，針對GEV和GPD兩種不同的分布模型，研究采用不同的適應性抽樣方法。對于GEV，結合“最大值”原則，本研究中選取最大1日徑流量、年最大3日徑流量、年最大5日徑流量和年最大洪峰流量為GEV 分布擬合的基礎數據；對于GPD，本研究，以年均超過定量發生次數2 倍為原則，針對突變前和突變后分別以58 和52 為樣本容量進行樣本抽取，針對澄碧河水庫的匯流時間，以3 d 作為抽樣間隔，以保證每個樣本的獨立性。

3.3 分布模型擬合及優選

研究選取了GEV和GPD分布模型進行擬合與對比分析，其突變前、突變后時段擬合的PP 圖和QQ 圖為圖2和圖3，同時將擬合模型的各個參數值和K-S 模型檢驗結果匯總在表2。圖表中am 和pot 分別表示年最大取樣法和超定量取樣法，pre 和pos分別表示氣候突變前、后時段，1 d、3 d、5 d 和peak 分別對應最大1日徑流量（萬m3）、年最大3日徑流量（萬m3）、年最大5日徑流量（萬m3）和年最大洪峰流量（m3/s）。

從GEV分布和GPD分布模型的PP圖和QQ圖可知，兩個分布均能較好地反映經驗觀測值與模擬值之間的關系，從散點圖與45°對角線的關系可以看出，兩者也都緊密地環繞在對角線周圍。通過K-S 檢驗量化的統計結果（表2）發現：①所有情況下的統計量值均小于對應的顯著性標準值。由此表明，GEV 分布和GPD 分布均能通過α=0.05 的顯著性檢驗標準。②除了天然時段的年最大5日徑流量（pre-5d）和氣候變化影響時段的年最大1日徑流量（pos-1d）之外，其他所有情況下的K-S 檢驗結果中GPD 分布對應的統計量均小于（或等于）GEV 分布對應的統計量。因此綜合上述比較分析可知，雖然GEV分布和GPD分布均能很好地反映澄碧河流域的極端徑流序列，但相比之下，GPD分布在擬合澄碧河流域的極端徑流方面更具優勢。

表2 GEV和GPD分布模型參數及K-S檢驗結果

從GPD 分布的參數變化來看，突變后時段相比突變前時段，年最大1日徑流量和年最大洪峰的形狀參數均有所減小，而年最大3日徑流和年最大5日徑流則表現為形狀參數增大。尺度參數方面，除最大1日徑流增大外，其余3 個指標均有所減小。而位置參數方面，只有年最大洪峰流量表現為增加，其余則都是減小的趨勢。由此表明，分布參數的變化在一定程度上受到氣候變化的影響。

3.3 氣候變化對各重現期下極端徑流的影響

通過分布模型的對比分析，確定澄碧河流域的極端徑流頻率分析模型為GPD 模型。因此，以POT 抽取的樣本觀測數據為基礎，采用GPD 分布進行數據擬合，采用極大似然估計進行參數估計得到的線型可以作為重現期計算的依據。根據重現期與累積概率分布函數的關系可以計算天然時段和氣候變化影響時段的各重現期及其對應的極端徑流，進而研究在氣候變化條件下各重現期對應的極端徑流變化情況。

對計算的兩個時段的重現期進行分析，結果見圖4。

從圖4和表3變化可知，突變后時段的各個極端徑流指標在不同重現期水平下，較突變前均有不同程度的減小。其中，年最大洪峰變化最為明顯，10年一遇重現期水平下，其值由550.54 m3/s 減小至448.09 m3/s，減小比例為18.61%；100年一遇重現期水平下，其值由1 141.02 m3/s 減少至602.63 m3/s，減小比例為47.18%；1 000年一遇重現期水平下，其值由2 130.04 m3/s減少至718.52 m3/s，減小比例為66.27%。其次是最大1日徑流量，其減小比例為23.03%～65.51%。相比之下，年最大5日徑流量變化最小，其減小比例為12.95%～21.94%；由此可見，氣候變化在一定程度上確實影響了極端徑流的頻率，使其發生了變化。

表3 澄碧河流域氣候變化前后極端徑流重現期

4 討論與結論

4.1 討 論

受氣候變化的影響，處于巖溶區的澄碧河流域的極端徑流序列平穩性發生了變化，作為極端徑流指標的洪峰和洪量序列的平穩性遭到了破壞，都由原本的“平穩序列”轉變為了“非平穩序列”。綜合分析，序列的突變年發生在20 世紀90年代初。研究表明亞洲區域，尤其是中國范圍內，氣候變化在90年代表現最為劇烈，這在一定程度上佐證了本研究的正確性。此外，據廣西百色市《凌云縣縣志》考證，澄碧河流域自20 世紀70年代之后就沒有再有過大規模的人類活動，因此造成極端徑流頻率變化的原因應該還是來自于氣候變化。氣候變化影響著流域的氣溫、降水、蒸發等水文循環過程，最終表現在極端徑流上［11］。

從結果來看，氣候變化影響最大的是年最大洪峰流量和年最大1日徑流，這可能與澄碧河流域的山區地形和巖溶區地貌有關。澄碧河流域弄塘以上流域多為巖溶區，巖石裂隙、洼地、溶溝、溶槽和地下暗河等地質形態發育而成了交錯縱橫的地下河網體系。加上該地區土壤相對疏松，水流垂向運動頻繁，暴雨轉化為徑流的效率較高，一旦發生極端降水形成大洪水，也能很快被削峰坦化。從流域的匯流面積及匯流時間分析，即便是大洪水，基本也能在3 d 左右完成匯流過程。所以該流域對于短時降水引發的極端徑流，表現是較為明顯的。這與研究得到的年最大洪峰流量和年最大1日徑流變化率較大，而年最大5日徑流量變化率最小的結論相符。通過巖溶區釋水過程及徑流組分的相關研究［12，13］可知，由于巖溶洼地、溶溝、溶槽和地下暗河等地質構造的存在，使得巖溶區流域釋水過程要明顯快于非巖溶區流域，其衰減系數也相對較大。此外，巖溶區徑流組分隨洪峰的變化也較非巖溶區劇烈，這也是在氣候變化影響下，極端徑流頻率變化明顯的重要原因。

4.2 結 論

本文以澄碧河流域為實例，研究了氣候變化對巖溶區流域極端徑流頻率分析的影響，發現如下結論。

（1）澄碧河流域的氣候突變時間發生在1991年，可將其極端徑流序列分成突變前（天然時段）和突變后（氣候變化影響時段）兩個部分進行研究。

（2）GEV 分布和GPD 分布均能較好地擬合澄碧河流域的4個極端徑流序列，但相比之下，GPD分布表現更為突出。

（3）氣候變化導致該流域4 個極端徑流指標在不同重現期時都出現了不同程度的減小。其中，年最大洪峰流量減小比例最大，為18.61%～66.27%，年最大1日徑流量次之，年最大5日徑流量減小比例最小。

（4）巖溶區流域獨特的地質構造可能加劇了氣候變化對極端徑流頻率的影響。 □