黃河源區(qū)水文要素變化及驅(qū)動(dòng)因素分析

張 鳴，朱 奎，魯 帆，戴雁宇，宋昕熠，于嵩彬

（1.中國(guó)礦業(yè)大學(xué) 資源與地球科學(xué)學(xué)院，江蘇 徐州 221008；2.中國(guó)水利水電科學(xué)研究院 水資源研究所，北京 100038；3.長(zhǎng)沙理工大學(xué) 水利與環(huán)境工程學(xué)院，湖南 長(zhǎng)沙 410114）

0 引言

黃河是我國(guó)僅次于長(zhǎng)江的第二大河，多年平均水資源總量為647 億m3，不到長(zhǎng)江的7%，但水資源開(kāi)發(fā)利用率高達(dá)80%，遠(yuǎn)超40%的生態(tài)警戒線，水資源短缺是黃河流域最大的矛盾［1－3］。黃河源區(qū)唐乃亥水文站的多年平均徑流量占黃河水資源總量的34.5%［4］，黃河源區(qū)的水資源狀況不僅關(guān)系當(dāng)?shù)氐慕?jīng)濟(jì)社會(huì)發(fā)展和生態(tài)環(huán)境保護(hù)，而且對(duì)整個(gè)黃河流域水資源開(kāi)發(fā)利用有著重要影響。近年來(lái)，已有不少學(xué)者對(duì)黃河源區(qū)徑流變化進(jìn)行了研究，蘇中海等［5］、劉希勝等［6］認(rèn)為黃河源區(qū)上游徑流量呈現(xiàn)不顯著增大趨勢(shì)，而劉晶等［7］、Yan 等［8］認(rèn)為黃河源區(qū)徑流量整體呈減小趨勢(shì)。各學(xué)者研究選取的時(shí)間段和數(shù)據(jù)不同，結(jié)論也有所不同。對(duì)于徑流驅(qū)動(dòng)因素，商瀅等［9］采用雙累積曲線得出非降水因素是徑流變化的主要驅(qū)動(dòng)因素的結(jié)論，張士鋒等［10］通過(guò)建立驅(qū)動(dòng)模型得出降水、潛在蒸散發(fā)對(duì)徑流變化分別起正向、負(fù)向驅(qū)動(dòng)作用的結(jié)論。上述學(xué)者研究徑流時(shí)多將研究期看作一個(gè)整體，較少考慮研究期內(nèi)徑流變化過(guò)程，驅(qū)動(dòng)因素分析大都只考慮部分氣象要素和人類活動(dòng)的影響，較少綜合考慮氣象要素和下墊面的影響并定量分析。黃河源區(qū)位于素有“亞洲水塔”和“第三極”之稱的青藏高原，對(duì)氣候變化較敏感，在全球變暖的環(huán)境下以往的研究結(jié)論已不完全適用，需補(bǔ)充研究論證。

目前采用Budyko 法、水文模型、雙累積曲線法均能定量揭示氣候變化和人類活動(dòng)對(duì)徑流的影響。采用Budyko 法便于從氣候變化中分離出不同影響因子。近年來(lái)，國(guó)內(nèi)外學(xué)者采用Budyko 法建立多種經(jīng)驗(yàn)公式并在長(zhǎng)江和黃河等流域廣泛應(yīng)用［11－15］，但傳統(tǒng)Budyko法下墊面參數(shù)在不同時(shí)段會(huì)出現(xiàn)跳躍式的變化，這與實(shí)際情況存在一定偏差。本文基于1959—2019 年黃河源區(qū)水文氣象數(shù)據(jù)，分析黃河源區(qū)水文要素演變規(guī)律，運(yùn)用累積距平曲線和Mann－Kendall 突變檢驗(yàn)法確定徑流突變年份，采用變參數(shù)的Budyko 公式定量分離降水、潛在蒸散發(fā)和下墊面變化對(duì)徑流量變化的貢獻(xiàn)率，以期為黃河流域合理開(kāi)發(fā)利用水資源、建設(shè)水利工程以及防汛抗旱提供參考。

1 研究區(qū)概況和數(shù)據(jù)來(lái)源

1.1 研究區(qū)概況

黃河流域位于東經(jīng)96°—119°、北緯32°—42°，發(fā)源于青藏高原巴顏喀拉山北麓的約古宗列盆地，流域地勢(shì)西高東低，西部河源地區(qū)平均高程在4 000 m 以上，源區(qū)涉及青海、四川、甘肅3 個(gè)省的6 個(gè)州18 個(gè)縣。黃河源區(qū)唐乃亥水文站控制集水面積為12.2 萬(wàn)km2，黃河源區(qū)是世界上高原水體與濕地集中分布的地區(qū)之一。研究區(qū)內(nèi)水系、湖泊、10 個(gè)氣象站和唐乃亥水文站分布見(jiàn)圖1，為增加站網(wǎng)密度，選取的部分氣象站位于集水面積之外。

圖1 黃河源區(qū)范圍

1.2 數(shù)據(jù)來(lái)源

使用的唐乃亥水文站數(shù)據(jù)源自《黃河水文年鑒》《黃河水資源公報(bào)》《黃河泥沙公報(bào)》，10 個(gè)氣象站數(shù)據(jù)源自中國(guó)氣象數(shù)據(jù)網(wǎng)，1982—2019 年歸一化植被指數(shù)（精度為5 km）源自國(guó)家地球系統(tǒng)科學(xué)數(shù)據(jù)中心，1980—2020 年土地利用類型遙感監(jiān)測(cè)數(shù)據(jù)（精度為1 km）源自中國(guó)科學(xué)院資源環(huán)境科學(xué)與數(shù)據(jù)中心。

2 研究方法

2.1 趨勢(shì)檢驗(yàn)和突變檢驗(yàn)

趨勢(shì)檢驗(yàn)采用線性回歸分析法和Mann－Kendall檢驗(yàn)法（簡(jiǎn)稱M－K 檢驗(yàn)）［16］。線性回歸分析法可以反映變量隨時(shí)間的變化趨勢(shì)。M－K 檢驗(yàn)是一種非參數(shù)檢驗(yàn)法，其優(yōu)點(diǎn)是不需要樣本遵從一定的分布規(guī)律，也不受少數(shù)異常值的干擾。

突變檢驗(yàn)采用M－K 檢驗(yàn)法和累積距平法，累積距平法是根據(jù)時(shí)間序列累積值與均值的差值繪制累積距平曲線，能比較直觀地判斷樣本的變化趨勢(shì)和突變年份。

2.2 潛在蒸散發(fā)量與植被覆蓋度計(jì)算

采用FAO Penman－Monteith 方法計(jì)算流域潛在蒸散發(fā)量：

式中：E0為潛在蒸散發(fā)量，Tmean為日均氣溫，Rn為凈輻射，G為土壤熱通量，γ為干濕表常數(shù)，u2為2 m 高風(fēng)速，es、ea分別為飽和、實(shí)際水氣壓，Δ為飽和水氣壓曲線斜率。

李苗苗等［17］基于像元二分模型建立植被覆蓋度計(jì)算公式：

式中：fc為植被覆蓋度，N為歸一化植被指數(shù)，Nsoil為裸土或者無(wú)植被覆蓋區(qū)域的歸一化植被指數(shù)，Nveg為完全被植被覆蓋區(qū)域的歸一化植被指數(shù)。

由于存在噪聲等其他因素的影響，因此取5%的N值為Nsoil，取95%的N值為Nveg。為更好地反映研究區(qū)內(nèi)植被變化情況，歸一化植被指數(shù)可選取年內(nèi)最大值。

2.3 徑流歸因分析

基于Budyko 假設(shè)的長(zhǎng)時(shí)間序列中，水分輸入可看作降水量P，能量輸入可看作潛在蒸散發(fā)量E0，考慮水量平衡和能量平衡的水熱耦合方程為

式中：E為流域蒸散發(fā)量。

對(duì)于閉合流域來(lái)說(shuō)，長(zhǎng)期水量平衡方程為

式中：Q為徑流量，ΔS為流域水分蓄變量。

考慮到區(qū)域氣候、地理環(huán)境等的差異，國(guó)內(nèi)外學(xué)者基于Budyko 理論陸續(xù)建立了各種考慮下墊面參數(shù)的經(jīng)驗(yàn)公式［18－21］，本文選用Zhang 等［18］建立的經(jīng)驗(yàn)公式：

式中：ε為干旱指數(shù)，ε＝E0/P；n為下墊面參數(shù)。

采用式（5）分析徑流變化時(shí)，下墊面參數(shù)和水文要素會(huì)在不同時(shí)段呈現(xiàn)跳躍性變化，這與實(shí)際情況存在偏差。Jiang 等［15］提出了變參數(shù)的Budyko 公式，選取氣象和人類活動(dòng)等要素為協(xié)變量構(gòu)建下墊面參數(shù)的相關(guān)函數(shù)，假定選取m個(gè)自變量，滑動(dòng)平均序列為y1，y2，…，yt，第i（i＝1，2，…，t）年的m個(gè)自變量為，，…，，可選用一次函數(shù)、二次函數(shù)、指數(shù)函數(shù)等擬合下墊面參數(shù)和m個(gè)自變量之間的關(guān)系，第i年的徑流量Qi為

式中：Pi為第i年降水量，εi為第i年干旱指數(shù)，βi為第i年回歸系數(shù)，ei為第i年殘差。

假定ei符合正態(tài)分布N（0，σ2），則Qi符合正態(tài)分布，根據(jù)參數(shù)β和σ的似然函數(shù)求出極大似然估計(jì)量和，對(duì)應(yīng)有ni＝Y(jié)i，ni為第i年的下墊面參數(shù)。

為定量描述氣象因子和下墊面參數(shù)對(duì)徑流量變化的影響，引入彈性系數(shù)φ來(lái)表示徑流量對(duì)某因子x變化的敏感程度：

假設(shè)變化期和基準(zhǔn)期內(nèi)降水量、潛在蒸散發(fā)量、下墊面參數(shù)變化分別為ΔP、ΔE0、Δn，根據(jù)基準(zhǔn)期信息推求變化期內(nèi)徑流變化量：

3 水文要素變化

基于唐乃亥水文站實(shí)測(cè)徑流數(shù)據(jù)，分析1959—2019 年黃河源區(qū)年徑流量變化情況，見(jiàn)圖2（a）。M－K趨勢(shì)檢驗(yàn)結(jié)果顯示，整體上年徑流量呈現(xiàn)減小趨勢(shì)，線性減小速率約為2.36 億m3/10 a，多年實(shí)測(cè)年均徑流量為203.9 億m3，1959 年年徑流量（為328.08 億m3）和2002 年年徑流量（為105.74 億m3）分別為61 a樣本系列中徑流量的最大值和最小值。利用M－K 檢驗(yàn)方法對(duì)徑流序列進(jìn)行突變檢驗(yàn)，見(jiàn)圖2（b），在0.05顯著性水平下，UBK曲線與UFK曲線在1987—1990 年置信區(qū)間內(nèi)多次交匯。為明確徑流突變年份，引入年徑流量累積距平曲線，見(jiàn)圖2（c），年徑流量在1959—1989 年整體呈現(xiàn)上升趨勢(shì)，在1989—2019 年整體呈現(xiàn)下降趨勢(shì)，在1989 年發(fā)生一次明顯突變，綜合判斷1989 年為徑流突變年份。為進(jìn)一步反映近年來(lái)全球變暖環(huán)境下黃河源區(qū)的徑流變化，繪制1990—2019 年年徑流量累積距平曲線，見(jiàn)圖2（d），2004 年為第2 個(gè)突變年份。總體上，從1959—1989 年、1990—2004 年和2005—2019 年3 個(gè)時(shí)段來(lái)看，年徑流量呈現(xiàn)先下降后上升趨勢(shì)。

圖2 1959—2019 年黃河源區(qū)徑流量變化

對(duì)1959—2019 年黃河源區(qū)年內(nèi)徑流量進(jìn)行分析，1981 年9 月月徑流量最大（為92.02 億m3），2003 年1月月徑流量最?。?.33 億m3），每年7 月月均徑流量最大（為35.04 億m3），每年1 月月均徑流量最?。?.58 億m3），7—9 月徑流量占全年總徑流量的46.97%，而1—3 月徑流量?jī)H占全年總徑流量的7.32%，由此可見(jiàn)，黃河源區(qū)徑流量季節(jié)分配不均。徑流集中度指月徑流量以矢量形式累積，其各分矢量之和占年徑流量的百分比，反映年內(nèi)徑流量的集中程度。集中期反映全年徑流量集中的月份［22］。唐乃亥水文站年內(nèi)徑流集中度和集中期計(jì)算結(jié)果見(jiàn)圖3，集中度為35%～55%，整體上呈現(xiàn)下降趨勢(shì)，下降率為0.716%/10 a。集中期為7—9 月，整體上呈現(xiàn)略微提前的趨勢(shì)，對(duì)于唐乃亥水文站和下游水庫(kù)來(lái)說(shuō)，洪峰有提前的趨勢(shì)。

圖3 唐乃亥水文站年內(nèi)徑流集中度和集中期計(jì)算結(jié)果

基于10 個(gè)氣象站降水?dāng)?shù)據(jù)分析黃河源區(qū)1959—2019 年降水量變化趨勢(shì)，見(jiàn)圖4。1959—2019 年多年平均降水量為504.83 mm，1962 年降水量（為398.24 mm）和2018 年降水量（為674.21 mm）分別為61 a 樣本系列中的最小值和最大值。M－K 趨勢(shì)檢驗(yàn)結(jié)果顯示，統(tǒng)計(jì)量Z為3.72＞2.56，降水量呈現(xiàn)顯著增加趨勢(shì)（顯著性水平α＝0.01），平均線性增加速率為9.1 mm/10 a。

圖4 1959—2019 年降水量變化

徑流系數(shù)代表流域產(chǎn)流能力，綜合反映流域內(nèi)自然地理要素對(duì)降水—徑流關(guān)系的影響。1959—2019年徑流系數(shù)變化情況見(jiàn)圖5，1959—1989 年、1990—2004 年和2005—2019 年對(duì)應(yīng)的徑流系數(shù)均值分別為0.36、0.28 和0.32，徑流系數(shù)呈現(xiàn)先減小后增大的趨勢(shì)，與徑流量變化趨勢(shì)相同。

圖5 1959—2019 年徑流系數(shù)變化情況

4 驅(qū)動(dòng)因素

4.1 降水、氣溫和潛在蒸散發(fā)分析

基于黃河源區(qū)10 個(gè)氣象站數(shù)據(jù)對(duì)降水量進(jìn)行Kriging 插值分析，得出多年平均降水量的空間分布，見(jiàn)圖6（a），黃河源區(qū)1959—2019 年多年平均降水量在空間上呈現(xiàn)由西北向東南遞增的趨勢(shì)。M－K 突變檢驗(yàn)結(jié)果［見(jiàn)圖6（b）］和降水量累積距平曲線［見(jiàn)圖6（c）］顯示，2002 年和2016 年為降水序列的突變點(diǎn)。1959—2002 年、2003—2016 年和2017—2019 年多年平均降水量分別為490.67、525.97、613.84 mm，尤其是近年來(lái)降水量增加趨勢(shì)愈發(fā)明顯。

圖6 降水量變化

基于10 個(gè)氣象站實(shí)測(cè)數(shù)據(jù)分析黃河源區(qū)1959—2019 年年均氣溫變化趨勢(shì)，見(jiàn)圖7（a）。多年平均氣溫為－0.66 ℃，1965 年年均氣溫最低（為－2.02 ℃），2016 年年均氣溫最高（為0.73 ℃），M－K 趨勢(shì)檢驗(yàn)統(tǒng)計(jì)量Z為11.40＞2.56，氣溫呈現(xiàn)顯著上升趨勢(shì)（α＝0.01），平均線性上升速率為0.3 ℃/10 a。由氣溫累積距平曲線［見(jiàn)圖7（b）］可以看出，氣溫累積距平值在1997 年發(fā)生轉(zhuǎn)折，結(jié)合趨勢(shì)檢驗(yàn)結(jié)果，氣溫在1959—1997 年呈現(xiàn)不顯著上升趨勢(shì)，在1997 年之后開(kāi)始顯著上升。

圖7 1959—2019 年氣溫變化趨勢(shì)與氣溫累積距平曲線

基于10 個(gè)氣象站實(shí)測(cè)氣溫、風(fēng)速、日照和相對(duì)濕度數(shù)據(jù)，利用FAO Penman－Monteith 方法計(jì)算各氣象站潛在蒸散發(fā)量，采用Thiessen 多邊形方法計(jì)算出黃河源區(qū)年平均潛在蒸散發(fā)量，見(jiàn)圖8（a）。1959—2019年多年平均潛在蒸散發(fā)量為724.8 mm，1965 年潛在蒸散發(fā)量最小（為672.4 mm），2010 年潛在蒸散發(fā)量最大（為784.5 mm）。M－K 趨勢(shì)檢驗(yàn)統(tǒng)計(jì)量Z為7.89＞2.56，潛在蒸散發(fā)量呈現(xiàn)顯著增加趨勢(shì)（α＝0.01），平均線性增加率為8.2 mm/10 a。潛在蒸散發(fā)量主要受日照時(shí)數(shù)、風(fēng)速、氣溫和相對(duì)濕度的影響，因此利用多元線性回歸法對(duì)每個(gè)氣象站的潛在蒸散發(fā)量進(jìn)行歸因分析。結(jié)果表明，氣溫是90%氣象站潛在蒸散發(fā)量的主要影響因素，風(fēng)速是10%氣象站潛在蒸散發(fā)量的主要影響因素。M－K 突變檢驗(yàn)結(jié)果［見(jiàn)圖8（b）］顯示，潛在蒸散發(fā)量在1997 年發(fā)生一次突變，突變時(shí)間點(diǎn)與氣溫的相同，說(shuō)明氣溫與潛在蒸散發(fā)量之間具有較強(qiáng)的相關(guān)性。

圖8 1959—2019 年潛在蒸散發(fā)量變化趨勢(shì)與M－K 突變檢驗(yàn)

為分析1959—2019 年黃河源區(qū)降水量、氣溫、潛在蒸散發(fā)量與徑流量的相關(guān)性，對(duì)其進(jìn)行單因素相關(guān)性檢驗(yàn)，結(jié)果顯示，徑流量與降水量之間的皮爾遜系數(shù)為0.778，具有極強(qiáng)的相關(guān)性，遠(yuǎn)遠(yuǎn)大于氣溫的和潛在蒸散發(fā)量的。

4.2 下墊面變化分析

下墊面主要指與大氣直接接觸的地球表面，本文主要分析黃河源區(qū)植被、冰雪和凍土的變化情況。植被覆蓋度是衡量地表植被變化的重要指標(biāo)，黃河源區(qū)植被覆蓋度變化趨勢(shì)與累積距平曲線見(jiàn)圖9。M－K 趨勢(shì)檢驗(yàn)統(tǒng)計(jì)量Z為2.01＜2.56，植被覆蓋度呈現(xiàn)不顯著上升趨勢(shì)（α＝0.01），這與管曉祥等［23］、劉啟興等［24］、韓思淇等［25］的研究結(jié)果一致。在空間上植被覆蓋度呈現(xiàn)由西北向東南遞增的趨勢(shì)，與降水量的空間分布一致，這與高思琦等［26］認(rèn)為黃河源區(qū)植被的空間分布主要受降水影響的觀點(diǎn)一致。

為分析黃河源區(qū)內(nèi)不同植被類型的變化情況，對(duì)1980—2020 年土地利用類型數(shù)據(jù)進(jìn)行分析，見(jiàn)圖10。源區(qū)內(nèi)土地利用類型主要為草地，城鎮(zhèn)用地和湖泊面積變化不大，而林地、草地分別在1990—1995 年、1996—2000 年呈現(xiàn)略微退化趨勢(shì)。2000—2013 年未利用土地面積呈現(xiàn)顯著減小趨勢(shì)，而林地面積和草地面積呈現(xiàn)顯著增大趨勢(shì)。植被變化是人類活動(dòng)和氣候因子共同作用的結(jié)果，在生態(tài)環(huán)境較為敏感的青藏高原地區(qū)植被變化更為顯著。氣候變化對(duì)植被的影響主要通過(guò)降水量和氣溫的變化改變土壤濕度和溫度，調(diào)節(jié)植物生長(zhǎng)過(guò)程中酶的活性，影響植物光合作用和呼吸作用［27－29］。降水量增加會(huì)促使高寒草原和荒漠植被生長(zhǎng)，增加植被面積，氣溫升高則促使植被覆蓋度提高［30］，植被的改善會(huì)增強(qiáng)區(qū)域蒸散發(fā)，進(jìn)而增加降水量。

圖10 1980—2020 年土地利用類型變化

歸一化植被指數(shù)是表征區(qū)域植被變化的重要指數(shù)，1982—2019 年黃河源區(qū)降水量與歸一化植被指數(shù)的相關(guān)性分析及顯著性檢驗(yàn)結(jié)果見(jiàn)圖11。78%的區(qū)域降水量對(duì)歸一化植被指數(shù)起正向促進(jìn)作用，兩者成顯著正相關(guān)的區(qū)域主要集中在降水相對(duì)貧乏的西北部地區(qū)和東北部地區(qū)；21%的區(qū)域降水量與歸一化植被指數(shù)負(fù)相關(guān)，兩者成顯著負(fù)相關(guān)的區(qū)域很少，植被類型以林地和耕地為主，降水較充沛，這些地區(qū)的植被變化可能受人類活動(dòng)影響較大，這與史丹丹等［31］的研究結(jié)果一致。管曉祥等［23］研究發(fā)現(xiàn)黃河源區(qū)氣溫與歸一化植被指數(shù)顯著正相關(guān)，且春季的影響最大。黃河源區(qū)多年平均降水量為504.83 mm，降水較充沛，因此水量不是限制植被生長(zhǎng)的主要條件。黃河源區(qū)高程普遍在4 000 m以上，熱量不足限制植被生長(zhǎng)，因此氣溫是影響黃河源區(qū)植被變化的主要因子［23，27－28，32］。此外，黃河源區(qū)的人類活動(dòng)也對(duì)植被變化產(chǎn)生一定影響，如過(guò)度砍伐和放牧、退耕還林、減畜工程和成立三江源國(guó)家自然保護(hù)區(qū)等。

圖11 1982—2019 年降水量與歸一化植被指數(shù)相關(guān)性分析

黃河源區(qū)冰雪面積變化情況見(jiàn)圖12，累積距平曲線顯示，1982—1997 年冰雪面積呈現(xiàn)增加趨勢(shì)，1997—2015 年冰雪呈現(xiàn)消融趨勢(shì)，冰雪面積變化的突變時(shí)間為1997 年，與氣溫的一致，表明氣溫是控制冰雪變化的主要因素［33］。

圖12 1982—2015 年冰雪面積占比累積距平曲線

關(guān)于黃河源區(qū)凍土的變化情況，曹慧宇［34］的研究表明2003—2019 年多年凍土面積由3.59 萬(wàn)km2減小到3.42 萬(wàn)km2，4.82%的多年凍土轉(zhuǎn)化為季節(jié)性凍土，季節(jié)性凍土最大凍結(jié)深度減小16 cm，活動(dòng)層厚度增加9.28 cm。馮雨晴［35］的研究表明1961—2015 年黃河源區(qū)凍土融化夾層厚度和活動(dòng)層厚度呈現(xiàn)增加趨勢(shì)，1997 年之后呈現(xiàn)顯著增大趨勢(shì)。多年凍土的退化以及活動(dòng)層厚度的增加導(dǎo)致土壤滲透性增強(qiáng)，地表徑流更多入滲成為地下水，土壤含水量整體增加，蒸散發(fā)量增加，流域地下水儲(chǔ)量增加，從而減少夏季徑流量，增加冬季徑流量和基流量，調(diào)節(jié)徑流量的季節(jié)性分布［36－40］。此外，活動(dòng)層底部的凍結(jié)層上水可能側(cè)向補(bǔ)給江河湖泊，黃河源區(qū)地下冰融水對(duì)地表徑流的補(bǔ)給比例為13.2%～16.7%，但其融水對(duì)徑流的改變較小［41－44］。

4.3 Budyko 假設(shè)定量分離

根據(jù)徑流突變結(jié)果，以1989 年和2004 年為突變年份將研究期分為3 個(gè)時(shí)段，3 個(gè)時(shí)段的多年平均徑流深、降水量、潛在蒸散發(fā)量、氣溫見(jiàn)表1。

表1 黃河源區(qū)各時(shí)段徑流深和協(xié)變量

為消除氣候變化及地下水儲(chǔ)量對(duì)水量平衡的影響，構(gòu)建下墊面參數(shù)時(shí)滑動(dòng)平均步長(zhǎng)取9 a，利用式（5）計(jì)算滑動(dòng)平均后的下墊面參數(shù)n?；贘iang等［15］提出的變參數(shù)Budyko 公式，選取降水量、潛在蒸散發(fā)量和氣溫作為協(xié)變量，選用線性函數(shù)進(jìn)行擬合，下墊面參數(shù)n與降水量、潛在蒸散發(fā)量和氣溫的關(guān)系式為

式（9）標(biāo)準(zhǔn)化系數(shù)顯示降水量是下墊面參數(shù)的主要影響因素，其次是氣溫。根據(jù)式（5）和式（4）利用氣象數(shù)據(jù)和徑流數(shù)據(jù)推求氣候變化造成的徑流變化量dQc，再根據(jù)dQh＝dQ－dQc計(jì)算分析人類活動(dòng)對(duì)徑流的影響（dQh為人類活動(dòng)對(duì)徑流的改變量，dQ為徑流總改變量）。結(jié)果顯示，在第1 個(gè)突變期（1959—2004年）氣候變化對(duì)徑流量減小的貢獻(xiàn)率為86.63%、人類活動(dòng)對(duì)徑流量減小的貢獻(xiàn)率為13.37%，在第2 個(gè)突變期（1990—2019 年）氣候變化對(duì)徑流量增加的貢獻(xiàn)率為101.18%、人類活動(dòng)對(duì)徑流量增加的貢獻(xiàn)率為－1.18%。在整個(gè)研究期內(nèi)氣候變化是徑流量變化的主導(dǎo)因素，人類活動(dòng)對(duì)徑流量的增大一直起負(fù)向作用，但對(duì)徑流量變化的影響程度明顯降低。

為定量分離不同氣象因子和下墊面變化的貢獻(xiàn)率，分別計(jì)算降水量、潛在蒸散發(fā)量和下墊面參數(shù)在不同時(shí)期對(duì)應(yīng)的彈性系數(shù)，結(jié)果見(jiàn)表2?；鶞?zhǔn)期和變化期（分別以1959—1989 年為基準(zhǔn)期、1990—2004 年為變化期，1990—2004 年為基準(zhǔn)期、2005—2019 年為變化期）降水量的彈性系數(shù)始終大于潛在蒸散發(fā)量的和下墊面參數(shù)的，說(shuō)明徑流對(duì)降水的響應(yīng)最為敏感。

表2 黃河源區(qū)下墊面參數(shù)和各協(xié)變量對(duì)應(yīng)的彈性系數(shù)

根據(jù)式（8）計(jì)算突變期各因子對(duì)徑流量變化的貢獻(xiàn)率，在第1 個(gè)突變期內(nèi)下墊面參數(shù)n對(duì)徑流量變化的貢獻(xiàn)率達(dá)到60.69%，是徑流量減小的主要驅(qū)動(dòng)因素，降水其次（為31.75%），潛在蒸散發(fā)的貢獻(xiàn)率最?。?.56%）；第2 個(gè)突變期內(nèi)降水對(duì)徑流量變化的貢獻(xiàn)率達(dá)到103.91%，是徑流量增加的主要驅(qū)動(dòng)因素，潛在蒸散發(fā)其次（為－10.37%），下墊面變化的貢獻(xiàn)率最小（為6.46%）。此外，計(jì)算得出2 個(gè)突變期內(nèi)各因子引起的模擬徑流變化量之和與實(shí)測(cè)徑流變化量的誤差分別為3.98%和0.1%，誤差較小，表明結(jié)果具有準(zhǔn)確性。

5 結(jié)論

1）黃河源區(qū)1959—2019 年徑流整體呈現(xiàn)先下降后上升的變化趨勢(shì)，多年實(shí)測(cè)年均徑流量為203.9 億m3，年內(nèi)徑流集中度為35%～55%，集中期位于7—9 月，年內(nèi)徑流量季節(jié)分配不均，1989 年和2004 年為徑流量突變年份。1959—2019 年年降水量和潛在蒸散發(fā)量都呈現(xiàn)顯著增加趨勢(shì)，氣溫在1997 年之后顯著上升。1982—2019 年植被覆蓋度呈現(xiàn)不顯著上升趨勢(shì)，在空間上植被覆蓋度呈現(xiàn)由西北向東南遞增的趨勢(shì)，與降水量的空間分布一致。

2）基于變參數(shù)的Budyko 方法定量分離不同因子對(duì)徑流量變化的影響，氣候變化是黃河源區(qū)1959—2019 年徑流量變化的主要影響因素，徑流對(duì)降水的響應(yīng)最為敏感。1959—2004 年下墊面參數(shù)n對(duì)徑流量變化的貢獻(xiàn)率達(dá)到60.69%，是徑流減小的主要驅(qū)動(dòng)因素，降水其次，潛在蒸散發(fā)的貢獻(xiàn)率最小。1990—2019年降水對(duì)徑流量變化的貢獻(xiàn)率達(dá)到103.91%，是徑流增大的主要驅(qū)動(dòng)因素，潛在蒸散發(fā)次之，下墊面的貢獻(xiàn)率最小。

3）徑流量與降水量之間的皮爾遜相關(guān)系數(shù)為0.778，遠(yuǎn)遠(yuǎn)大于氣溫的和潛在蒸散發(fā)量的，因此1959—2019 年黃河源區(qū)徑流量改變最主要的影響因素是降水，其次為氣溫。