劉洋, 陳菡, 譚學志*

(1. 中山大學 土木工程學院水資源與環(huán)境研究中心，廣州510275; 2. 湖北省水利水電規(guī)劃勘測設計院, 武漢 430064)

0 引言

流域河川徑流的歸因識別是制定氣候變化應對策略和水資源可持續(xù)利用的基礎[1]。在全球氣候變化與人類活動對流域徑流的共同作用下，全球諸多流域發(fā)生了明顯的徑流變化[2]。越來越多的學者開始著眼于氣候變化與人類活動對流域徑流變化影響及歸因的研究。其中對徑流變化的歸因方法有許多。部分研究采用基于Budyko理論的流域水熱耦合平衡方程對流域徑流變化進行歸因的方法，并取得了較好的研究結果[3-8]。該方法與通常采用的分布式流域水文模型相比較，不僅方法簡單，且避免了模型參數不確定性對歸因結果的影響[8]。

韓江流域河流開發(fā)利用較早，流域內修建有多處水利工程，取水引水等人類活動較為頻繁。流域內人口稠密，經濟發(fā)達，經濟社會的日益發(fā)展對區(qū)域水資源的需求日益增加。分析韓江流域徑流變化及對其成因的量化對該流域水資源合理開發(fā)與可持續(xù)利用等方面有著重要意義。在以往的研究中，已有學者對韓江流域及其子流域的河川徑流變化進行了分析與研究[9-11]，但關于徑流變化歸因的量化分析研究相對較少，也缺乏對各影響因子貢獻的量化分析。本研究從氣候變化和人類活動兩方面，基于Budyko理論的水熱耦合平衡方程與彈性系數理論，對韓江流域及其子流域即梅江流域、汀江流域的徑流變化進行敏感性分析，并對各因子對徑流變化造成的貢獻進行定量歸因分析，旨在為韓江流域在河流開發(fā)策略、水災害防控、水資源合理規(guī)劃與可持續(xù)利用等方面提供科學依據。

1 研究區(qū)域與數據來源

1.1 研究區(qū)域概況

韓江流域位于粵東和閩西南地區(qū),地理位置處于東經115°13′～117°09′，北緯 23°17′～ 26°05′。韓江流域水系主要由流域上游的梅江子流域與汀江子流域組成，梅江與汀江于大埔三河壩匯合后，形成韓江干流，而后向南延伸，最終于汕頭匯入南海。其中，韓江干流全長470 km，總流域面積30 112 km2;梅江子流域面積為13 929 km2;汀江子流域面積為11 802 km2。韓江流域屬亞熱帶季風氣候，多年平均氣溫介于20～21.5 ℃之間，多年平均降雨量約為1 620 mm,氣候溫和，雨量充沛，流域內大部分地區(qū)植被覆蓋率高，生態(tài)資源豐富[12]。

1.2 數據來源

表 1 韓江流域水文站Table 1 Hydrological stations in Han River Basin

圖 1 韓江流域水文站點與氣象站點分布 Figure 1 Locations of hydro-climatic stations in Han River Basin

采用韓江流域潮安、溪口、橫山3個代表性水文站點的1960—2000年共41年逐日平均流量數據資料。其中，潮安水文站為韓江主要控制站，橫山水文站為梅江主要控制站，溪口水文站為汀江主要控制站，站點詳細信息見表1。此外，本研究采用韓江流域及其周邊區(qū)域的7個氣象站點1960—2000年共41年降水及氣象數據(圖1)。其中，用于計算日潛在蒸散發(fā)的氣象指標包括逐日氣壓、氣溫、風速、相對濕度和日照等數據。采用的氣象站數據均來自中國氣象科學數據共享服務網提供的《中國地面氣候資料日值數據集》。

2 研究方法

2.1 徑流序列變化突變分析

采用非參數Mann-Kendall突變點檢驗法與累積距平法來檢測年徑流系列的突變點，以尋找徑流序列變化趨勢發(fā)生轉變的年份。并將兩種方法得到的突變點進行對比，互相驗證，結合實際情況分析，以確定徑流序列變化趨勢發(fā)生突變的具體年份。其中，(1)非參數Mann-Kendall突變點檢驗法是一種廣泛應用于水文、氣象等領域的非參數檢驗法[13]。(2)累積距平法可以通過觀察數據序列差積曲線的變化來判斷數據點的離散程度和變化趨勢，并找到數據序列的突變位置[13]。在確定突變點后，將徑流系列劃分為突變點前與突變點后2個子序列，并對突變前后子系列的各水文氣象要素進行分析，進而對造成徑流突變的各影響因素進行歸因分析。

2.2 徑流變化敏感性分析與歸因識別

(a) 多年尺度下的流域水量平衡方程

流域水量平衡方程為：

P=E+R+?S

(1)

式(1)中：R為徑流深，P為降水量，E為實際蒸散發(fā)量，ΔS為水量變化量。在較長時間尺度的計算中，ΔS可忽略不計，則式(1)簡化為下式：

R=P-E

(2)

(b) Choudhury-Yang公式

Budyko假設[14]認為在一定氣候和植被條件下，流域水分和能量之間存在著水熱耦合平衡關系，流域多年平均降水量、潛在蒸散發(fā)量和實際蒸散發(fā)量的關系可以用經驗曲線來描述[15]。在此基礎上，諸多學者通過數學分析得到了一些衍生公式。其中，Choudhury[16]與Yang[17]等基于經驗法或解析法等方法分別得到了形式一致的Choudhury-Yang公式，本研究采用該公式描述流域水量與能量的平衡關系，具體形式如下：

(3)

式(3)中：ET0—長期平均的年潛在蒸散發(fā)量，n—反映流域下墊面特征的參數。ET0采用世界糧農組織提出的Penman-Monteith公式[18]計算：

(4)

式(4)中：ET0—可能蒸散量(mm/d)，Δ—飽和水氣壓曲線斜率(kPa/℃)，Rn—地表凈輻射(MJ/(m×d))，G—土壤熱通量(MJ/(m2×d))，γ—干濕表常數(kPa/℃)，Tmean—日平均溫度(℃)，u2—2 m高處風速(m·s-1)，es—飽和水氣壓(kPa)，ea—實際水氣壓(kPa)。

(c) 徑流影響因素的彈性系數

(5)

(6)

(7)

在敏感性分析中，為判斷各影響因子的彈性系數年際間的變化趨勢，本研究采用M-K趨勢檢驗方法來對其趨勢進行顯著性檢驗[8]，顯著性水平設為α=0.05。

(d) 徑流影響因素對徑流變化的貢獻率

由上述突變點檢驗方法將徑流序列分為突變點前、突變點后兩個子序列，記突變點前子序列為R1，突變點后子序列為R2，則△R=R2-R1，降水、潛在蒸散發(fā)與下墊面特征指數的突變前后變化同理。將徑流的變化量分解為降水變化引起的徑流變化量△RCP，潛在蒸散發(fā)變化引起的徑流變化量△RCET0，下墊面特征指數變化引起的徑流變化量△RCn，則

(8)

(9)

(10)

將各影響因子造成的徑流變化量相加為△RC=△RCP+△RCET+△RCn，則流域氣象要素與下墊面變化對徑流變化的貢獻率為：

(11)

(12)

(13)

式(11～13)中:CP為降水變化對流域徑流變化的貢獻率，CET0為潛在蒸散發(fā)變化對流域徑流變化的貢獻率，Cn為下墊面特征指數變化對流域徑流變化的貢獻率。其中，CP與CET0之和可看作氣候變化對流域徑流變化的貢獻率，Cn可看作人類活動對流域徑流變化的貢獻率。

3 結果與分析

3.1 徑流變化趨勢及突變點檢驗

采用非參數Mann-Kendall突變點檢驗法對潮安、橫山、溪口的流域徑流過程發(fā)生突變的年份進行識別，得到突變檢驗結果如圖2所示，并以此作為依據分析韓江流域徑流突變點發(fā)生的年份。綜合考慮3個水文站的累積距平趨勢，其中3個水文站的累積距平曲線均于1960—1982年呈快速下降趨勢，于1982—2000年呈平緩上升趨勢。韓江流域3個水文站點均于1961與1973年出現(xiàn)劇烈變化，尤其以溪口站于1973年的起伏最為明顯。此外，溪口站還于1990年出現(xiàn)劇烈變化。由累積距平法得到的結果分析，潮安站和橫山站徑流發(fā)生突變的年份為1972、1982年前后，溪口站徑流發(fā)生突變的年份為1972、1982、1991年前后。由M-K檢驗法得到的結果分析可知，潮安站徑流發(fā)生突變的年份為1982、1984、1988年前后，橫山站徑流發(fā)生突變的年份為1975、1979年前后，溪口站徑流發(fā)生突變的年份為1989、1991年前后。根據3個站點的累積距平曲線與M-K突變檢驗法結果，以韓江干流的突變時間作為全流域徑流發(fā)生突變的主要依據，綜合實際情況分析，于1960—2000年時間序列中，取1982年為韓江全流域徑流發(fā)生突變的時間點。

為驗證結論的可靠性，與已有的研究進行對比?？娺B華[9]等以1951—2010年潮安水文站徑流時間序列進行分析，發(fā)現(xiàn)其突變點為1980年，與本研究選取1982年為潮安站發(fā)生突變年份的結論較為接近，考慮到序列時間長度選取的不同，故結論略有差異，認為可以接受。有關臨近流域的突變年份，賴天锃[5]等對珠江流域下東江流域1959—2000年多個水文站徑流量進行突變分析，發(fā)現(xiàn)其突變時間點為1972年。同時，根據李艷[20]等得出的結論：珠江流域下的北江水系在1980年后的多年平均徑流量比1980年以前更大。以上各項研究雖然研究區(qū)域與選取的時間序列不盡相同，用于判斷突變點發(fā)生的方法也不盡相同。但綜合來看，珠江流域下各大水系基本于在二十世紀七八十年代年前后發(fā)生突變，且突變后流量大致呈上升趨勢，故認為本研究取1982年為韓江流域徑流突變年份大致合理。

圖 2 潮安、橫山、溪口水文站年平均徑流深累積距平法與M-K法檢驗結果Figure 2 Cumulative departure curves and M-K test curves of annual runoff depth at three stations

3.2 徑流變化影響因素敏感性分析

由公式(5～7)計算韓江流域、梅江流域和汀江流域的徑流對降水、潛在蒸散發(fā)、下墊面特征指數的彈性系數(1960—2000年)，其結果如圖3所示。各影響因子的彈性系數年際波動較大，曲線起伏較為明顯?？傮w而言，徑流變化與降水正相關，與潛在蒸散發(fā)和下墊面特征指數均呈負相關。降水彈性系數介于1.23～2.22之間，潛在蒸散發(fā)彈性系數介于-1.22～-0.23之間，下墊面特征指數彈性系數介于-1.07～-0.18之間。各彈性系數取值表明，降水每增加1%，將導致徑流量增加1.23%～2.22%，潛在蒸散發(fā)每增加1%，將導致徑流量減少0.23%～1.22%，下墊面特征指數每增加1%，將導致徑流量減少0.18%～1.07%。

表 2 徑流對各影響因子彈性系數的M-K趨勢檢驗結果Table 2 The M-K trend test results of elastic coefficient of runoff on various impact factors

表 3 徑流對各影響因子的多年彈性系數Table 3 Multi-year elastic coefficient of runoff on various impact factors

圖 3 1960—2000年潮安、橫山和溪口水文站實測徑流對(a)降水、(b)潛在蒸散發(fā)和(c)下墊面特征指數的彈性系數Figure 3 Elastic coefficient of runoff depth on precipitation (a), potential evapotranspiration (b) and underlying surface index (c) at three stations during 1960—2000

表 4 韓江流域徑流變化歸因分析Table 4 Attribution analysis of runoff changes in Han River Basin

對比3個站點3種彈性系數可得，梅江子流域的絕對值普遍較高，汀江子流域的絕對值普遍較低，韓江流域則普遍位于二者之間。即梅江流域對各影響因子的敏感性最高，汀江流域對各影響因子的敏感性最低，韓江流域總體的敏感性介于梅江子流域與汀江子流域之間。由表2的M-K趨勢檢驗結果(顯著性水平α=0.05)可知，降水彈性系數呈不顯著下降趨勢，潛在蒸散發(fā)與下墊面特征指數的彈性系數呈不顯著上升趨勢，表明在韓江及其子流域內，降水對徑流變化的正相關性有略微減弱的趨勢，潛在蒸散發(fā)與下墊面特征指數與徑流變化的負相關性有略微減弱的趨勢，但3種彈性系數的年際波動均較大，總體變化趨勢不明顯。

由3個流域徑流與氣候數據的多年平均值，通過公式(5～7)計算得到3個研究區(qū)域多年尺度的彈性系數，如表3所示。其結果與年際變化彈性系數的結果分析一致，梅江子流域對3個影響因素的敏感性最高，汀江子流域的敏感性最低，韓江流域處于二者之間。且韓江流域徑流變化均與降水呈正相關，與潛在蒸散發(fā)與下墊面特征指數均為負相關。

3.3 徑流變化歸因分析

由公式(11～13)計算出韓江流域、梅江子流域、汀江子流域于1982年前后氣候變化與人類活動對徑流變化的貢獻率(表4)。各影響因子造成的徑流變化量之和(△RC=△RCP+△RCET0+△RCn)與實測徑流變化量△R差異很小(0.13～0.4 mm)，因此說明本研究采用的歸因方法是有效可行的。1960—2000年，氣候變化是引起韓江各流域徑流變化的主要因素，其中，又以降水變化影響為主，降水變化對韓江流域、梅江子流域和汀江子流域的徑流變化貢獻率分別為58.24%、52.72%和51.70%，潛在蒸散發(fā)變化對韓江流域、梅江子流域和汀江子流域的徑流變化貢獻率分別為10.80%、5.96%和25.75%，其中汀江子流域的潛在蒸散發(fā)變化的貢獻較高，其原因是汀江子流域的潛在蒸散發(fā)變化較大，達-53.59 mm(表5)。下墊面特征指數變化也對韓江各流域徑流變化起到了很大的影響，其中對梅江子流域的影響最大，高達41.32%，對汀江子流域的影響最小，為22.55%，對韓江總流域的徑流變化影響為30.96%，其中梅江子流域的下墊面特征指數變化的影響較高，其原因是梅江子流域的下墊面特征指數變化較大，達-0.15(表5)。由此可見，在韓江流域，氣候變化是徑流變化的主要因素，約占60%～70%，其中大部分貢獻由降水變化造成，而人類活動對徑流變化的影響也同樣重要，約占30%～40%。

表 5 韓江流域各影響要素1982年前后變化Table 5 Changes in the impact factors of the runoff changes in Han River Basin before and after 1982

4 總結與討論

4.1 氣候變化對韓江流域徑流變化的影響

韓江流域及其子流域梅江、汀江流域對氣候變化均呈現(xiàn)相似響應。降水變化是氣候變化對流域徑流變化產生影響的關鍵因素，其變化往往與流域徑流量的變化有直接關系。在本研究中，根據敏感性分析結果，降水是韓江流域徑流變化的主要影響因子，對徑流的貢獻呈正相關，且在3個影響因子中其影響程度最高。在本研究中，潛在蒸散發(fā)是氣候因素中的次要因素，對徑流變化的影響相對較小,且與徑流呈負相關關系，這與目前許多研究的觀點一致[1,4,6-8]，即潛在蒸散發(fā)對流域徑流變化的彈性系數為負。由表5可知，于1982年前后，韓江流域及其子流域降水量均上漲，潛在蒸散發(fā)量則在一定程度上降低，故二者變化皆導致了韓江流域徑流量的增加。

4.2 人類活動對韓江流域徑流變化的影響

依據表5中n值在兩個時間段的變化，韓江流域的n值于1982年前后有所減小，由于敏感性分析的結果，n值與徑流變化負相關，因此n值的減小也導致了徑流量的增加。根據現(xiàn)有研究，通常認為流域下墊面條件中的地形、地質和土壤條件相對穩(wěn)定，因此下墊面參數n的變化主要由土地利用變化/植被覆蓋變化導致。一般認為n值的增加是由流域內植被覆蓋的提高導致的[21]。在本研究中，韓江流域及其子流域的下墊面特征指數均呈現(xiàn)5%～10%左右的下降，可以認為這主要是流域內土地利用變化所致。根據黃苑君的研究[22]，1986—2006年韓江流域未利用地面積減少62 076.0 hm2，建設用地增加48 559.3 hm2，建設用地增加導致地表硬化，流域下墊面不透水面積所占比例增加，致使流域內徑流增加。

5 結論

本研究采用累積距平法與M-K突變點檢驗法分析了韓江流域及其子流域梅江與汀江流域在1960—2000年徑流量變化，根據基于Budyko理論的水熱耦合平衡方程，結合彈性系數法，對韓江流域干支流對各影響因子的敏感性及其年際變化進行了分析，并從氣候變化與人類活動兩方面對韓江流域徑流變化產生的影響，進行了定量歸因分析。主要結論如下：

(1)結合累積距平法與M-K突變點檢驗法的結果分析，韓江流域干支流于1982年前后發(fā)生突變，1982年后韓江流域3個水文站點的多年平均徑流量均顯著增加。

(2)基于Budyko理論的水熱耦合平衡方程，通過彈性系數法對韓江流域干支流進行敏感性分析。在彈性系數的年際變化中，各影響因子的彈性系數年際波動較大，降水彈性系數呈不顯著下降趨勢，潛在蒸散發(fā)與下墊面特征指數的彈性系數呈不顯著上升趨勢，但總體變化趨勢不明顯。總體而言，徑流變化與降水呈正相關，與潛在蒸散發(fā)和下墊面特征指數均呈負相關。

(3)對韓江流域徑流變化進行定量歸因分析，結果表明：氣候變化對韓江流域徑流變化的貢獻率為59%～77%，其中氣候變化的貢獻主要體現(xiàn)在降水的變化，對徑流變化的貢獻率約為50%左右。人類活動對韓江流域徑流變化的貢獻率為22.6%～41.3%，其中梅江子流域的人類活動貢獻率較高，汀江子流域的人類活動貢獻率相對較低。