劉金平， 任艷群， 張萬昌， 陶 輝， 易 路

（1.華北水利水電大學(xué)測(cè)繪與地理信息學(xué)院，河南鄭州 450046； 2.中國科學(xué)院空天信息創(chuàng)新研究院，北京 100094； 3.中國科學(xué)院新疆生態(tài)與地理研究所荒漠與綠洲生態(tài)國家重點(diǎn)實(shí)驗(yàn)室，新疆烏魯木齊 830011； 4.西湖大學(xué)工學(xué)院浙江省海岸帶環(huán)境與資源研究重點(diǎn)實(shí)驗(yàn)室，浙江杭州 310024）

0 引言

水是一種重要的自然資源和戰(zhàn)略經(jīng)濟(jì)資源，在維系和促進(jìn)經(jīng)濟(jì)社會(huì)可持續(xù)發(fā)展方面發(fā)揮著不可替代的作用，也是大氣環(huán)流和水文循環(huán)過程中的重要因素，最易受到氣候變化的直接影響[1-4］，這種影響因土地利用等條件的不同而存在著巨大的區(qū)域差異和不確定性[5］，導(dǎo)致水資源在時(shí)間和空間上的重新分配以及水資源總量的改變，增加洪澇、干旱等極端災(zāi)害發(fā)生的頻率和強(qiáng)度[6-10］。政府間氣候變化專門委員會(huì)（Intergovernmental Panel on Climate Change，IPCC）第五次評(píng)估報(bào)告（AR5）證實(shí)了20 世紀(jì)以來全球氣候變暖的事實(shí)[11］。在全球氣候變化的驅(qū)動(dòng)下，雅魯藏布江流域以冰川強(qiáng)烈退縮、季節(jié)性積雪變化明顯、部分湖泊擴(kuò)張（如然烏湖）、土地覆被變化為主要特征的下墊面條件發(fā)生著劇烈的變化[12-13］，導(dǎo)致了流域徑流的變化[14-17］。

國內(nèi)外學(xué)者在氣候變化和下墊面變化對(duì)徑流影響方面的研究，可以概括為以下三種方法[18］：第一種是流域?qū)Ρ仍囼?yàn)法，適用于較小流域中消除氣候變化的影響[19］；第二種是水文特征參數(shù)法，利用簡單的數(shù)理統(tǒng)計(jì)方法對(duì)長時(shí)間序列的水文特征參數(shù)進(jìn)行統(tǒng)計(jì)分析，未考慮水文過程的空間異質(zhì)性和氣候變化的機(jī)制，適用于下墊面條件比較均勻和降水量空間差異不大的流域[20-21］；第三種是流域水文模型模擬法，尤其是隨著遙感技術(shù)的快速發(fā)展，利用分布式水文模型研究氣候變化或下墊面變化對(duì)徑流的影響已經(jīng)成熟[22-24］。與前兩種方法相比，流域水文模型模擬法不僅考慮了空間異質(zhì)性，而且具有更好地描述水文物理過程的優(yōu)勢(shì)，總體思路是先建立水文模型，制定不同的氣候和土地利用模擬策略進(jìn)行模擬，最后根據(jù)獲得的模擬徑流厘清氣候變化和下墊面變化對(duì)徑流的貢獻(xiàn)。盡管在這一領(lǐng)域已經(jīng)做了大量研究，由于對(duì)深層機(jī)制認(rèn)識(shí)不足，導(dǎo)致在不同流域其氣候變化和下墊面變化對(duì)徑流的影響結(jié)論未能達(dá)成一致[18］。此外，現(xiàn)有的一些研究主要集中在單一環(huán)境因子變化對(duì)徑流的影響，其中包括下墊面變化或氣候變化。針對(duì)雅魯藏布江流域的徑流變化及其影響因素研究已存在大量的工作，主要集中在對(duì)徑流演變規(guī)律進(jìn)行分析[15，25-26］和氣候變化對(duì)徑流的影響研究[6，16，27-28］。但是在雅魯藏布江流域厘清氣候和下墊面變化對(duì)徑流貢獻(xiàn)方面的研究較為缺乏，在這方面僅有的研究主要位于中游[29］和奴下站集水范圍內(nèi)[13］，未能進(jìn)行全流域研究，這一定程度上是由于該流域基礎(chǔ)數(shù)據(jù)資料匱乏、水文氣象站點(diǎn)布設(shè)不足且分配不均等[30-32］，如在流域出口處無水文觀測(cè)站點(diǎn)。另外，雅魯藏布江流域雪冰分布廣泛，積雪和冰川作為特殊的下墊面類型，易受到氣候變化的影響，有研究估計(jì)雪冰融水約占總徑流的30%以上[33］，但其作為特殊的下墊面對(duì)徑流的定量影響尚未有相關(guān)報(bào)道。因此，如何改善雅魯藏布江流域徑流影響機(jī)制的分析手段成為流域合理開發(fā)利用的一大難題，這也是目前高原寒區(qū)流域研究中亟待解決的問題之一[34］。

綜上，流域水文模型結(jié)合衛(wèi)星數(shù)據(jù)的模擬系統(tǒng)無疑為資料匱乏的高原寒區(qū)流域研究提供了一種理想的重要手段[35］。水文模型通過對(duì)水文過程的物理概化，利用計(jì)算機(jī)數(shù)值模擬方式，實(shí)現(xiàn)流域產(chǎn)匯流過程及地表徑流和地下徑流等水文過程變量的計(jì)算。尤其是近幾十年快速發(fā)展起來的集成3S技術(shù)的分布式水文模型將最新獲取的遙感下墊面信息和降水空間數(shù)據(jù)作為模型參數(shù)輸入，通過數(shù)值模擬再現(xiàn)流域內(nèi)部空間降雨-徑流形成直至最后匯至流域出口的時(shí)空過程，發(fā)展相對(duì)成熟的模型如SWAT（Soil and Water Assessment Tool）模型、TOP?MODEL 模型和SHE（System Hydrologic European）模型等。但是以上提到的模型在現(xiàn)階段由于未全面考慮寒區(qū)流域水文的冰川和凍土兩大特性，以至于不能準(zhǔn)確模擬高原寒區(qū)流域的水文過程。Zhao等[36］基于VIC（Variable Infiltration Capacity）模型框架改進(jìn)的VIC-CAS 模型雖然能夠有效的針對(duì)單條冰川進(jìn)行模擬，但是VIC 模型具有的大尺度的特性不能滿足本研究要求高空間分辨率的需求（1 km×1 km）。而滲蓄一體化動(dòng)態(tài)產(chǎn)流機(jī)制的分布式水文模 型ESSI（infiltration Excess and Saturation excess Soil-water dynamic Integration model for hydrology）由于考慮了高寒山區(qū)的下墊面特征具有較高的模擬精度，在寒旱區(qū)的黑河流域模擬月徑流的精度達(dá)0.92（Nash-Sutcliffe 效率系數(shù)）[37］。ESSI 模型是由張萬昌團(tuán)隊(duì)研建而成，以柵格為基本單元，考慮并設(shè)計(jì)了水文模擬時(shí)段內(nèi)水文年的豐、平和枯的不同情況，也具有積雪模塊，但該模型在降水聚集狀態(tài)判定及融雪量上存在一定的不足，為了提高在雅魯藏布江流域的模擬精度，需要對(duì)ESSI 模型進(jìn)行改進(jìn)。

因此，本研究以典型的高原寒區(qū)雅魯藏布江流域?yàn)檠芯繀^(qū)，基于ESSI模型并改進(jìn)其降水聚集狀態(tài)判定及融雪量計(jì)算模塊，采用動(dòng)態(tài)下墊面數(shù)據(jù)，通過制定不同的模擬策略，實(shí)現(xiàn)在資料稀缺的雅魯藏布江流域厘清氣候和下墊面變化對(duì)徑流的影響，將在三個(gè)方面展開分析：（1）氣候和下墊面變化對(duì)出水口徑流變化的貢獻(xiàn)率，（2）氣候和下墊面變化對(duì)空間產(chǎn)流變化的貢獻(xiàn)率，（3）典型下墊面類型（積雪和冰川）的徑流變化特征以及對(duì)流域總徑流的貢獻(xiàn)。本研究的開展既是區(qū)域水文水資源基礎(chǔ)性研究的需要，具有重要的理論研究意義，同時(shí)對(duì)于更科學(xué)地利用與分配該流域的水資源，以及對(duì)我國乃至亞洲的可持續(xù)發(fā)展亦具有重要的現(xiàn)實(shí)意義。

1 研究區(qū)概況

雅魯藏布江流域（圖1）在中國境內(nèi)的部分全長約2 057 km，地理范圍介于82°01′~97°06′E 和27°49′~31°17′N，總面積約為25.8×104km2。它是世界上最高的流域，平均海拔在4 600 m 以上[38］，但是流域內(nèi)海拔差異極大（149~7 159 m）。從西至東被拉孜水文站和羊村水文站劃分為上游、中游和下游，出口處無水文觀測(cè)站點(diǎn)，但在流域下游存在奴下水文站，其集水面積約占流域面積的80%。受孟加拉灣暖濕氣流、西風(fēng)環(huán)流以及高原地理環(huán)境的影響，流域內(nèi)降水分布非常不均勻，從流域下游至上游呈梯度減少趨勢(shì)[39］。流域內(nèi)氣溫具有東南—西北梯度減少的趨勢(shì)，且與海拔呈顯著負(fù)相關(guān)關(guān)系。由于獨(dú)特的地理環(huán)境，其生態(tài)環(huán)境非常復(fù)雜，植被類型多樣，但主要以高原草地為主，約占流域面積的64%。流域內(nèi)冰川條數(shù)占青藏高原冰川數(shù)量的29%，占比中國冰川數(shù)量、面積、儲(chǔ)量分別為23%、24%和23%[12］；冬季積雪大約相當(dāng)于7.4×1010m3的水量，冰川和積雪季節(jié)性融水對(duì)于下游水資源利用具有重要作用，對(duì)高原寒區(qū)生態(tài)系統(tǒng)有重要影響[40-42］。

圖1 雅魯藏布江流域基礎(chǔ)地理概況Fig.1 The fundamental geographic information of the Yarlung Zangbo River basin

2 數(shù)據(jù)與方法

2.1 數(shù)據(jù)來源

2.1.1 氣象水文數(shù)據(jù)

雅魯藏布江流域氣象觀測(cè)站點(diǎn)較少并且分布非常不均勻，特別是流域內(nèi)上游地區(qū)無可用觀測(cè)資料[39］，因此需要空間分布的陸面數(shù)據(jù)支持。全球陸面同化系統(tǒng)GLDAS（Global Land Data Assimilation Systems）是全球變化與水循環(huán)研究的重要數(shù)據(jù)源之一，分發(fā)了V1.0和V2.0兩個(gè)版本，但GLDAS V2.0在數(shù)據(jù)一致性、數(shù)據(jù)質(zhì)量季節(jié)穩(wěn)定性及對(duì)趨勢(shì)性描述能力方面則明顯優(yōu)于GLDAS V1.0 數(shù)據(jù)[43-44］。GLDAS V2.0 包括Noah、CLSM（Catchment Land Surface Model）兩個(gè)陸面過程模型輸出數(shù)據(jù)集，其中CLSM 是由美國航空航天局（NASA）全球建模和同化辦公室（GMAO）構(gòu)思并不斷發(fā)展而成的陸面模型，采用以地形獲取的流域替代傳統(tǒng)陸面模型的網(wǎng)格作為陸面單元[45］，比Noah陸面模型更適用于流域的氣候要素分析。因此，本研究采用GLDAS V2.0的CLSM 陸面模型輸出的氣溫、降水和實(shí)際蒸散發(fā)數(shù)據(jù)集，時(shí)間序列為1986—2010 年，空間分辨率為0.25°×0.25°。為了提高數(shù)據(jù)精度以及滿足驅(qū)動(dòng)分布式水文模型ESSI的要求，首先基于地面氣象站點(diǎn)觀測(cè)數(shù)據(jù)對(duì)CLSM 陸面數(shù)據(jù)進(jìn)行校正，再利用梯度距離平方反比法[46-47］（GIDS，Gradient plus Inverse Distance Squared）將數(shù)據(jù)插值成1 km×1 km 空間分辨率。

本研究采用奴下水文站1986—2010 年月尺度徑流觀測(cè)數(shù)據(jù)進(jìn)行水文模型的率定和驗(yàn)證，以獲得雅魯藏布江流域出口處的徑流量。

2.1.2 土地利用數(shù)據(jù)

土地利用數(shù)據(jù)采用中國科學(xué)院資源環(huán)境科學(xué)數(shù)據(jù)中心（http：//www. resdc. cn/）提供的中國1：10 萬比例尺土地利用現(xiàn)狀遙感監(jiān)測(cè)數(shù)據(jù)庫。該數(shù)據(jù)集已經(jīng)在國土資源調(diào)查、水文、生態(tài)研究中被廣泛使用，被認(rèn)為是我國精度最高的土地利用產(chǎn)品，包括1970年代末期（1970s，相當(dāng)于1980 年）、1980 年代末期（1980s，相當(dāng)于1990 年）、1995 年、2000 年、2005 年、2010 年和2015 年共7 期土地利用空間分布數(shù)據(jù)，本底數(shù)據(jù)是以美國陸地衛(wèi)星Landsat 遙感影像數(shù)據(jù)作為主信息源，通過人工目視解譯獲取，該數(shù)據(jù)集土地利用一級(jí)類型綜合評(píng)價(jià)精度達(dá)到94.3%以上，二級(jí)類型分類綜合精度達(dá)到91.2%以上[48］。目前，1 000 m 柵格數(shù)據(jù)可免費(fèi)下載。本研究選用1990 年、1995年、2000年和2010年的土地利用數(shù)據(jù)。

2.1.3 地形數(shù)據(jù)

本研究采用SRTM V4（Shuttle Radar Topogra?phy Mission，Version 4）DEM 數(shù)據(jù)表征流域內(nèi)地形高度。SRTM 是由美國航空航天局（NASA）、美國圖像測(cè)繪局（NIMA）、德國及意大利航天局共同實(shí)施的航天飛機(jī)雷達(dá)地形測(cè)量任務(wù)。SRTM DEM 數(shù)據(jù)包括1 000 m、90 m、30 m 三個(gè)級(jí)別的分辨率。本文采用1 000 m 空間分辨率的DEM 數(shù)據(jù)。在此基礎(chǔ)上，利用IDL編程提取DEM 各個(gè)像元對(duì)應(yīng)的經(jīng)度和緯度，用作GLDAS數(shù)據(jù)插值的輔助數(shù)據(jù)。

2.2 研究方法

2.2.1 ESSI水文模型改進(jìn)方法

由于原ESSI中融雪模塊設(shè)計(jì)較為薄弱，僅利用一個(gè)溫度閾值作為雨雪分離和融雪臨界點(diǎn)，不能反映自然界的實(shí)際觀測(cè)現(xiàn)象，包括雪、雨或雨加雪，因此，揭示了模型中應(yīng)該存在兩個(gè)溫度閾值，低于某個(gè)溫度閾值時(shí)，降水全部以降雪的形式發(fā)生；高于某個(gè)溫度閾值時(shí)，降水則全部以降雨的形式發(fā)生；介于這兩個(gè)溫度閾值中間時(shí)，則以雨、雪混合狀態(tài)發(fā)生。為了更加準(zhǔn)確地判斷降水的聚集狀態(tài)，借助遙感手段可以避免在計(jì)算融雪量時(shí)對(duì)無雪區(qū)的錯(cuò)誤估算。因此，本研究通過耦合降水聚集狀態(tài)判定、融雪量修正等方面改進(jìn)ESSI模型。

（1）降水聚集狀態(tài)判定

原ESSI模型中，降水聚集狀態(tài)是以一個(gè)簡單的溫度閾值決定的，這種設(shè)計(jì)不符合實(shí)際觀測(cè)到的自然現(xiàn)象，因此Klok等[49］提出了一種有效的模糊轉(zhuǎn)換函數(shù)（Fuzzy transition function）[50］，當(dāng)處于兩閾值之間時(shí)，降雪百分比可以用式（1）表示：

式中：Psnow為降雪占總降水的百分比；TR/S為雨、雪各占50%時(shí)的溫度；T為當(dāng)前的溫度；2Ttrans為發(fā)生雨雪混合狀態(tài)時(shí)的溫度范圍。當(dāng)前溫度小于TR/S-Ttrans時(shí)，降水全部以降雪的形式發(fā)生，當(dāng)溫度大于TR/S+Ttrans時(shí)，降水全部以降雨的形式發(fā)生。

（2）融雪量修正相關(guān)算法

根據(jù)度-日模型，融雪量Ms可用公式（2）表示：

式中：DDFsnow為積雪度日因子；T為實(shí)際溫度；T0為融雪溫度閾值。

如果實(shí)際積雪量小于理論融雪量，則正積溫先將積雪融化完后，判斷積雪層下是否有冰川覆蓋，如果有積雪覆蓋，則剩余的正積溫再利用冰川度日因子進(jìn)行計(jì)算融水量，其假設(shè)是冰川在短時(shí)間內(nèi)不會(huì)發(fā)生面積的變化。由于原ESSI 中雨雪分離閾值和融雪閾值使用同一個(gè)溫度值，實(shí)際上兩者應(yīng)該設(shè)置獨(dú)立的閾值，因此，本研究在融雪模塊中設(shè)置獨(dú)立的融雪閾值。另一方面，采用MODIS 積雪產(chǎn)品進(jìn)行修正計(jì)算融雪量時(shí)對(duì)無雪區(qū)的錯(cuò)誤估算。因此，計(jì)算融雪量分兩種情況：

1）當(dāng)遙感數(shù)據(jù)觀測(cè)到有積雪覆蓋時(shí)，融雪按照

式（2）進(jìn)行計(jì)算，然后按照公式（3）進(jìn)行判斷實(shí)際融雪量：

式中：M0為雪冰總?cè)诨?；Ms為理論融雪量；Hs為累積積雪量，需要注意的是，當(dāng)Ms≥Hs時(shí)，計(jì)算M0后要將Hs校正為0；Mice為理論融冰量，由于假設(shè)短時(shí)間內(nèi)積雪不會(huì)發(fā)生面積的變化，所以冰川足夠融化，Mice也為實(shí)際融冰量，其計(jì)算見式（4）：

2）當(dāng)遙感數(shù)據(jù)觀測(cè)為無積雪覆蓋時(shí)，融雪量為0，修正積雪累積量HS為0，正積溫如果大于0 且有冰川覆蓋時(shí)，則正積溫全部用于融冰，此時(shí)，M0=Mice，見式（5）：

根據(jù)氣候水文和土地利用變化特征，將研究時(shí)段（1986—2010 年）劃分為一個(gè)率定期（1986—1990年）和4 個(gè)驗(yàn)證期（1991—1995 年、1996—2000 年、2001—2005 年和2006—2010 年）?；谂抡?986—2010 年的逐月徑流數(shù)據(jù)，對(duì)模擬的徑流量進(jìn)行率定和驗(yàn)證。在各模擬時(shí)段內(nèi)ESSI 模型逐月徑流模擬精度及模擬效率如表1 所示，以及模擬結(jié)果和實(shí)測(cè)值對(duì)比情況見圖2。從表1 和圖2 可以看到，無論率定期還是驗(yàn)證期都取得了令人滿意的結(jié)果，R2都在0.84 以上，NSE（Nash-Sutcliffe 效率系數(shù)）都達(dá)到了0.81以上，滿足本研究的需要。

圖2 雅魯藏布江流域率定期和驗(yàn)證期月徑流量模擬值與實(shí)測(cè)值對(duì)比Fig.2 Comparison between simulated and measured monthly runoff of Yarlung Zangbo River basin during calibration（a）and validation periods（b）~（e）

表1 雅魯藏布江流域月徑流量模擬結(jié)果的精度評(píng)價(jià)Table 1 Accuracy evaluation of monthly runoff simulation results in Yarlung Zangbo River basin

2.2.2 模擬策略

雅魯藏布江流域徑流變化主要受氣候和下墊面變化影響，為了定量區(qū)分這兩種因素的變化對(duì)流域徑流變化的影響，假設(shè)它們對(duì)徑流的影響是相互獨(dú)立的，本研究設(shè)定了不同的模擬情景（表2），基于表2所示各種情景方案的設(shè)計(jì)。

表2 氣候變化和下墊面變化對(duì)徑流變化影響評(píng)估情景方案設(shè)計(jì)Table 2 Scenarios design for quantifying the impacts of climate change and human activities on runoff

2.2.3 貢獻(xiàn)率計(jì)算方法

氣候變化和下墊面變化對(duì)出水口徑流變化的貢獻(xiàn)率，是指氣候變化或下墊面變化對(duì)徑流造成的影響與兩者共同對(duì)徑流造成的影響的比值，可用式（6）和（7）表示：

式中：Cc和Cs分別為時(shí)段內(nèi)氣候和下墊面變化對(duì)徑流的貢獻(xiàn)率；Rc和Rs分別為時(shí)段內(nèi)因氣候和下墊面變化造成的徑流變化量；Rc+s為時(shí)段內(nèi)徑流的總變化量。以時(shí)段1 和時(shí)段2 為例，其技術(shù)路線可由圖3表示。

圖3 氣候變化和下墊面變化對(duì)徑流變化影響的研究路線（以時(shí)段1和時(shí)段2為例）Fig.3 Route of quantifying the contribution rates of climate change and underlying surface change（taking period 1 and period 2 for example）

3 結(jié)果分析

3.1 土地利用變化特征分析

表3 列出了雅魯藏布江流域1990 年、1995 年、2000 年和2010 年不同土地利用類型的面積以及在不同時(shí)段間土地利用類型變化情況。由此可知，在1990—1995 年時(shí)段間草地變化幅度最大且面積減少1.24×104km2，其中高覆蓋度草地減少4.87×104km2，而中覆蓋度草地和低覆蓋草地都呈擴(kuò)張趨勢(shì)且面積分別增加1.85×104km2和1.77×104km2。其次，未利用地和林地面積增加明顯。耕地（包括水田和旱地）和部分水域（河渠、湖泊和灘地）面積存在不同程度的減少。但是，冰川面積增加了847.9 km2。在1995—2000 年時(shí)段各類土地利用變化方向總體上與1990—1995年時(shí)段相反，但草地仍然是變化幅度最大的土地利用類型，具體表現(xiàn)為：草地面積增加1.2×104km2，其中高覆蓋度草地增加4.85×104km2，而中覆蓋度草地和低覆蓋度草地分別減少1.84×104km2和1.78×104km2。其次，未利用地和林地面積減少明顯。耕地和部分水域（河渠、湖泊和灘地）存在不同程度的擴(kuò)張。但冰川面積減少811.7 km2。

表3 1990年、1995年、2000年和2010年土地利用類型面積及不同時(shí)期面積變化Table 3 Area of different land use types in 1990，1995，2000 and 2010，and the associated change in different periods

在2000—2010 年間，在各類土地利用類型中，面積變化最大的仍然是草地（-26.9 km2），高覆蓋度草地、中覆蓋度草地和低覆蓋度草地面積都呈減少趨勢(shì)。耕地減少13.4 km2，而城鄉(xiāng)、工礦、居民用地和水域分別增加了25.8 km2和11.8 km2，其中冰川面積沒有發(fā)生變化。在整個(gè)土地利用變化時(shí)段內(nèi)（1990—2010 年），各類土地利用類型中變化幅度最大的是草地，在20 年間減少了138.9 km2，其中高覆蓋度草地和低覆蓋度草地分別減少了191.4 km2和92.6 km2，而中覆蓋度草地增加了145.2 km2。其次，林地、城鄉(xiāng)、工礦、居民用地和未利用地分別增加了108.1 km2、49.3 km2和12.8 km2。而耕地和水域在1990—2010年分別減少了20.1 km2和19.1 km2。

3.2 氣候和下墊面變化對(duì)出水口徑流變化的貢獻(xiàn)率

不同時(shí)段氣候變化和下墊面對(duì)徑流變化的影響以及貢獻(xiàn)率如表4所示?？偝鏊诘膹搅髯兓沂玖苏麄€(gè)雅魯藏布江流域氣候變化和下墊面變化對(duì)徑流變化的綜合影響。結(jié)果可知：從時(shí)段1 至?xí)r段2，氣候變化導(dǎo)致雅魯藏布江流域總出水口徑流增加了31.5%，時(shí)段2 至?xí)r段3 氣候變化導(dǎo)致徑流減少了13.1%，時(shí)段3 至?xí)r段4 氣候變化導(dǎo)致徑流增加20.3%，時(shí)段1 至?xí)r段4 氣候變化導(dǎo)致徑流增加36.3%；下墊面變化引起的徑流變化比例從時(shí)段1 到時(shí)段2、時(shí)段2 到時(shí)段3 逐漸增加，分別為2.3%、7.7%，而時(shí)段3 至?xí)r段4 下墊面變化引起徑流減少12.3%，從時(shí)段1 至?xí)r段4，下墊面變化引起的徑流變化比例僅為0.2%。因此，可以看出氣候變化引起徑流變化的比例（絕對(duì)值都大于13%）都高于下墊面變化引起徑流變化的比例（絕對(duì)值都小于13%）。從貢獻(xiàn)率來看，氣候變化對(duì)徑流變化的貢獻(xiàn)在4 個(gè)時(shí)段中一直處于主導(dǎo)地位，這是由于雅魯藏布江流域獨(dú)特的地理環(huán)境和復(fù)雜的地形，人類活動(dòng)較弱，下墊面變化較小，而高海拔更易受氣候變化的影響，因此，與氣候變化對(duì)徑流的影響相比，下墊面變化對(duì)徑流變化的貢獻(xiàn)非常小。

表4 不同時(shí)段間氣候和下墊面變化對(duì)徑流變化的貢獻(xiàn)率Table 4 Contribution rates of climate change and underlying surface change to the runoff variation between different time periods

3.3 氣候和下墊面變化對(duì)空間產(chǎn)流變化的貢獻(xiàn)率

為了定量分析氣候變化和下墊面變化對(duì)空間產(chǎn)流變化的影響，將分別研究時(shí)段1 至?xí)r段2、時(shí)段2 至?xí)r段3、時(shí)段3 至?xí)r段4、時(shí)段1 至?xí)r段4 的氣候變化和下墊面變化對(duì)空間產(chǎn)流變化的貢獻(xiàn)率。

（1）時(shí)段1至?xí)r段2

圖4 顯示了氣候變化和下墊面變化在時(shí)段1（1991—1995年）至?xí)r段2（1996—2000年）對(duì)空間產(chǎn)流變化的貢獻(xiàn)率。從整個(gè)流域來看，氣候變化對(duì)產(chǎn)流變化的貢獻(xiàn)率最大，氣候變化貢獻(xiàn)率[圖4（a）］大于75%的區(qū)域面積占流域總面積的96.9%，流域中游貢獻(xiàn)率最大，下游東北部區(qū)域貢獻(xiàn)率較小，而土地利用變化貢獻(xiàn)率[圖4（b）］在該區(qū)域貢獻(xiàn)率較大，表明從時(shí)段1 至?xí)r段2 該區(qū)主要受土地利用影響，而其他區(qū)域土地利用貢獻(xiàn)率相對(duì)較小。總體來講，從時(shí)段1至?xí)r段2，雅魯藏布江流域產(chǎn)流變化的主要貢獻(xiàn)來源于氣候變化。

圖4 在時(shí)段1（1991—1995年）和時(shí)段2（1996—2000年）期間氣候變化和下墊面變化對(duì)雅魯藏布江流域產(chǎn)流變化的貢獻(xiàn)率Fig.4 Contribution rate of climate change and underlying surface change to runoff change in the Yarlung Zangbo River basin during period 1（1991—1995）and 2（1996—2000）：climate change（a）；underlying surface change（b）

（2）時(shí)段2至?xí)r段3

圖5（a）和5（b）分別顯示了氣候變化和下墊面變化在時(shí)段2（1996—2000年）至?xí)r段3（2001—2005年）對(duì)空間產(chǎn)流變化的貢獻(xiàn)率。氣候變化對(duì)產(chǎn)流變化的貢獻(xiàn)率較大，氣候變化貢獻(xiàn)率大于75%的區(qū)域面積占流域總面積的92.8%，流域上游的氣候變化貢獻(xiàn)率最大，中游其次，下游最低。下墊面變化對(duì)產(chǎn)流的貢獻(xiàn)率主要表現(xiàn)為對(duì)上游產(chǎn)流量的負(fù)貢獻(xiàn)最大。

圖5 在時(shí)段2（1996—2000年）和時(shí)段3（2001—2005年）期間氣候變化和下墊面變化對(duì)雅魯藏布江流域產(chǎn)流變化的貢獻(xiàn)率Fig.5 Contribution rate of climate change and underlying surface change to runoff change in the Yarlung Zangbo River basin during period 2（1996—2000）and 3（2001—2005）：climate change（a）；underlying surface change（b）

（3）時(shí)段3至?xí)r段4

圖6 顯示了氣候變化和下墊面變化在時(shí)段3（2001—2005年）至?xí)r段4（2006—2010年）對(duì)空間產(chǎn)流變化的貢獻(xiàn)率。由圖可知，總體上，與下墊面變化相比，氣候變化對(duì)產(chǎn)流的貢獻(xiàn)率仍然是較大，與前兩個(gè)時(shí)段相比，時(shí)段3至?xí)r段4氣候變化的貢獻(xiàn)率在上、中、下游分布相對(duì)均勻。氣候變化貢獻(xiàn)率[圖6（a）］大于75%的區(qū)域面積占流域總面積的97.5%，下墊面變化貢獻(xiàn)率[圖6（b）］非常小，集中分布在-25%~25%范圍內(nèi)，這是由于從2001—2010 年下墊面變化（2000年和2010年土地利用數(shù)據(jù)）非常小。

圖6 在時(shí)段3（2001—2005年）和時(shí)段4（2006—2010年）期間氣候變化和下墊面變化對(duì)雅魯藏布江流域產(chǎn)流變化的貢獻(xiàn)率Fig.6 Contribution rate of climate change and underlying surface change to runoff change in the Yarlung Zangbo River basin during period 3（2001—2005）and 4（2006—2010）：climate change（a）；underlying surface change（b）

（4）時(shí)段1至?xí)r段4

圖7 顯示了氣候變化和下墊面變化在時(shí)段1（1991—1995年）至?xí)r段4（2006—2010年）對(duì)空間產(chǎn)流變化的貢獻(xiàn)率。從圖7（a）來看，在1991—2010年期間，除下游東北部區(qū)域外，氣候變化對(duì)徑流變化的貢獻(xiàn)率在流域內(nèi)更加均勻，貢獻(xiàn)率大于75%的區(qū)域面積占流域總面積的97.9%。從圖7（b）來看，下墊面變化對(duì)徑流的影響主要是流域下游東北部區(qū)域（多為冰川積雪發(fā)育區(qū)域），其他區(qū)域?qū)搅髯兓呢暙I(xiàn)非常小?？傮w上講，在1991—2010 年，雅魯藏布江流域氣候變化對(duì)產(chǎn)流的貢獻(xiàn)占主導(dǎo)地位，而下墊面變化對(duì)產(chǎn)流的貢獻(xiàn)主要位于下游雪冰發(fā)育區(qū)域，但遠(yuǎn)小于氣候變化對(duì)徑流的影響。這是由于該流域人類活動(dòng)較小，下墊面變化較弱，尤其在2000年之后，下墊面變化幅度更?。ū?）。

圖7 在時(shí)段1（1991—1995年）和時(shí)段4（2006—2010年）期間氣候變化和下墊面變化對(duì)雅魯藏布江流域產(chǎn)流變化的貢獻(xiàn)率Fig.7 Contribution rate of climate change and underlying surface change to runoff change in the Yarlung Zangbo River basin during period 1（1991—1995）and 4（2006—2010）：climate change（a）；underlying surface change（b）

3.4 典型下墊面類型的徑流變化特征以及對(duì)流域總徑流的貢獻(xiàn)

雅魯藏布江流域典型下墊面即為積雪和冰川，研究和分析雪冰融水徑流對(duì)流域徑流的貢獻(xiàn)具有十分重要的實(shí)際意義和科研價(jià)值。改進(jìn)后的ESSI具有積雪模塊和冰川、積雪融化模塊，能夠依據(jù)氣溫、降水和蒸散發(fā)等氣候條件的變化對(duì)季節(jié)性積雪和常年積雪、冰川進(jìn)行模擬。因此，冰雪融水徑流可以由設(shè)定融雪和融冰模塊是否參與計(jì)算得到，即首先按照率定好的參數(shù)，在融雪和融冰模塊參與條件下模擬雅魯藏布江流域的徑流（記為R1），然后在融雪和融冰模塊不參與條件下模擬該流域的徑流（記為R2），則R1 與R2 的差值表示冰雪融水徑流量。

圖8 為雅魯藏布江流域1991—2010 年模擬年徑流量與雪冰融水徑流量的年際變化趨勢(shì)。從圖8可知，隨著氣溫的升高，雪冰融水徑流呈增加的變化趨勢(shì)，對(duì)年徑流的平均貢獻(xiàn)率在21.1%~48.6%，最大年份為2005 年，最小年份為1999 年，多年平均為33.6%。在年內(nèi)變化上（圖9），流域雪冰融水一般從4 月開始增大，8 月達(dá)到最大為4 621.1 m3·s-1，10 月消融期結(jié)束。11 月至次年3 月為積累期，消融量較小。從雪冰融水徑流占總徑流的百分比來看，全年變化相對(duì)平穩(wěn)（10.4%~38.9%），總體來講，冷季時(shí)雪冰融水徑流占總徑流百分比較低，暖季時(shí)所占百分比較高。

圖8 雅魯藏布江流域總徑流、雪冰融水徑流以及雪冰融水徑流占總徑流百分比的年際變化趨勢(shì)Fig.8 Interannual trends of total runoff，snow-ice melt runoff，and snow-ice melt runoff as a percentage of total runoff in the Yarlung Zangbo River basin

圖9 雅魯藏布江流域總徑流、雪冰融水徑流以及雪冰融水徑流占總徑流百分比的年內(nèi)變化情況Fig.9 Annual variations of total runoff，snow-ice melt runoff，and snow-ice melt runoff as a percentage of total runoff in the Yarlung Zangbo River basin

4 討論

通過對(duì)氣候和下墊面變化對(duì)出水口徑流變化的貢獻(xiàn)率、氣候和下墊面變化對(duì)空間產(chǎn)流變化的貢獻(xiàn)率、典型下墊面類型（積雪和冰川）的徑流變化特征以及對(duì)流域總徑流的貢獻(xiàn)這三個(gè)方面展開研究，本研究較為全面地分析了整個(gè)雅魯藏布江流域氣候和下墊面變化對(duì)徑流的影響。由于雅魯藏布江流域在流域上游和下游氣象水文觀測(cè)資料稀缺，尤其是流域出口處無水文觀測(cè)站，導(dǎo)致針對(duì)該流域與徑流變化有關(guān)的研究主要集中于中游和奴下站以上的區(qū)域，如Liu 等[13］利用水量平衡法分析了1974—2000 年雅魯藏布江流域奴下站以上4 個(gè)區(qū)域（分別以拉孜、奴各沙、羊村和奴下水文站以分界點(diǎn)）的氣候和下墊面變化對(duì)徑流的影響，結(jié)果表明在冷季和暖季，氣候和下墊面對(duì)徑流的影響具有差異，但總體來說，氣候變化較下墊面變化對(duì)徑流的影響更為顯著；雪冰融水增加了流域徑流量，這與本研究的結(jié)論是一致的。王蕊[29］利用SWAT水文模型以雅魯藏布江流域中游（拉孜和羊村水文站之間的集水區(qū)域）為研究區(qū)，分析了該區(qū)域1990—2010 年氣候和下墊面變化對(duì)徑流的影響，結(jié)果表明雪冰融水對(duì)流域中游總徑流的平均貢獻(xiàn)率約為6.2%。而本研究以整個(gè)雅魯藏布江流域?yàn)檠芯繀^(qū)，得到雪冰融水對(duì)流域總徑流的多年平均貢獻(xiàn)率為33.6%，這與許朋琨等[33］的描述是一致的。其原因在于積雪和冰川主要集中在雅魯藏布江流域下游，而中游分布較少[39，51］，導(dǎo)致流域中游雪冰融水對(duì)總徑流的貢獻(xiàn)較低。

雖然本研究針對(duì)雅魯藏布江整個(gè)流域首次定量給出了氣候和下墊面對(duì)徑流（包括總徑流、空間產(chǎn)流）的影響，但是也存在一些限制，如未對(duì)上、中、下游進(jìn)行具體的分析，以及未具體地闡明氣候因子和除雪冰類型以外的下墊面因子在不同時(shí)段對(duì)徑流的影響，這將會(huì)在下一步的工作中進(jìn)一步分析。盡管在研究實(shí)施過程中力爭水文模擬精確，但是由于氣象水文資料稀缺和模型不確定性都會(huì)導(dǎo)致分析結(jié)果存在一定的不確性，即模型物理結(jié)構(gòu)、參數(shù)、驅(qū)動(dòng)和初始條件的不確定性都可能通過復(fù)雜的誤差傳遞反映到模擬結(jié)果中。在研究中發(fā)現(xiàn)土地利用數(shù)據(jù)在該流域可能存在一定的誤差，如在1990—2010 年間冰川這一類型的面積增加了36.2 km2，可能導(dǎo)致在一定程度上增加了雪冰徑流計(jì)算誤差。然而本研究盡可能地將誤差減少到最小，如改進(jìn)水文模型、利用GLDAS 和站點(diǎn)數(shù)據(jù)等，在典型高原寒區(qū)雅魯藏布江流域通過制定不同模擬策略厘清氣候和下墊面變化對(duì)徑流的影響，可為資料稀缺的高原寒區(qū)流域的徑流變化機(jī)制分析提供一個(gè)新的途徑，這些分析結(jié)果為了解雅魯藏布江流域徑流變化主要影響因素提供了一個(gè)獨(dú)特的機(jī)會(huì)，尤其是定量評(píng)估了雪冰徑流在總徑流的占比，為山區(qū)流域春、夏季融雪和融冰洪水預(yù)估提供了技術(shù)支撐。

5 結(jié)論

本研究以位于青藏高原東南部的典型高原寒區(qū)雅魯藏布江流域?yàn)檠芯繀^(qū)，通過耦合降水聚集狀態(tài)判定、融雪量修正等方面改進(jìn)ESSI 模型，對(duì)雅魯藏布江流域進(jìn)行水文模擬，制定不同的模擬策略以定量評(píng)估流域1991—2010 年不同時(shí)段間氣候和下墊面變化對(duì)總徑流和空間產(chǎn)流的貢獻(xiàn)率，尤其定量分析了作為流域典型下墊面（積雪和冰川）的融水徑流的變化特征以及對(duì)總徑流的貢獻(xiàn)率。主要研究工作及結(jié)論如下：

（1）為了定量分析氣候和下墊面變化對(duì)流域徑流變化的影響，在1991—2010 年設(shè)計(jì)了8 個(gè)模擬情景，得到在不同的時(shí)間段內(nèi)徑流對(duì)氣候變化和下墊面變化的響應(yīng)，結(jié)果表明，在1991—2010年期間，由于下墊面的變化，總體上引起徑流具有增加的趨勢(shì)。在全球變暖的背景下，雅魯藏布江流域氣候變化引起徑流的變化量遠(yuǎn)大于下墊面引起徑流的變化量，除在2001—2005年期間由于氣候的變化使徑流量減少外，其他時(shí)間段由于氣候變化都使徑流量增加，總體上，氣候變化使雅魯藏布江流域徑流量增加，且遠(yuǎn)高于下墊面變化對(duì)徑流的影響。

（2）4 個(gè)時(shí)段中（時(shí)段1 至?xí)r段2、時(shí)段2 至?xí)r段3、時(shí)段3 至?xí)r段4、時(shí)段1 至?xí)r段4），不同時(shí)段間氣候和下墊面變化對(duì)徑流變化的貢獻(xiàn)率差異較大，總體上，氣候變化對(duì)徑流變化的貢獻(xiàn)率都遠(yuǎn)高于下墊面變化對(duì)徑流變化的貢獻(xiàn)率，即氣候變化主導(dǎo)了徑流變化的方向。從空間上看，在4個(gè)時(shí)段中，氣候變化對(duì)流域產(chǎn)流的貢獻(xiàn)率在上游和中游都較大，在下游東北部的貢獻(xiàn)率較小，而在該區(qū)域下墊面變化的貢獻(xiàn)率較大，主要是由于雪冰的融化使該區(qū)的下墊面類型發(fā)生了變化。

（3）1991—2010 年期間，隨著氣溫的升高，雪冰融水徑流呈增加的變化趨勢(shì)，最小和最大年份分別出現(xiàn)在1999 年和2005 年，對(duì)年徑流的平均貢獻(xiàn)率在21.1%~48.6% 范圍內(nèi)，多年平均貢獻(xiàn)率為33.6%。在年內(nèi)變化上，雪冰融水徑流一般從4 月份開始增大，8 月份達(dá)到最大，10 月份達(dá)到消融末期；從年內(nèi)來看，從1 月至12 月份，雪冰融水徑流對(duì)徑流的貢獻(xiàn)率在全年變化相對(duì)平緩，冷季時(shí)相對(duì)較低，而暖季時(shí)達(dá)到最高。