極端氣象氣候條件下土地利用變化的水文水資源效應研究
——以東江流域為例

秦美茵，董春雨，孫秀峰

(1.中山大學 土木工程學院，廣東 珠海 519082；2. 廣東省水利水電科學研究院，廣州 510635)

1 概述

土地利用變化作為影響流域水文循環(huán)的主要驅動因素之一，能夠顯著影響蒸散發(fā)、降雨截留、地表徑流及土壤下滲等水文過程，進而改變河道的徑流量，并引發(fā)洪旱災害、水資源短缺等問題[1]。目前國內外已有較多學者對土地利用變化對徑流的影響展開研究。Matens等[2]采用試驗流域法對多個土地利用條件不同流域的年徑流量進行分析，研究表明園地和裸地年徑流系數最大，而自然植被的徑流系數最小。歐陽欽從等[1]對晉江流域土地利用變化對水文過程影響進行分析，認為降雨量越小，蒸散發(fā)和地表徑流的效應越顯著，而土壤剖面下滲和地下徑流效應則相反。呂樂婷等[3]對土地利用變化對東江流域產流過程影響進行定量評估，認為地表徑流、側向流對土地利用變化的響應在雨季更為敏感，實際蒸發(fā)、潛在蒸發(fā)及徑流總量對土地利用變化的響應在旱季更為敏感。目前相關研究大多集中于研究土地利用變化對年徑流、月徑流及產流過程的響應分析，對極端氣象氣候條件下土地利用變化的水文效應研究相對較少。

近百年來，全球氣候持續(xù)變暖，地表氣溫升高，在一定程度上加劇了蒸散發(fā)等水循環(huán)過程，導致暴雨和干旱等災害事件發(fā)生的頻率和強度增加[4-5]。極端氣象氣候條件下，土地利用對水文過程的影響可能會有放大效應，但目前學界仍沒有明確結論，亟需深入研究[6-8]。因此，研究極端氣象氣候條件下土地利用變化的水文效應對于地方防災減災、生態(tài)保護和國土空間規(guī)劃等具有重要意義。

東江流域位于經濟發(fā)達的東部沿海地區(qū)，1990年以來，其土地利用結構變化主要表現為城市化過程顯著，非城鎮(zhèn)用地結構變化明顯，以及林地和灌草地減少的總體特征；隨著粵港澳大灣區(qū)的建設開發(fā)，未來該流域內土地開發(fā)利用活動將更加頻繁[9-10]。近半個多世紀以來，珠江流域極端高溫及短時間極端降水事件持續(xù)增多，極端高溫事件的增多，降水的集中以及干期的變長使珠江流域面臨著高溫干旱和暴雨洪水的威脅[11]。目前關于東江流域土地利用變化對水文過程的影響分析主要基于實際發(fā)生的土地利用變化，少有針對未來不同土地利用情景對水文過程和水資源的影響進行相關探討。因此，本文以東江流域為案例研究區(qū)，基于SWAT水文模型，模擬1990年、2018年以及另外六種城鎮(zhèn)用地及林地面積遞增的土地利用情景下，不同程度強降水和干旱條件下的流域水文過程和水資源量變化規(guī)律，探討極端氣象氣候條件下土地利用變化的水文效應和水資源效應，旨在為東江流域的土地利用規(guī)劃、水資源配置和防洪減災提供科學依據。

2 研究區(qū)概況與數據來源

2.1 研究區(qū)概況

東江流域位于珠江三角洲的東北端，地理坐標為113°52′～115°52′E，22°38′～25°14′N。流域多年平均氣溫為20℃～22℃，多年平均降雨量為1 500～2 400 mm，降雨時空分布不均，西南多、東北少；降水年際變化大，年內分配也極不均勻，主要集中在4—9月[12]。本研究選取東江干流下游主要控制站博羅站以上區(qū)域作為研究區(qū)(見圖1)，流域控制面積為25 325 km2，占東江流域總面積的71.7%。土壤類型主要有紅壤、赤紅壤和水稻土；土地利用類型以林地為主，其次為耕地、草地、建筑用地等。

圖1 研究區(qū)概況示意

2.2 數據來源

本研究采用的數據包括研究區(qū)數字高程模型(Digital Elevation Model，DEM)、氣象、土壤、土地利用及實測徑流量數據。DEM分辨率為30 m，來源于2020年2月NASA發(fā)布的NASADEM數據。氣象數據來自國家青藏高原科學數據中心孟現勇建立的CMADS V1.1數據集[13]，可直接用于驅動SWAT模型，數據水平空間分辨率為 0.1°，時間分辨率為3 h，數據覆蓋時段為1979—2018年，數據要素包含流域內日均降雨量、日最高和最低氣溫、風速、相對濕度和太陽輻射值等。土壤數據采用聯(lián)合國農糧組織和維也納國際應用系統(tǒng)研究所于2009年發(fā)布的世界土壤數據庫(Harmonized World Soil Database，HWSD)，分辨率為1 km。土地利用數據來源于中國科學院資源環(huán)境科學數據中心，分辨率為30 m。

3 研究方法

3.1 SWAT模型

SWAT(Soil and Water Assessment Tool)水文模型是由美國農業(yè)部農業(yè)研究中心于1994年開發(fā)的，近年來在水文水資源領域已得到廣泛應用[14]。其模擬的流域水文過程分為產流、坡面匯流、河道匯流過程[15]，產流和坡面匯流過程描述的是子流域內水流、泥沙、營養(yǎng)物等匯入主河道的過程，河道匯流過程描述的是流域河網中水流、泥沙等的輸移過程[16]。SWAT模型根據地形及河網將流域劃分為多個子流域，再將子流域中土地利用和土壤類型的唯一組合定義為水文響應單元(HRU)，HRU是SWAT水文模擬的基本單位，1個子流域可劃分為1個或多個HRU[16]。SWAT模型中的水量平衡方程如下：

(1)

式中:

SWt——土壤最終含水量；

SW0——第i天土壤初始含水量；

t——時間；

Rday——第i天的降水量；

Qsurf——第i天地表徑流；

Ea——第i天的蒸散發(fā)量；

Wseep——第i天離開土壤剖面底部的滲透水流和旁通水流水量；

Qgw——第i天回歸流的水量。

SWAT模型可在控制其它輸入數據不變的前提下輸入不同類型的土地利用數據，以獲得不同土地利用情景下不同時間尺度的徑流量以及藍綠水的分量(子流域產水量、深層含水層補給量和蒸散發(fā)、土壤含水量)。

3.2 藍綠水定義及計算

近年來，藍綠水概念已被廣泛應用于水資源管理與評估。藍水是指地表徑流和地下徑流，即可直接被人類利用的水資源，包含深層含水層補給、湖泊、濕地儲水量；綠水由綠水流和綠水儲組成，綠水流指實際蒸散量，綠水儲指儲存在土壤中的水分[17]。根據藍水、綠色的定義以及SWAT模型的輸出結果可知，藍水可用SWAT模型輸出的子流域產水量(WYLD)以及深層含水層補給量(DA_RCHG)之和表示；綠水流為實際蒸散發(fā)量(ET)；綠水儲為土壤含水量(SW)，綠水即為實際蒸散發(fā)量(ET)與土壤含水量(SW)之和[18-19]。

3.3 帕爾默干旱指數(PDSI)

帕爾默指數(Palmer Drought Severity Index, PDSI)是表征一段時間內某地區(qū)實際水分供應持續(xù)少于當地氣候適宜水分供應的水分虧缺情況[20]，其數值落于區(qū)間(-∞，-4]時為極度干旱，落于區(qū)間(-4，-3]時為嚴重干旱，落于區(qū)間(-3，-2]時為中度干旱，落于區(qū)間(-2，-1]時為輕度干旱。本文應用遙感云平臺Google Earth Engine提取1979—2018年間東江流域PDSI數據集，統(tǒng)計不同干旱等級情景下的藍綠水變化情況，以研究東江流域在不同土地利用變化情景下藍綠水資源對不同干旱程度的響應特征。

3.4 土地利用情景構建

根據SWAT模型的分類標準，將東江流域土地利用重分類為林地(FRST)、耕地(AGRL)、草地(PAST)、城鎮(zhèn)用地(URHD)、村莊(URLD)、工業(yè)用地(UIDU)、裸地(BALD)和水域(WATR)8類。為研究徑流量對不同土地利用情景的響應、不同土地利用變化情景下徑流量對不同降雨等級的響應，以及不同土地利用變化情景下藍綠水對不同干旱等級的響應，本文將對2個實際土地利用情景及6個假設土地利用情景下的水文過程和水資源量進行對比分析。各土地利用情景分類見表1。

表1 研究區(qū)土地利用情景分類

3.5 SWAT模型構建及評價

本文首先利用SWAT模型進行子流域劃分，通過DEM設置、河網定義、修改閾值、流域總出口選擇以及子流域參數計算等過程，完成子流域劃分，共劃分子流域(subbasin)101個，水文響應單元(HRU)1 943個。隨后依次輸入處理后的土壤數據、土地利用數據(模型參數率定使用1990年土地利用數據)，設置流域坡度分級對水文響應單元進行定義；輸入處理后的CMADS氣象數據集，數據包含日均降雨量、日最高/最低氣溫、風速、相對濕度和太陽輻射值；定義流域內三大水庫(新豐江、白盆珠、楓樹壩水庫)出口，輸入水庫相關參數[21]。徑流模擬選用SCS曲線法，潛在蒸散發(fā)模擬選用Penman-Monteith方法，河道驗算選用Muskingum公式法。配置好SWAT模型后，月徑流模型分別選取博羅站1979—1998年、1999—2007年實測徑流數據作為率定期和驗證期，其中模型預熱期2 a；日徑流模型分別選取博羅站1989—2003年、2004—2007年實測徑流數據作為率定期和驗證期，其中模型預熱期2 a。

利用SWAT-CUP中SUFI-2算法對SWAT模型輸出結果進行參數敏感性分析、率定及驗證。選用確定性系數(R2)和Nash-Suttclife系數(Ens)，評估模型模擬值和實測值的擬合效果[22]。經過率定，流域出口博羅站月徑流(見圖2)率定期和驗證期確定性系數R2和Nash-Suttclife系數Ens分別為0.87、0.87和0.89、0.88，P-factor和R-factor分別為0.64和0.68；日徑流(見圖3)率定期和驗證期確定性系數R2和Nash-Suttclife系數Ens分別為0.75、0.74和0.78、0.78，P-factor和R-factor分別為0.85和0.84；表明SWAT模型模擬東江流域水文過程適用性良好，精度較高，可以開展不同土地利用情景的水文水資源效應分析。

圖2 博羅站月徑流率定期(1981—1998年)、驗證期(1999—2007年)模擬結果示意

圖3 博羅站日徑流率定期(1989—2003年)、驗證期(2004—2007年)模擬結果示意

4 結果分析與討論

4.1 土地利用變化分析

由表2、表3可知，林地為東江流域最主要的用地類型(見圖4)，面積占比分別為75.4%(1990年)、74.5%(2018年)，1990—2018年間，林地面積共減少231.6 km2，變化率為-1.2%，減少的林地主要轉變?yōu)楦?、草地、工業(yè)用地及城鎮(zhèn)用地。城鎮(zhèn)用地為東江流域增長量最大的用地類型，1990—2018年，流域內城鎮(zhèn)用地面積共增長321.3 km2，增長率為272.7%，增加的城鎮(zhèn)用地主要源于耕地、林地及草地。

表2 東江流域土地利用結構及變化趨勢 km2，%

表3 1990—2018年東江流域土地利用轉移矩陣 km2，%

1990年

4.2 徑流量對不同土地利用情景的響應分析

基于表1中所示8種土地利用情景，本文利用SWAT水文模型模擬了東江流域下游博羅站逐月徑流量。以1990年土地利用情景下的徑流量為基準，本文計算了其他7種土地利用變化情景引起的徑流量變化。由圖5顯示結果可知，流域內城鎮(zhèn)用地及林地面積的增加均對流域出口博羅站月均徑流量產生不同程度的影響，表現為擴大的城鎮(zhèn)用地面積加劇了徑流量年內分配不均。具體而言，在豐水季(4—9月)，城鎮(zhèn)用地面積的增加導致徑流量增加，且徑流量變化量隨城鎮(zhèn)用地面積增加而增加；在枯水季(10月至次年3月)，城鎮(zhèn)用地面積的增加導致徑流量減少；林地面積的增加在全年各季節(jié)均導致流域出口博羅站月均徑流量減少。

圖5 月均徑流量變化量示意

以上土地利用變化引起的水文響應，原因可能在于：流域城鎮(zhèn)用地面積增加導致不透水面增加，減小了流域總體下滲能力和滯水能力，在豐水季，使得流域徑流量增加；在枯水季，徑流主要由基流組成，而基流一般來源于地下水或其他延遲部分的徑流[23-24]，枯水季徑流量減少，可能是由于不透明水面的增加導致土壤含水量減少，枯水季水量補給減少，從而使得流域徑流量減少。流域林地的增加，一方面導致陸面蒸散發(fā)增加，另一方面伴隨著冠層截留、枯落物層和土壤層持水量增加，從而降低土壤水分和地表徑流量，導致豐水季徑流量明顯減少；而枯水季因降水量較少且植被蒸散發(fā)較弱，林地面積變化引起的徑流量變化相對不明顯。

4.3 土地利用變化情景下徑流量對強降雨的響應分析

根據流域降雨特性及子流域降雨計算結果，本文對流域強降雨事件按24 h降雨量為100～150 mm，150～200 mm和大于200 mm進行劃分，以研究城鎮(zhèn)用地及林地變化情景下徑流量對不同強度降雨事件的響應規(guī)律。據統(tǒng)計，1981—2018年，24 h降雨量在100～150 mm、150～200 mm之間和大于200 mm的強降雨事件分別發(fā)生了127次、26次和5次。

由圖6～圖8可知，發(fā)生強降雨時，林地面積的增加，使得徑流量減少，變化量隨降雨強度增加而增加，表明林地面積的增加對發(fā)生強降雨時的地表徑流的減弱效應顯著，可減緩和削弱洪峰。由圖9～圖11可知，城鎮(zhèn)用地面積的增加導致徑流量產生明顯變化，變化率最大達38%，變化量隨降雨強度及城鎮(zhèn)用地面積增加而增加，表明城鎮(zhèn)用地面積的增加對發(fā)生強降雨時的地表徑流的增強效應顯著，能加速和加劇流域內洪峰的形成，擴大洪澇災害風險。

圖6 遭遇200 mm/d以上強降雨時林地面積擴張引起的徑流量變化示意

圖7 遭遇150～200 mm/d強降雨時林地面積擴張引起的徑流量變化示意

圖8 遭遇100～150 mm/d強降雨時林地面積擴張引起的徑流量變化示意

圖9 遭遇200 mm/d以上強降雨時城鎮(zhèn)用地面積擴張引起的徑流量變化示意

圖10 遭遇150～200 mm/d強降雨時城鎮(zhèn)用地面積擴張引起的徑流量變化示意

圖11 遭遇100～150 mm/d強降雨時城鎮(zhèn)用地面積擴張引起的徑流量變化示意

4.4 土地利用變化情景下藍綠水對干旱的響應分析

為探討不同土地利用情景在發(fā)生干旱時的水資源效應變化規(guī)律，本文進一步分析了表1所示各土地利用情景下，東江流域藍綠水資源量在不同強度干旱條件下的變化情況。據統(tǒng)計，1981—2018年，東江流域分別發(fā)生重度干旱(-3≥PDSI指數>-4)事件3次、中度干旱(-2≥PDSI指數>-3)事件4次、輕度干旱(-1≥PDSI指數>-2)事件3次。

由圖12可知，流域林地面積增加，使得藍水減少、綠水增加，且變化量隨林地面積的增加而增加。在林地面積擴張1.1倍和1.13倍時，綠水儲對林地擴張的響應不敏感，綠水的變化以綠水流的變化為主。林地擴張引發(fā)的綠水流的增加在中度干旱時取得最大值，表明林地擴張的蒸散發(fā)增加效應在中度干旱時達到峰值。發(fā)生重度干旱時林地面積擴張引起的綠水變化最少，可能由于部分森林植被因強烈干旱脅迫而生理受損，使得蒸散發(fā)增加量明顯減小。以上結果與前人研究發(fā)現相一致，即一定程度的干旱會導致植被蒸散發(fā)增強，尤其干旱期間伴隨異常高溫時，植被需要更多水分進行降溫；而當干旱繼續(xù)加劇時，植被葉片氣孔關閉，光合作用降低，甚至因干旱脅迫死亡，蒸散發(fā)下降[25-26]。此外，圖12結果顯示，林地面積擴張導致藍水的減少量也在中度干旱時取得最大值，這與蒸散發(fā)在中度干旱時達到峰值的效應相符，根據水量平衡原理，蒸散發(fā)增大，徑流量減小。

圖12 不同干旱等級情景下林地面積擴張引起的藍綠水變化示意

由圖13可知，流域城鎮(zhèn)用地面積增加，使得藍水增加、綠水減少，且變化量隨城鎮(zhèn)用地面積的增加而增加。在城鎮(zhèn)用地向外擴張2.5倍和5倍時，綠水流對城鎮(zhèn)用地擴張的響應相對不敏感，綠水的變化以綠水儲的變化為主，綠水儲的減少在中度干旱時取得最大值，這可能是由于中度干旱時土壤含水量已接近其下限值，干旱等級繼續(xù)增強后，土壤水分變化量開始減小。藍水變化量在豐水季隨干旱等級增強而減小，反映豐水季的降水隨干旱增強而減少；藍水變化量在枯水季為中度干旱時最大，這可能是由于枯水季藍水有較大比例為地下徑流，而地下徑流又與土壤含水量相關，中度干旱時土壤含水量已接近下限值，干旱等級增強后，土壤水下降速度變慢，變化量開始減小。

圖13 不同干旱等級情景下城鎮(zhèn)用地面積擴張引起的藍綠水變化示意

5 結語

以上結果分析顯示，本文建立的基于SWAT的東江流域水文模型能夠較好地模擬該流域的日徑流、月徑流過程，適用于評估東江流域土地利用變化的水文水資源效應。通過以上分析，本文主要得出以下幾點結論：

1) 通過分析不同土地利用情景下東江流域月均徑流的響應特征，表明城鎮(zhèn)用地面積的增加會加劇徑流年內分配不均，在豐水季，城鎮(zhèn)用地面積的增加導致徑流量增加，且徑流量變化量隨城鎮(zhèn)用地面積增加而增加；在枯水季，城鎮(zhèn)用地面積的增加導致徑流量減少；林地面積的增加在全年均導致流域出口月均徑流量減少。

2) 對不同土地利用變化情景下日徑流對不同強度降雨的響應分析結果表明，城鎮(zhèn)用地及林地面積的增加對發(fā)生強降雨時的地表徑流的增強、減弱效應顯著，且變化效應隨土地利用結構變化強度增加而增加，因此，城鎮(zhèn)用地的持續(xù)擴張及植樹造林導致的土地利用結構變化都會對流域洪水過程產生顯著影響。

3) 對土地利用變化情景下藍綠水資源對不同干旱等級的響應分析結果表明，城鎮(zhèn)用地及林地面積的增加對發(fā)生干旱時藍綠水資源量總體影響較小，綠水儲對林地擴張的響應不敏感，林地擴張引起的綠水的變化以綠水流的變化為主，綠水流的增加為中度干旱時最大，藍水的減少量為中度干旱時最大。綠水流對城鎮(zhèn)用地擴張的響應較不敏感，城鎮(zhèn)用地擴張引起的綠水的變化以綠水儲的變化為主，綠水儲的減少量為中度干旱時最大；城鎮(zhèn)用地擴張引起的藍水變化量在豐水季隨干旱等級增強而減小，在枯水季為中度干旱時最大。