近60 a渭干河上游與黑孜河交匯處水沙變化特征及影響因素分析

張明旺，侍克斌，朱躍亮，孔選昭，姚海坤

(1.新疆農(nóng)業(yè)大學(xué)水利與土木工程學(xué)院，新疆 烏魯木齊 830052；2.新疆克孜爾水庫(kù)管理局，新疆 阿克蘇 842313)

在氣候變化和人類活動(dòng)的共同影響下，流域的整個(gè)水循環(huán)正在發(fā)生重大變化[1]。水沙是河流系統(tǒng)的最基本要素，水沙的年內(nèi)分配及年際變化規(guī)律是水利部門掌握河流特征的重要基礎(chǔ)。一方面，水沙的年內(nèi)分配決定了水利工程的運(yùn)行和調(diào)度方式，以及河流上游和下游的河槽演變過程;另一方面，水沙的年際變化在一定程度上揭示了長(zhǎng)期以來自然及人類活動(dòng)影響下的河流變化規(guī)律。因此，研究水沙變化特征具有十分重要的科學(xué)意義[2-3]。

近年來，對(duì)流域內(nèi)的水沙問題研究逐漸深入，如：姚文藝等[4]對(duì)黃河百年水沙變化進(jìn)行分析，發(fā)現(xiàn)近30 a泥沙主要減于中游，而徑流主要減于上游，水沙搭配關(guān)系趨好；楊永利[5]對(duì)大洋河沙里站近50 a的泥沙徑流資料進(jìn)行研究，發(fā)現(xiàn)輸沙量存在明顯的下降趨勢(shì)，徑流量的下降趨勢(shì)不顯著；張輝等[6]分析淮河干流水沙變化特征，淮河干流年徑流量沒有明顯增加或減小趨勢(shì)，但來沙量表現(xiàn)出明顯的減少趨勢(shì)；葉晨等[7]以五華河30 a的徑流、泥沙和降雨為研究對(duì)象，發(fā)現(xiàn)降雨變化對(duì)徑流、泥沙變化的影響為73%和21%，人類活動(dòng)對(duì)徑流、泥沙變化的影響為27%和79%；冉大川等[8]對(duì)黃河頭道拐多年平均徑流、輸沙量的增減幅度進(jìn)行研究，發(fā)現(xiàn)同時(shí)期輸沙量減少趨勢(shì)更明顯，灌區(qū)引水和支流水土保持綜合治理減沙作用明顯；胡春宏等[9]分析了中國(guó)主要河流代表站總徑流量和總輸沙量的變化趨勢(shì)，發(fā)現(xiàn)自然條件以及人類活動(dòng)是影響中國(guó)主要河流水沙變化的關(guān)鍵性因素；趙廣舉等[10]分析皇甫川流域水沙變化特征，結(jié)果表明人類活動(dòng)對(duì)流域水沙變化的貢獻(xiàn)占主導(dǎo)因素；趙秀蘭等[11]分析祖厲河流域的水沙及周期變化,為當(dāng)?shù)氐乃Y源管理及水土流失治理提供科學(xué)依據(jù)。

水沙變化特性及其影響因素的研究一直是國(guó)內(nèi)外眾多學(xué)者關(guān)注的重點(diǎn)。在前人的研究中，確定趨勢(shì)性及突變點(diǎn)的方法大多為單一方法，缺乏相互印證，相互對(duì)比。因此，本文采用線性趨勢(shì)估計(jì)法、Mann-Kendall秩次相關(guān)檢驗(yàn)法分析渭干河上游與黑孜河交匯處水沙序列的趨勢(shì)性；采用有序類聚法、Mann-Kendall非參數(shù)檢驗(yàn)方法和滑動(dòng)t檢驗(yàn)法分析其水沙序列的突變性；采用雙累積曲線法定量分析降水及人類活動(dòng)對(duì)其徑流量和輸沙量的影響程度。分析渭干河上游與黑孜河交匯處多年來的水沙變化規(guī)律及其影響因素，可為水庫(kù)的優(yōu)化運(yùn)行調(diào)度、水沙預(yù)測(cè)及清淤工作提供一定的基礎(chǔ)數(shù)據(jù)，具有重要的理論意義和參考價(jià)值。

1 材料與方法

1.1 研究區(qū)概況

克孜爾水庫(kù)是渭干河灌區(qū)上游最大的水庫(kù)控制工程，它承擔(dān)著下游農(nóng)業(yè)灌溉和生活用水安全保障的重要責(zé)任。入庫(kù)水沙主要來自于哈爾克他烏山脈的干流木扎提河及支流克孜河，干流木札提河由卡木斯浪河、臺(tái)勒外丘克河及卡拉蘇河匯合而成。本文研究數(shù)據(jù)來自于干流木札提河及支流卡木斯浪河、臺(tái)勒外丘克河、卡拉蘇河、克孜河上的各水文站點(diǎn)，研究區(qū)域水系及站點(diǎn)分布見圖1。文中選擇克孜爾水庫(kù)1960年至2017年入庫(kù)逐年實(shí)測(cè)徑流量與輸沙量資料作為基礎(chǔ)數(shù)據(jù)進(jìn)行水沙變化分析。

圖1 渭干河流域水系及水文站網(wǎng)分布圖

1.2 研究方法

1.2.1 趨勢(shì)性分析 研究首先使用線性趨勢(shì)估計(jì)法[12]分析克孜爾水庫(kù)上游徑流和輸沙序列的長(zhǎng)期變化趨勢(shì)，然后應(yīng)用Mann-Kendall秩次相關(guān)檢驗(yàn)法[13]識(shí)別水文時(shí)間序列趨勢(shì)[14]。

1.2.2 突變性分析 首先采用Mann-Kendall參數(shù)檢驗(yàn)法分析水庫(kù)上游徑流和輸沙序列的突變性，接著采用有序聚類法和滑動(dòng)t檢驗(yàn)法[15]檢驗(yàn)水沙序列的突變性。非單一的確定方法避免了結(jié)果的偶然性，提高了結(jié)果的精度，相互印證對(duì)比，綜合3種判定方法最終確定突變點(diǎn)。

1.2.3 周期性分析 采用基于小波分析的EMD周期分析方法，EMD方法可以平滑復(fù)雜的序列，從原始序列中提取出一系列具有不同層次的波動(dòng)，并獲得具有不同尺度的多個(gè)IMF分量。這是一種有效的篩選方法，可以逐步分解序列中不同尺度的波動(dòng)和趨勢(shì)分量[16]。

2 結(jié)果與分析

2.1 水庫(kù)水沙特性分析

2.1.1 趨勢(shì)性分析 分析圖2-A可知，多年徑流量總體上呈上升趨勢(shì)，年徑流量上升的趨向率為1.165×108m3/10a，用Mann-kendall秩次相關(guān)檢驗(yàn)法，確定徑流量變化趨勢(shì)的顯著程度，徑流量上升趨勢(shì)顯著(表1)。該流域的徑流來源主要為冰川積雪融水，隨著全球氣候變暖現(xiàn)象逐漸加劇，冰川每年的消融量也逐漸變大，導(dǎo)致河流的徑流量呈上升趨勢(shì)，與分析結(jié)果相統(tǒng)一。

分析圖2-B可知，多年輸沙量總體呈上升趨勢(shì)，年輸沙量上升的趨向率為104.3×104t/10 a，用Mann-kendall秩次相關(guān)檢驗(yàn)法，確定輸沙量變化趨勢(shì)的顯著程度，輸沙量上升趨勢(shì)不顯著(表1)。河流的含沙量與當(dāng)?shù)氐耐恋乩妙愋陀嘘P(guān)，由于該地區(qū)的土地利用類型逐漸增多，土壤撓動(dòng)問題變得尤為突出，導(dǎo)致進(jìn)入河道的泥沙量增多。同時(shí)，洪水的發(fā)生頻率在近些年也增加，洪水沖刷地表，攜帶大量的疏松土壤進(jìn)入河流，使河流的輸沙量也呈上升趨勢(shì)。

表1 克孜爾水庫(kù)上游流域年徑流量與年輸沙量趨勢(shì)性分析結(jié)果

圖2 徑流量與輸沙量年際變化過程線

2.1.2 突變性分析 使用有序類聚法分析年徑流量序列的突變，推估年徑流量序列的分割點(diǎn)，并計(jì)算年徑流序列的總離差平方和Sn(τ)與分割點(diǎn)τ之間的關(guān)系曲線(圖3-B)。初步分析年徑流量的變異點(diǎn)，最小的年份是1995年，第二小的年份是1993年；為進(jìn)一步分析變異點(diǎn)的顯著性，使用非參數(shù)檢驗(yàn)方法Mann-kendall對(duì)變異點(diǎn)的顯著性進(jìn)行檢驗(yàn)(3-A)，突變的時(shí)間是1994年；采用滑動(dòng)t檢驗(yàn)法對(duì)突變點(diǎn)進(jìn)行檢驗(yàn)(圖3-C)，突變的時(shí)間是1994年；結(jié)合以上的分析數(shù)據(jù)，確定克孜爾水庫(kù)年徑流量序列的突變年份為1994年。

圖3 徑流量突變性分析

利用有序類聚法對(duì)年輸沙量序列進(jìn)行突變分析，推估年輸沙量序列的分割點(diǎn)，計(jì)算得到年輸沙序列總離差平方和Sn(τ)與分割點(diǎn)τ的關(guān)系曲線(圖4-B)，初步分析年徑流量的變異點(diǎn)，最小值為1977年，第二小值為1976年；為進(jìn)一步分析變異點(diǎn)的顯著性，使用非參數(shù)檢驗(yàn)方法Mann-Kendall對(duì)變異點(diǎn)的顯著性進(jìn)行檢驗(yàn)(圖4-A)，曲線的交點(diǎn)在1966年、2008年與2009年之間；采用滑動(dòng)t檢驗(yàn)法對(duì)突變點(diǎn)進(jìn)行檢驗(yàn)(圖4-C)，突變的時(shí)間是1978年、1991年和1992年；輸沙序列的突變性檢驗(yàn)經(jīng)多種方法計(jì)算未能達(dá)到相統(tǒng)一的結(jié)論，這是因?yàn)榭俗螤査畮?kù)的建成時(shí)間為1992年，當(dāng)水庫(kù)建成以后，泥沙由于建筑物的阻攔，改變了其原來的運(yùn)移狀態(tài)，輸沙量的突變年份與克孜爾水庫(kù)的修建存在一定的關(guān)聯(lián)性，在多種方法判定的基礎(chǔ)上，最終確定克孜爾水庫(kù)年輸沙量序列的突變年份為1992年。

圖4 輸沙量突變性分析

2.1.3 周期性分析 1960～2017年克孜爾水庫(kù)上游流域徑流量序列為非平穩(wěn)序列，可以分解為3個(gè)具有不同波動(dòng)周期的IMF分量和1個(gè)Res趨勢(shì)分量(圖5)，反映了克孜爾水庫(kù)上游流域徑流量變化復(fù)雜的時(shí)間尺度特征。

圖5 克孜爾水庫(kù)上游流域年徑流量序列EMD分解

第一個(gè)本征模態(tài)函數(shù)IMF1具有最大振幅，最高頻率和最短波長(zhǎng)，IMF2和IMF3分量的幅度逐漸減小，頻率逐漸減小，波長(zhǎng)逐漸增大。

Res分量顯示了年徑流量的總體變化趨勢(shì)，并且可以看出多年來有明顯的上升趨勢(shì)。

為進(jìn)一步獲得各IMF分量準(zhǔn)確的周期特點(diǎn)，在EMD分解得到IMF分量的基礎(chǔ)上，對(duì)各分量進(jìn)行Morlet小波分析。利用matlab軟件編程，得到各IMF分量的小波方差圖，見圖6。

由圖6可知，IMF1分量具有兩個(gè)明顯的峰值，時(shí)間尺度為6～9 a，準(zhǔn)周期為7 a；IMF2分量具有較高的峰值，時(shí)間尺度為15～18 a，準(zhǔn)周期為17 a；而IMF3分量有一個(gè)峰值，時(shí)間尺度為26～32 a，準(zhǔn)周期均為30 a。分析顯示了克孜爾水庫(kù)上游年徑流量變化呈復(fù)雜的多時(shí)間尺度性，本結(jié)論可為克孜爾水庫(kù)未來的水資源配置、防洪減災(zāi)和抗旱提供科學(xué)指導(dǎo)。

圖6 年徑流量EMD分量的小波方差

克孜爾水庫(kù)上游流域1960～2017年間的輸沙量序列為非平穩(wěn)序列，可以分解為4個(gè)具有不同波動(dòng)周期的IMF分量和1個(gè)Res趨勢(shì)分量(圖7)，反映了克孜爾水庫(kù)上游流域輸沙量變化的復(fù)雜的多時(shí)間尺度性；

圖7 克孜爾水庫(kù)上游流域年輸沙量序列EMD分解

第一個(gè)本征模函數(shù)IMF1具有最大振幅，最高頻率和最短波長(zhǎng)。IMF2，IMF3和IMF4的分量幅度逐漸減小，頻率逐漸減小，波長(zhǎng)逐漸增大；

Res分量顯示了年輸沙量的總體變化趨勢(shì)，可以看出多年來有明顯的上升趨勢(shì)。

為進(jìn)一步獲得各IMF分量準(zhǔn)確的周期特點(diǎn)，在EMD分解得到IMF分量的基礎(chǔ)上，對(duì)各分量進(jìn)行Morlet小波分析。利用matlab軟件編程，得到各IMF分量的小波方差圖(圖8)。

由圖8可知，IMF1分量具有兩個(gè)明顯的峰值，時(shí)間尺度為5～7 a，其準(zhǔn)周期為6 a；IMF2分量存在一個(gè)明顯的峰值，時(shí)間尺度為10～13 a，其準(zhǔn)周期為11 a；IMF3分量存在一個(gè)明顯的峰值，時(shí)間尺度為25～28 a，其準(zhǔn)周期為26 a；IMF4分量無峰值。分析顯示了克孜爾水庫(kù)上游年輸沙量變化呈復(fù)雜的多時(shí)間尺度性，本結(jié)論可為克孜爾水庫(kù)的優(yōu)化調(diào)度運(yùn)行及排沙清淤延長(zhǎng)水庫(kù)壽命提供幫助。

圖8 年輸沙量EMD分量的小波方差

2.2 影響因素分析

2.2.1 氣候活動(dòng)和人類活動(dòng)對(duì)徑流量及輸沙量影響程度 為了區(qū)分人類活動(dòng)和降水對(duì)克孜爾水庫(kù)上游水沙變化的影響，首先對(duì)克孜爾水庫(kù)上游1960～2017年的年徑流量和輸沙量進(jìn)行突變性檢驗(yàn)，將突變年份之前的階段視為僅有降雨影響的基準(zhǔn)期，突變年份之后的時(shí)期視為有人類活動(dòng)和降雨共同作用的時(shí)期(圖9)。由前部分計(jì)算結(jié)果所示，1960～2017年克孜爾水庫(kù)上游徑流量的突變年份為1994年，輸沙量的突變年份為1992年，其次，用降水-徑流量和降水-輸沙量的雙累積曲線對(duì)人類活動(dòng)和降雨量對(duì)于徑流量和輸沙量的貢獻(xiàn)進(jìn)行模擬回歸分析(表2)。通過降水-徑流量和降水-輸沙量雙累積曲線擬合方程得出一個(gè)計(jì)算值，計(jì)算值與實(shí)測(cè)值的差值即為人類活動(dòng)貢獻(xiàn)(表3)。

表2 克孜爾水庫(kù)上游降水-徑流量和降水-輸沙量雙累積曲線線性模擬評(píng)估

表3 克孜爾水庫(kù)上游突變年份前后時(shí)期降水和人類活動(dòng)的貢獻(xiàn)率

圖9 累積降水量與累積徑流量、累積輸沙量雙累積曲線

從1960年到2012年，與計(jì)算出的累計(jì)徑流量相比，實(shí)測(cè)累積徑流量減少了72.73億m3。徑流量在突變年份之后，實(shí)測(cè)累積徑流量增加了0.82億m3。人類活動(dòng)對(duì)上游徑流量增加的貢獻(xiàn)率占-385.37%，而降水量對(duì)徑流量增加的貢獻(xiàn)率占485.37%。結(jié)果表明：在克孜爾水庫(kù)上游徑流量突變年份之后，降水對(duì)徑流量的影響比人類活動(dòng)大?？俗螤査畮?kù)上游流域內(nèi)的人口密度不是很大，靠近河流的盡頭是冰川的發(fā)源地，人類活動(dòng)很少，此處的徑流量受降水的影響很明顯，降水是此區(qū)域的主要影響因素。

從1960年到2012年，與計(jì)算出的累計(jì)輸沙量相比，實(shí)測(cè)累積輸沙量增加了819.47萬(wàn)t。輸沙量在突變年份之后，實(shí)測(cè)累積輸沙增加了585.9萬(wàn)t。人類活動(dòng)對(duì)上游輸沙量增加的貢獻(xiàn)率占26.43%，而降雨量對(duì)徑流量的增加貢獻(xiàn)率占73.57%。結(jié)果表明：在克孜爾水庫(kù)上游輸沙量在突變年份之后，降水是影響輸沙量的主要因素。此區(qū)域地表植被的覆蓋率不是很高，氣溫較高，蒸發(fā)量很大，一定程度上導(dǎo)致地表土壤疏松，在降水時(shí)候會(huì)沖刷地表，帶走疏松的土壤進(jìn)入河道；此區(qū)域的水土保持措施不是很完善，隨著時(shí)間的推移，原有的水土保持措施已漸漸失去了當(dāng)初的效果，不能滿足當(dāng)前的需求，需結(jié)合現(xiàn)狀進(jìn)一步推進(jìn)水土保持工作的開展。

2.2.2 太陽(yáng)活動(dòng)對(duì)徑流量與輸沙量的影響 采用前文的基于小波分析法的EMD分解周期法，求解出太陽(yáng)活動(dòng)的小波方差分析圖(圖10)。由圖10可知，IMF1分量具有兩個(gè)明顯的峰值，其準(zhǔn)周期為16 a和31 a；IMF2分量存在一個(gè)明顯的峰值，其準(zhǔn)周期為31 a；IMF2分量無明顯峰值；IMF4分量存在一個(gè)明顯的峰值，其準(zhǔn)周期為23 a。

圖10 太陽(yáng)黑子活動(dòng)EMD分量的小波方差圖

克孜爾水庫(kù)上游年徑流量存在7 a、17 a、30 a左右的主周期，第一主周期是30 a。在15～18 a和26～32 a尺度上，年徑流量和太陽(yáng)黑子相對(duì)數(shù)序列具有相似的主周期，但兩個(gè)周期在較短的時(shí)間尺度上是不同的。以年徑流量變化為代表的7 a周期與黑子活動(dòng)無關(guān)，而17 a和30 a周期應(yīng)主要受黑子活動(dòng)影響，并且與黑子活動(dòng)密切相關(guān)。徑流量的增大會(huì)影響河流中攜帶的泥沙，進(jìn)一步影響克孜爾水庫(kù)的入庫(kù)輸沙量，輸沙量的周期波動(dòng)性與徑流量緊密相關(guān)，徑流量的變化受到太陽(yáng)活動(dòng)的影響，輸沙量同時(shí)也因徑流量的變化而受到影響。

2.2.3 ENSO對(duì)徑流量及輸沙量的影響分析 厄爾尼諾-南方濤動(dòng)(El Nio-Southern Oscillation，ENSO)是熱帶太平洋最顯著的海氣耦合信號(hào)，暖、冷位相分別表現(xiàn)為厄爾尼諾事件和拉尼娜事件[17]。ENSO事件在南疆地區(qū)表現(xiàn)的周期性為：20世紀(jì)1980年之前和20世紀(jì)90年代末期表現(xiàn)出6 a的周期性，20世紀(jì)80年代到90年代中期表現(xiàn)出準(zhǔn)4 a的周期變化，1995年后這種周期有變長(zhǎng)的趨勢(shì)。這與克孜爾水庫(kù)水沙變化的6～9 a的短周期一致，可以推測(cè)出克孜爾水庫(kù)水沙受ENSO信號(hào)的調(diào)節(jié)。分析二者之間的關(guān)系，可以發(fā)現(xiàn)水沙變化周期相比較ENSO事件具有3～4 a的滯后，這是因?yàn)镋NSO對(duì)新疆雨季降水的影響具有明顯的滯后性，新疆降水和ENSO之間存在著一種物理的遙相關(guān)關(guān)系[18]。

3 討論

水沙變化趨勢(shì)性分析方法主要包括雙累積曲線法、滑動(dòng)平均法以及Spearman秩相關(guān)分析法等[19]；水沙變化突變性分析主要采用Pettitt檢驗(yàn)法、Mann-Kendall非參數(shù)秩檢驗(yàn)法、水文情勢(shì)突變指標(biāo)法及有序聚類法等[20-21]；水沙周期性分析主要采用灰譜分析、EMD分解、極大熵譜分析以及小波分析法等[22-23]。本文對(duì)趨勢(shì)性判斷采用的是線性趨勢(shì)法及Mann-Kendall秩次檢驗(yàn)法，得到的結(jié)果與實(shí)際情況相符合，對(duì)于水沙的未來趨勢(shì)預(yù)測(cè)還有待深入研究，依據(jù)其過去的數(shù)據(jù)分析未來趨勢(shì)；對(duì)突變點(diǎn)采用有序類聚法、Mann-kendall非參數(shù)檢驗(yàn)方法和滑動(dòng)t檢驗(yàn)法分析其突變的年份，在文中對(duì)徑流量的突變檢驗(yàn)結(jié)果可以得到統(tǒng)一結(jié)論，對(duì)輸沙量的突變檢驗(yàn)結(jié)果偏差很大，需要結(jié)合水庫(kù)的實(shí)際情況進(jìn)行判斷。本文采取的3種突變檢驗(yàn)方法是否適合輸沙序列的突變檢驗(yàn)，還需要深入研究。水沙序列周期的分析采用基于小波分析的EMD分解方法，其重要的環(huán)節(jié)為數(shù)據(jù)系列的分解，經(jīng)典模態(tài)分解在本文已做闡述，對(duì)水沙序列分解的方法不同，其精度是否會(huì)提高，仍需繼續(xù)研究。本文主要分析了降水、人類活動(dòng)、太陽(yáng)活動(dòng)及ENSO現(xiàn)象對(duì)水沙變化的影響，至于其他影響水沙變化的因素還需進(jìn)一步研究。對(duì)于克孜爾水庫(kù)上游流域而言，大氣環(huán)境變化導(dǎo)致的冰川融水也會(huì)對(duì)其流域的水沙變化造成一定的影響，靠冰川積雪融水補(bǔ)給河流，是該流域的一大特色，需要深入研究。

4 結(jié)論

1) 克孜爾水庫(kù)壩址處入庫(kù)徑流量及輸沙量整體上均呈增加趨勢(shì)。壩址處入庫(kù)徑流量及輸沙量發(fā)生由少至多的突變年份分別為1994年、1992年。入庫(kù)徑流量存在3種時(shí)間長(zhǎng)度，表現(xiàn)為6～9 a、15～18 a、26～32 a，周期具體表現(xiàn)為7 a、17 a、30 a；輸沙量存在3種時(shí)間長(zhǎng)度，表現(xiàn)為5～7 a、10～13 a、25～28 a，周期具體表現(xiàn)為6 a、11 a、26 a；徑流量及輸沙量的第一主周期分別為30 a、26 a。

2) 徑流量在突變年份之后，實(shí)測(cè)累積徑流量增加0.82×108m3，人類活動(dòng)對(duì)徑流量增加貢獻(xiàn)率占-385.37%，而降水對(duì)徑流量的增加貢獻(xiàn)率占485.37%；輸沙量在突變年份之后，實(shí)測(cè)累積輸沙增加585.9×104t。人類活動(dòng)對(duì)輸沙量增加的貢獻(xiàn)率占26.43%，而降水對(duì)輸沙量的增加貢獻(xiàn)率占73.57%。

3) 年徑流量變化表現(xiàn)的7 a周期與太陽(yáng)黑子活動(dòng)關(guān)系不大，而17 a和30 a周期應(yīng)主要是受太陽(yáng)黑子活動(dòng)的影響，與太陽(yáng)黑子活動(dòng)關(guān)系緊密?？俗螤査畮?kù)水沙受ENSO信號(hào)的調(diào)節(jié)，水沙變化周期相比較ENSO事件具有3～4 a的滯后，這是因?yàn)镋NSO對(duì)新疆雨季降水的影響具有明顯的滯后性。