姚春艷，劉洪鵠，劉 競(jìng)

（1.中國(guó)科學(xué)院教育部水土保持與生態(tài)環(huán)境研究中心，陜西 楊凌 712100；2.中國(guó)科學(xué)院水利部水土保持研究所，陜西 楊凌 712100；3.中國(guó)科學(xué)院大學(xué)，北京 100049；4.長(zhǎng)江水利委員會(huì)長(zhǎng)江科學(xué)院，湖北 武漢 430010）

近年來(lái)，隨著氣候變化和人類活動(dòng)的加劇，極端氣候事件頻發(fā)，使得流域徑流、泥沙等水文要素及其輸移過(guò)程發(fā)生改變［1］。國(guó)內(nèi)外許多學(xué)者對(duì)尼羅河［2］、密西西比河［3］、黃河［4］和長(zhǎng)江［5］等重要河流開(kāi)展了河流水沙情勢(shì)變化規(guī)律及歸因分析研究。其中，氣候要素中降水對(duì)河流水沙的影響更為顯著［5-6］，人類活動(dòng)則是通過(guò)水土保持工程等的實(shí)施減少流域水土流失，進(jìn)而影響水沙變化［5-7］。不可忽略的是，流域水土保持工程和水庫(kù)建設(shè)是近年來(lái)造成河流輸沙量下降的重要原因［1］，因此，開(kāi)展流域水土流失與水沙規(guī)律的綜合研究對(duì)流域水土保持和泥沙治理工作十分重要。

青藏高原是亞洲眾多河流的發(fā)源地，被稱為“亞洲水塔”［5］，對(duì)全球氣候變化的響應(yīng)極為敏感［8］。長(zhǎng)江源區(qū)地處青藏高原腹地，是我國(guó)重要的水源涵養(yǎng)地與生態(tài)安全屏障區(qū)［8］，區(qū)內(nèi)凍土、沼澤、湖泊和高寒濕地等廣泛分布，生態(tài)環(huán)境極為脆弱［9］。近年來(lái)，隨著全球氣候變暖的持續(xù)作用，源區(qū)內(nèi)水土流失量明顯增加，徑流量和泥沙等要素隨之發(fā)生變化［10］，進(jìn)而對(duì)下游河道和流域管理產(chǎn)生影響［5］。已有研究表明，源區(qū)內(nèi)徑流的季節(jié)性變化主要受降水的影響［11］，土地利用/覆蓋和冰雪融水對(duì)徑流變化也有一定調(diào)節(jié)作用［12］。水土流失和入河泥沙量增長(zhǎng)的主要原因是降水的增加和頻繁的人類活動(dòng)［13］，而長(zhǎng)時(shí)間序列內(nèi)植被變化則與流域內(nèi)泥沙變化相關(guān)性不大［14］。Teng等［15］和Wang等［16］基于模型模擬結(jié)果發(fā)現(xiàn)降水對(duì)源區(qū)水土流失影響較大，且未來(lái)水土流失量將會(huì)繼續(xù)增加。由此可見(jiàn)，降水是引發(fā)源區(qū)土壤侵蝕和水沙變化的主要驅(qū)動(dòng)力［15］，一方面直接影響坡面土壤的侵蝕，進(jìn)而影響泥沙過(guò)程［16-17］，另一方面通過(guò)影響徑流變化間接影響河流攜沙能力［18］。盡管對(duì)長(zhǎng)江源區(qū)徑流和泥沙變化及影響因素進(jìn)行了多項(xiàng)研究［11,14］，但是大/中流域尺度中侵蝕產(chǎn)沙與輸移的綜合研究仍較少，同時(shí)，氣候變化是影響源區(qū)土壤侵蝕與河流水沙變化的重要因素。因此，有必要開(kāi)展長(zhǎng)江源區(qū)侵蝕與水沙變化及其與氣候變化耦合的響應(yīng)關(guān)系研究，以準(zhǔn)確把握長(zhǎng)江上游泥沙來(lái)源及流域侵蝕產(chǎn)沙對(duì)降雨在長(zhǎng)時(shí)間尺度下的響應(yīng)關(guān)系，為源區(qū)生態(tài)工程實(shí)施效果的研究與評(píng)價(jià)提供科學(xué)依據(jù)。

因此，本文基于長(zhǎng)江源區(qū)1980—2020 年降水、徑流和泥沙等資料，研究長(zhǎng)江源區(qū)年土壤侵蝕量、徑流量和泥沙變化特征，闡明長(zhǎng)江源區(qū)土壤侵蝕變化規(guī)律，河道水沙變化對(duì)地表侵蝕的響應(yīng)，以及長(zhǎng)時(shí)間尺度下泥沙變化對(duì)侵蝕和降水的響應(yīng)關(guān)系，以期為推動(dòng)小流域科學(xué)治理和長(zhǎng)江流域高質(zhì)量發(fā)展提供理論依據(jù)。

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

1.1 研究區(qū)概況

長(zhǎng)江源區(qū)為直門達(dá)水文站控制以上的流域（圖1），地理位置位于90°～98°E、31°～36°N，總面積約為15.7×104km2。源區(qū)地勢(shì)西高東低，平均海拔超過(guò)4000 m，太陽(yáng)輻射較強(qiáng)，日溫差大但年溫差?。?6］，多年降水量顯著增加（P<0.05），年均約為370 mm，水系包括通天河水系、沱沱河水系，呈扇狀分布。主要土壤類型為高山草原土和高山草甸土，其中高山草甸土分布面積最大，約占源區(qū)總面積63%，植被類型從東南到西北依次分布有高寒草甸、高寒草原、高寒荒漠，凍土面積約占青藏高原多年凍土面積的20%［18］。人口以藏族居民為主，人口密度不到1人·km-2，牧業(yè)是區(qū)域經(jīng)濟(jì)發(fā)展的主要生產(chǎn)方式。

圖1 長(zhǎng)江源區(qū)示意圖Fig.1 Basic information of the headwaters of the Yangtze River

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

長(zhǎng)江源區(qū)直門達(dá)水文站1980—2020 年逐年徑流及泥沙數(shù)據(jù)由青海省水文局提供。從國(guó)家氣象科學(xué)數(shù)據(jù)中心（http://data.cma.cn）下載長(zhǎng)江源區(qū)及其周圍共12 個(gè)氣象站1980—2020 年逐日降水資料。長(zhǎng)江源區(qū)高程數(shù)據(jù)（空間分辨率為1 km）、長(zhǎng)江源區(qū)歸一化植被指數(shù)（Normalized Difference Vegeta?tion Index，NDVI）數(shù)據(jù)（1998—2018年，空間分辨率為1 km）和1980年、1990年、1995年、2000年、2005年、2010年、2015年和2020年8期土地利用數(shù)據(jù)（圖2）來(lái)源于資源環(huán)境科學(xué)數(shù)據(jù)中心（https://www.resdc.cn）。

圖2 長(zhǎng)江源區(qū)1980—2020年土地利用變化Fig.2 Landuse changes in the headwaters of the Yangtze River during 1980-2020

1.3 研究方法

1.3.1 土壤侵蝕量 利用中國(guó)土壤流失方程（CSLE）［19］估算長(zhǎng)江源區(qū)1980—2020 年的土壤侵蝕量。在ArcGIS 10.2的支持下，對(duì)各個(gè)因子賦值疊加計(jì)算，得到土壤侵蝕量。模型如公式（1）所示：

式中：A是土壤侵蝕模數(shù)（t·hm-2·a-1）；R是降雨侵蝕力因子（MJ·mm·hm-2·h-1·a-1）；K是土壤可蝕性（t·hm2·h·hm-2·MJ-1·mm-1）；L是坡長(zhǎng)因子，S是坡度因子，二者統(tǒng)稱為地形因子，無(wú)量綱；B是生物措施因子；E是工程措施因子；T是耕作措施因子，B、E、T因子均無(wú)量綱（取值范圍為0～1）。

（1）降雨侵蝕力因子（R）

基于章文波等［20］建立的降雨侵蝕力算法模型計(jì)算年降雨侵蝕力因子值（圖3），

圖3 長(zhǎng)江源區(qū)1980—2020年年降水量(a)和降雨侵蝕力(b)變化Fig.3 Variations of annual average precipitation(a)and rainfall erosivity factor(b)in the headwaters of Yangtze River during 1980-2020

式中：R為半月周期的降雨侵蝕力（MJ·mm·hm-2·h-1·a-1）；k為一個(gè)半月周期（每個(gè)月的前15 d 為一個(gè)半月周期，該月的其余天數(shù)為另一個(gè)半月周期）；Pj是一個(gè)半月周期內(nèi)超過(guò)12 mm 的日降水量（mm）；α和β分別為模型的待定參數(shù)，計(jì)算方程式如公式（3）和公式（4）：

式中：Pd12代表日雨量超過(guò)12 mm 的日降雨量（mm）；Py12代表一年內(nèi)超過(guò)12 mm 的日雨量的和（mm）。

（2）土壤可蝕性因子（K）

土壤可蝕性因子（K）是評(píng)價(jià)土壤被降雨分離、沖蝕和搬運(yùn)難易程度的一項(xiàng)指標(biāo)［21］。源區(qū)土壤可蝕性因子?xùn)鸥駭?shù)據(jù)集是基于1∶1 000 000 中國(guó)土壤數(shù)據(jù)庫(kù)、傳統(tǒng)的“土壤發(fā)生分類”系統(tǒng)和Nomo 圖法建立的K值與土壤性質(zhì)之間的Nomo 模型得出。長(zhǎng)江源區(qū)1000 m 分辨率土壤可蝕性K值來(lái)源于國(guó)家地球系統(tǒng)科學(xué)數(shù)據(jù)共享服務(wù)平臺(tái)（http://www.geoda?ta.cn）。利用ArcGIS 10.2進(jìn)行重采樣和插值得到本研究土壤可蝕性因子分布（圖4a）。

圖4 長(zhǎng)江源區(qū)K值(a)及LS值(b)空間分布Fig.4 Spatial distribution of K(a)and LS value(b)

（3）地形因子（LS）

基于DEM 數(shù)據(jù)計(jì)算地形因子［22］，得到LS值分布圖（圖4b）。方程如下：

式中：LS是地形因子；l是累計(jì)坡長(zhǎng)（m）；β是坡度（°）；m是坡度變量，當(dāng)坡度≤1°時(shí)，m=0.2；當(dāng)坡度為1°～3°時(shí)，m=0.3；當(dāng)坡度為3.5°～4.5°時(shí)，m=0.4；當(dāng)坡度為≥5°時(shí)，m=0.5。

（4）生物措施因子（B）

生物措施因子（B）主要反映植被覆蓋對(duì)土壤侵蝕的影響，其賦值范圍為0～1，越接近0，表示防治越好，幾乎不發(fā)生侵蝕；越接近1，表示未采取任何措施。具體包括退耕還林、還草等。本研究基于8 期土地利用數(shù)據(jù)，結(jié)合Liu 等［19］和Fu 等［23］提出的賦值原則（表1），利用ArcGIS 10.2 重分類得到B值的空間分布（圖5）。

表1 不同土地利用類型的B值Tab.1 B value of different landuse types

圖5 不同年份生物措施因子值空間分布Fig.5 Spatial distribution of biological practice factor at the different years

（5）工程措施因子（E）和耕作措施因子（T）

工程措施因子（E）反映了實(shí)施的工程措施對(duì)于流域內(nèi)水土流失的影響。長(zhǎng)江源區(qū)居民將旱地改為石坎水平梯田，通過(guò)查表可得石坎水平梯田E值為0.121。耕作措施（E）通過(guò)提高土壤抗蝕性，進(jìn)而防治土壤侵蝕。長(zhǎng)江源區(qū)屬于北部中高原半干旱喜涼作物一熟區(qū)［24］，根據(jù)土地利用類型圖將研究區(qū)耕地賦值為0.121 和0.488。其他區(qū)域分配為1，得到工程措施因子和耕作措施因子分布。

由于長(zhǎng)江源區(qū)土地利用變化較?。?5］，因此，各期土地利用可代表相鄰5～10 a 變化，基于CSLE 模型，估算得到長(zhǎng)江源區(qū)1980—2020 年土壤侵蝕模數(shù)，根據(jù)《土壤侵蝕分類分級(jí)標(biāo)準(zhǔn)》（SL190-2007）將長(zhǎng)江源區(qū)土壤侵蝕模數(shù)劃分為6個(gè)等級(jí)（0～5 t·hm-2·a-1為微度侵蝕，5～25 t·hm-2·a-1為輕度侵蝕，25～50 t·hm-2·a-1為中度侵蝕，50～80 t·hm-2·a-1為強(qiáng)烈侵蝕，80～150 t·hm-2·a-1為極強(qiáng)烈侵蝕，>150 t·hm-2·a-1為劇烈侵蝕），由于微度侵蝕被定義為小于土壤容許流失量，因此本文僅統(tǒng)計(jì)輕度及以上侵蝕強(qiáng)度等級(jí)的面積占比［19］。

1.3.2 偏最小二乘結(jié)構(gòu)方程（PLS-SEM）偏最小二乘結(jié)構(gòu)方程（PLS-SEM）是一種結(jié)合因子分析和回歸的非參數(shù)方法［26］，通過(guò)迭代估計(jì)最大限度地解釋獨(dú)立潛變量的變化。本研究利用PLS-SEM分析，通過(guò)與多個(gè)潛在變量的相關(guān)性來(lái)解釋降水、土壤侵蝕與水沙變化耦合模式。本文所提出的結(jié)構(gòu)模型由降水量（P）、徑流量（Q）、輸沙量（S）和土壤侵蝕量（M）組成，使用SmartPLS 3.0軟件繪制4個(gè)潛變量的影響路徑，通過(guò)300 次迭代執(zhí)行PLS 算法以獲得路徑系數(shù)（β），并使用5000次抽樣進(jìn)行bootstrapping以測(cè)試路徑系數(shù)的顯著性。當(dāng)模型擬合中NFI>0.9，SRMR<0.12 時(shí)，表明PLS-SEM 模型具有預(yù)測(cè)有效性，能較強(qiáng)地解釋因子間相互影響的變化。潛變量間路徑是否通過(guò)檢驗(yàn)的主要依據(jù)為顯著性P值，當(dāng)P<0.05時(shí)路徑檢驗(yàn)通過(guò)。

1.3.3 突變分析方法 利用累積距平法和滑動(dòng)t檢驗(yàn)法［27］判斷長(zhǎng)江源區(qū)1980—2020年徑流量、含沙量和輸沙量序列的突變情況。

2 結(jié)果與分析

2.1 土壤侵蝕的變化特征

圖6 展示了長(zhǎng)江源區(qū)1980—2020 年各侵蝕強(qiáng)度面積占比變化。如圖所示，輕度侵蝕面積約占總侵蝕面積的70.27%，強(qiáng)烈及以上強(qiáng)度侵蝕面積約占總侵蝕面積的10.27%。長(zhǎng)江源區(qū)1980—2020 年土壤侵蝕模數(shù)呈顯著增加趨勢(shì)（P<0.05）（圖7），多年平均土壤侵蝕模數(shù)為4.71 t·hm-2·a-1。最大值為8.16 t·hm-2·a-1，出現(xiàn)在1985年；最小值為2.13 t·hm-2·a-1，出現(xiàn)在2015 年。土壤侵蝕模數(shù)變化規(guī)律與降雨侵蝕力變化規(guī)律一致（圖3b）。

圖6 長(zhǎng)江源區(qū)1980—2020年土壤侵蝕強(qiáng)度面積占比變化Fig.6 Changes of soil loss intensity area proportion in the headwaters of Yangtze River during 1980-2020

圖7 長(zhǎng)江源區(qū)1980—2020年土壤侵蝕模數(shù)變化Fig.7 Change of annual soil loss modulus in the headwaters of Yangtze River during 1980-2020

圖8展示了長(zhǎng)江源區(qū)1980年、1990年、1995年、2000 年、2005 年、2010 年、2015 年和2020 年的土壤侵蝕強(qiáng)度的空間分布狀況。從中發(fā)現(xiàn)源區(qū)1980—2020 年土壤侵蝕強(qiáng)度整體較低，以輕度侵蝕為主，侵蝕強(qiáng)度較大的地方主要分布在東南部，如玉樹(shù)縣。東南部降水量相對(duì)較多，且地形多以高山為主，坡面處極容易發(fā)生水力侵蝕，導(dǎo)致水土流失情況相對(duì)其他區(qū)域更為嚴(yán)重。如圖7 和圖8 所示，長(zhǎng)江源區(qū)土壤侵蝕時(shí)空變化較為復(fù)雜。1980 年源區(qū)東南部土壤侵蝕面積較大且侵蝕強(qiáng)度較高，1990—1995 年侵蝕強(qiáng)度與侵蝕面積呈現(xiàn)減小趨勢(shì)。土壤侵蝕模數(shù)降低至3.82 t·hm-2·a-1，強(qiáng)烈及以上侵蝕所占面積明顯減少，土壤侵蝕得到一定遏制。然而，2000 年后侵蝕面積和侵蝕強(qiáng)度均較高且土壤侵蝕模數(shù)逐漸增加，輕度侵蝕面積下降，中度及以上侵蝕面積占比逐漸增加。2010 年后侵蝕量增加速度有所放緩，其中，2015 年土壤侵蝕面積最小且整體侵蝕強(qiáng)度最低，這與2015 年降雨侵蝕力較小有關(guān)，值僅為125.5 MJ·mm·hm-2·h-1·a-1（圖3b）。

圖8 長(zhǎng)江源區(qū)1980—2020年土壤侵蝕模數(shù)空間分布Fig.8 Spatial distribution of soil loss modulus in the headwaters of Yangtze River during 1980-2020

2.2 年徑流量、含沙量和輸沙量的變化規(guī)律

直門達(dá)水文站1980—2020年年徑流量、年含沙量和年輸沙量的變化趨勢(shì)如圖9 所示，直門達(dá)水文站多年徑流量顯著增加（P<0.05），多年平均徑流量為142.03×108m3；年徑流量最大值在2009 年，為245.14×108m3，最小值為81.95×108m3，出現(xiàn)在1997年。相比于90 年代，21 世紀(jì)以來(lái)徑流量明顯增加，增幅約為79.7%。多年含沙量和輸沙量均無(wú)顯著變化趨勢(shì)，多年平均含沙量為0.69 kg·m-3，最大值出現(xiàn)在1989年，為1.41 kg·m-3，最小值出現(xiàn)在1994年，為0.23 kg·m-3。多年平均輸沙量為1063.41×104t；輸沙量與含沙量變化較為一致，最大值均出現(xiàn)在1989年，為2980×104t，最小值均出現(xiàn)在1994年，為203×104t。

圖9 長(zhǎng)江源區(qū)1980—2020年徑流量(a)、含沙量(b)和輸沙量(c)的變化趨勢(shì)Fig.9 Variations of runoff(a),suspend sediment concentrate(b),and sediment load(c)in the headwaters of Yangtze River during 1980-2020

2.3 降水、土壤侵蝕和水沙的耦合關(guān)系

基于1980—2020 年年降水量（P）、年土壤侵蝕量（M）、年徑流量（Q）和年輸沙量（S）數(shù)據(jù)，通過(guò)PLS-SEM 分析，得到長(zhǎng)江源區(qū)1980—2020 年降水、土壤侵蝕與水沙變化耦合關(guān)系圖（圖10）。計(jì)算得知，模型NFI=0.982（NFI>0.9），SRMR=0.03（SRMR<0.12），表明模型對(duì)潛變量的擬合較好。1980—2020年長(zhǎng)江源區(qū)各潛變量之間的主要關(guān)系如下：（1）以降水為驅(qū)動(dòng)因素，降水對(duì)土壤侵蝕和河道徑流均具有直接的顯著正向效應(yīng)［β=0.727（P<0.01）和β=0.837（P<0.01）］，并進(jìn)一步通過(guò)徑流和侵蝕產(chǎn)沙間接影響輸沙量變化。（2）徑流和土壤侵蝕對(duì)輸沙量變化具有直接正向效應(yīng)［β=0.885（P<0.01）和β=0.062］，徑流對(duì)輸沙量變化的貢獻(xiàn)度明顯大于土壤侵蝕。通過(guò)間接效應(yīng)分析發(fā)現(xiàn)，降水對(duì)輸沙量具有顯著正向影響。降水、徑流和土壤侵蝕對(duì)輸沙量變化的解釋度達(dá)72.5%，其影響路徑主要為降水?徑流?輸沙（表2）。

表2 1980—2020年長(zhǎng)江源區(qū)降水、土壤侵蝕、徑流和輸沙量路徑效應(yīng)Tab.2 Route effects of precipitation,soil loss,runoff,and sediment load in the headwaters of Yangtze River during 1980-2020

圖10 長(zhǎng)江源區(qū)1980—2020年降水、土壤侵蝕與水沙變化耦合關(guān)系Fig.10 Coupling relationship between precipitation,soil loss,runoff,and sediment load in the headwaters of Yangtze River during 1980-2020

3 討論

3.1 長(zhǎng)江源區(qū)土壤侵蝕對(duì)泥沙變化的影響

本研究利用CSLE模型估算1980—2020年長(zhǎng)江源區(qū)多年平均土壤侵蝕模數(shù)的結(jié)果與Liu等［28］和魏夢(mèng)美等［29］的分析結(jié)果較為一致。但是，顯著高于Teng 等［15］估算的2.76 t·hm-2·a-1。產(chǎn)生差異的主要原因是計(jì)算相關(guān)因子值的數(shù)據(jù)來(lái)源及估算方法不同，特別是R因子和K因子。由于所選氣象站點(diǎn)數(shù)量和位置的不同，計(jì)算得到的降雨侵蝕力具有一定差異，本研究中年均R值（269.25 MJ·mm·hm-2·h-1·a-1）與于文竹［30］研究的降雨侵蝕力范圍（25.72～448.35 MJ·mm·hm-2·h-1·a-1）略有不同，但在可接受變化范圍內(nèi)。同時(shí)，區(qū)域土壤侵蝕量和降雨侵蝕力的空間變化均為由西北向東南增加。另外，本文對(duì)K值柵格數(shù)據(jù)集進(jìn)行重采樣及空間插值后，發(fā)現(xiàn)源區(qū)侵蝕主要集中在草甸土，與全國(guó)水土流失與生態(tài)安全考察結(jié)果一致，這也使得降雨入滲量降低進(jìn)而導(dǎo)致徑流量增加［31］。

源區(qū)多年平均土壤侵蝕模數(shù)呈顯著增加趨勢(shì)，然而，多年輸沙量和含沙量無(wú)明顯變化趨勢(shì)。一方面，可能是源區(qū)水力侵蝕主要分布在源區(qū)東南部，該區(qū)域以高山為主，且LS值較高，進(jìn)而使泥沙容易在坡面淤積；另一方面，雖然新的土壤侵蝕仍在增加，但源區(qū)水土保持工作取得了一定成效，源區(qū)輸沙量整體處于一個(gè)相對(duì)平衡的狀態(tài)［9］。另外，20 世紀(jì)80 年代源區(qū)人為無(wú)序開(kāi)采沙金導(dǎo)致源區(qū)的侵蝕面積和侵蝕強(qiáng)度都較大［32］，整體侵蝕情況較為嚴(yán)重；20 世紀(jì)90 年代源區(qū)植被生態(tài)發(fā)展趨勢(shì)良好，沒(méi)有出現(xiàn)退化區(qū)域［33］，生態(tài)環(huán)境總體趨勢(shì)良好，因此，水土流失狀況得到明顯改善；2000—2004 年，近70%的地區(qū)出現(xiàn)植被輕微退化的現(xiàn)象［32］，且長(zhǎng)江源區(qū)降雨侵蝕力增加，導(dǎo)致影響水土流失的外部動(dòng)力條件增強(qiáng)［34］，進(jìn)一步加劇了水土流失狀況。但是，得益于近些年三江源自然生態(tài)保護(hù)規(guī)劃的實(shí)施［35］，2010年后侵蝕增加速度有所放緩（圖7）。

盡管近年來(lái)源區(qū)進(jìn)行了大量的生態(tài)保護(hù)與恢復(fù)工作，但是對(duì)比前期水土流失遙感調(diào)查數(shù)據(jù)發(fā)現(xiàn)，土壤侵蝕仍在加劇，對(duì)生態(tài)環(huán)境的威脅依然存在，長(zhǎng)江源區(qū)的水土保護(hù)及生態(tài)修復(fù)工作仍需繼續(xù)開(kāi)展。

3.2 不同時(shí)期降水、土壤侵蝕與水沙的響應(yīng)關(guān)系的差異

由1980—2020 年長(zhǎng)江源區(qū)年均徑流量累積距平圖可知（圖11），2004年是年徑流量累積距平值發(fā)生顯著變化的一年，在2004 年達(dá)到最低，隨后升高。結(jié)合滑動(dòng)t檢驗(yàn)法，進(jìn)一步確定了2004 年為徑流量的突變點(diǎn)。年輸沙量累積距平值在1988 年和2006 年發(fā)生變化，但兩個(gè)突變點(diǎn)均未通過(guò)檢驗(yàn)，年含沙量在1994年發(fā)生的突變也未通過(guò)檢驗(yàn)。因此，根據(jù)徑流量突變點(diǎn)將水文序列劃分為1980—2004年、2005—2020年兩個(gè)階段。

圖11 長(zhǎng)江源區(qū)1980—2020年徑流量累積距平變化Fig.11 Cumulative anomaly changes of runoff in the headwaters of Yangtze River during 1980-2020

源區(qū)直門達(dá)水文站輸沙量較為穩(wěn)定，徑流量在2004年出現(xiàn)拐點(diǎn)，與劉彥等［11］和張永勇等［36］的研究結(jié)果一致。徑流量發(fā)生突變的原因主要是21 世紀(jì)初降水的增加［37］與三江源生態(tài)保護(hù)項(xiàng)目的實(shí)施［34］。因此，將2004 年作為分界點(diǎn)，對(duì)1980—2004 年、2005—2020年兩個(gè)時(shí)期進(jìn)行耦合關(guān)系分析可知（圖12），降水對(duì)徑流量變化的解釋度上升［β=0.690（P<0.01）和β=0.869（P<0.01）］，結(jié)合圖10分析，3個(gè)時(shí)期中模型均顯示降水與河道徑流呈顯著正相關(guān)，且降水對(duì)徑流變化貢獻(xiàn)度達(dá)70%左右，降水是徑流量變化的主導(dǎo)因素，與李凱等［37］的研究結(jié)果較為一致。

圖12 長(zhǎng)江源區(qū)1980—2004年和2005—2020年降水、土壤侵蝕與水沙變化耦合關(guān)系Fig.12 Coupling relationship between precipitation,soil loss,runoff,and sediment load in the headwaters of Yangtze River during 1980-2004 and 2005-2020

通過(guò)對(duì)圖10 和圖12 中3 個(gè)時(shí)期耦合關(guān)系結(jié)果對(duì)比發(fā)現(xiàn)降水對(duì)土壤侵蝕的正向影響也略有上升，且均達(dá)到顯著水平。結(jié)合降雨侵蝕力（圖3b）與土壤侵蝕變化，認(rèn)為降水是引起源區(qū)侵蝕的主要因素。但是，源區(qū)仍面臨由于人類活動(dòng)導(dǎo)致草場(chǎng)退化及土地沙化的嚴(yán)峻形勢(shì)［38］，未來(lái)有可能會(huì)加劇源區(qū)侵蝕并影響源區(qū)徑流與泥沙變化。

目前大多數(shù)研究認(rèn)為區(qū)域土壤侵蝕對(duì)輸沙量有較大影響［39］，然而，本文研究結(jié)果顯示土壤侵蝕與輸沙量變化沒(méi)有顯著相關(guān)關(guān)系，源區(qū)土壤侵蝕與含沙量、輸沙量變化趨勢(shì)也明顯不同。一方面，水文特性和流域的形態(tài)等因素也會(huì)引起輸沙量變化；另一方面由于源區(qū)侵蝕泥沙在坡面淤積且源區(qū)內(nèi)河湖、濕地等存儲(chǔ)效應(yīng)的影響［18］，從而導(dǎo)致輸沙量變化較為穩(wěn)定，所以源區(qū)輸沙量的變化與侵蝕并沒(méi)有表現(xiàn)出必然聯(lián)系［40］。

3.3 土地利用及植被覆蓋變化對(duì)流域水沙的影響

土地利用/覆被變化對(duì)流域水文活動(dòng)起著重要的調(diào)節(jié)作用［9］，其對(duì)流域產(chǎn)流產(chǎn)沙的影響是流域水文學(xué)和生態(tài)水文學(xué)研究的關(guān)鍵問(wèn)題之一。

21世紀(jì)初期長(zhǎng)江源區(qū)人類擾動(dòng)較顯著，2010年后較為平穩(wěn)［41］。在本研究時(shí)段內(nèi)（1980—2020 年）源區(qū)草地多年平均占比約在68%～70%，其次是未利用地和水域，近40 a 土地利用格局相對(duì)穩(wěn)定（圖2）。基于1998—2009年、2010—2018年土地利用和1 km 空間分辨率的歸一化植被指數(shù)（NDVI）數(shù)據(jù)與水沙進(jìn)行相關(guān)性分析發(fā)現(xiàn)（表3），土地利用和植被指數(shù)變化與水沙變化整體幾乎呈現(xiàn)不顯著相關(guān)關(guān)系，與郭帥等［42］和Zhang 等［5］的研究結(jié)果相同。同時(shí)，Naveed 等［12］和Ji 等［43］認(rèn)為長(zhǎng)江源區(qū)徑流量變化主要由氣候變化引起，與土地利用相比，氣候變化的影響要大得多。另外，在本研究中土地利用的變化已囊括于土壤侵蝕的計(jì)算中，因此，本文重點(diǎn)探討了降水對(duì)水沙變化的影響，不再單獨(dú)討論土地利用和植被覆蓋變化對(duì)流域侵蝕與水沙變化的影響。隨著人類活動(dòng)的增加，加之氣候變化的影響和長(zhǎng)江源區(qū)各項(xiàng)生態(tài)項(xiàng)目的不斷推進(jìn)，未來(lái)長(zhǎng)江源區(qū)土地利用和植被覆蓋因子的變化不可忽視，應(yīng)注重未來(lái)源區(qū)生態(tài)保護(hù)及水源涵養(yǎng)［37］。

表3 1998—2009年、2010—2018年土地利用、植被覆蓋變化與徑流、泥沙相關(guān)性分析Tab.3 Correlation analysis of landuse,NDVI,runoff,and sediment load from 1998-2009,2010-2018

4 結(jié)論

本文基于中國(guó)土壤流失方程和結(jié)構(gòu)方程模型分析長(zhǎng)江源區(qū)水土流失與水沙變化特征及其與降水的響應(yīng)關(guān)系，得出以下結(jié)論：

（1）長(zhǎng)江源區(qū)1980—2020 年土壤侵蝕強(qiáng)度整體較低，以輕度侵蝕為主，約占總侵蝕面積的70.27%，年均土壤侵蝕模數(shù)為4.71 t·hm-2·a-1，呈顯著增加趨勢(shì)（P<0.05），但在2000年左右增加趨勢(shì)最為明顯，主要是由于降水的影響。侵蝕強(qiáng)度分布具有明顯空間異質(zhì)性，強(qiáng)度較大的地方主要分布在源區(qū)東南部。

（2）1980—2020年直門達(dá)水文站多年徑流量呈顯著增加趨勢(shì)，年均徑流量為142.03×108m3，輸沙量和含沙量均無(wú)明顯變化趨勢(shì)，水沙變化具有較好的一致性，徑流量在2004 年發(fā)生突變，突變后呈現(xiàn)增加趨勢(shì)。

（3）降水對(duì)土壤侵蝕、徑流和輸沙均具有顯著正向效應(yīng)，且降水是土壤侵蝕發(fā)生的主導(dǎo)因素，降水主要通過(guò)影響徑流對(duì)輸沙變化產(chǎn)生正向影響。未來(lái)應(yīng)及時(shí)關(guān)注降水帶來(lái)的土壤侵蝕影響，加強(qiáng)未來(lái)源區(qū)生態(tài)保護(hù)和恢復(fù)工作。