牧馬河流域年徑流變化特征分析

夏煥清，霍秀秀，曹苗苗

（陜西省漢中水文水資源勘測局，陜西 漢中 732000）

1 研究區(qū)概況與數(shù)據(jù)來源

1.1 研究區(qū)概況

牧馬河為漢江干流一級支流，源于城固、南鄭、西鄉(xiāng)3縣交界處米倉山北麓的白熊山老鷹崖，流域位于大巴山北部，山地高在2500 m～3000 m之間，大致東向流經(jīng)城固縣、西鄉(xiāng)縣，在西鄉(xiāng)城南三花石鄉(xiāng)回龍灣匯入漢江。西鄉(xiāng)水文站是牧馬河重要控制站，集水面積1224 km2，多年平均降雨量為877.8 mm，7月～9月降水量占全年的50%以上，多年平均徑流量為9.04億m3。

1.2 數(shù)據(jù)來源及研究方法

采用的徑流量數(shù)據(jù)為西鄉(xiāng)水文站(1974年～2018年)的年徑流量數(shù)據(jù)，見表1。

表1 西鄉(xiāng)水文站(1974年～2018年)的年徑流量數(shù)據(jù) 單位：億m3

研究方法為：采用線性回歸、五年滑動平均對牧馬河流域徑流變化趨勢進行分析，M-K檢驗、小波分析分別進行年徑流量突變檢驗與周期分析。數(shù)據(jù)處理工具：小波分析和M-K檢驗用Matlab，顯著性檢驗用皮爾遜相關(guān)性分析。

1.3 研究方法

1.3.1 線性回歸

通過建立牧馬河流域西鄉(xiāng)站年徑流與時間序列之間的線性回歸方程，通過回歸方程的斜率以及變化趨勢的顯著程度從而判斷因變量趨勢變化的方法，表達式：

式中：Xt為t時刻年徑流值；t為時間；b為方程截距；a為回歸方程斜率。

用10 a表示研究對象時間序列趨勢變化程度，當a＞0時，表明年徑流具有上升趨勢；當a＜0時，具有下降趨勢。此外，通過皮爾遜相關(guān)分析法來判斷年徑流變化趨勢是否顯著，當p＜0.05時，年徑流在時間序列上具有顯著性上升或下降趨勢；當p＞0.05時，年徑流不具有顯著變化趨勢，本文通過SPSS17.0對年徑流變化趨勢進行顯著性檢驗。

1.3.2 滑動平均

滑動平均是可以直觀反映年徑流時間序列變化趨勢的一種方法，通過計算西鄉(xiāng)站年徑流系列（1974年～2018年）44年徑流項連續(xù)5項序列的平均值。從第1項至5項的算術(shù)平均為滑動平均第1項，最后1項是原來44年徑流序列的第40至第44項的算術(shù)平均。若西鄉(xiāng)站年徑流序列具有趨勢成分，則5年滑動平均能夠?qū)⑽鬣l(xiāng)站年徑流的趨勢直觀的顯示出來。

1.3.3 Mann-Kendall突變檢驗

對于西鄉(xiāng)站年徑流樣本時間序列 X1、X2、X3、X4…Xn，mi表示第i個樣本Xi大于Xj（1≤j≤i）的累計數(shù)，則構(gòu)造統(tǒng)計量：

在序列相互獨立的假設條件下，的均值和方差為：

將公式標準化得到以下計算式：

UFk計算結(jié)果組成UF曲線，同樣以此方法引用到反序列中，計算得到 UBk，使得 UBk=-UFk。

給定顯著性水平0.05，統(tǒng)計量臨界值±1.96。通過繪制UFk、曲線圖UBk，如果兩條曲線在置信度0.05水平內(nèi)相交，則定義為突變點。如果UFk或UBk的值大于0，則表明序列具有上升趨勢，反之，序列呈下降趨勢。當曲線超出置信度區(qū)間時，表明序列具有顯著的上升或下降趨勢，該過程通過Matlab工具調(diào)用編輯的西鄉(xiāng)站年徑流M-K突變檢驗程序，來實現(xiàn)相關(guān)過程計算及作圖。

1.3.4 小波分析

小波分析方法是一種在空間固定的基礎(chǔ)上反映時頻局部發(fā)生變化的分析方法。基本原理用一簇小波函數(shù)系來表示或逼近某一信號或函數(shù)。小波分析最關(guān)鍵的部分是小波函數(shù)，小波函數(shù)φ（t）∈L2（R），且有：

φ（t）為基小波函數(shù)，它可通過尺度的伸縮和時間軸上的平移構(gòu)成一簇函數(shù)系：

式中：φa,b（t）為子小波；a為尺度因子，反映小波的周期長度；b為平移因子，反應時間上的平移。

對西鄉(xiāng)站年徑流的研究選用的是Morlet小波函數(shù)，其表達式為：

式中：ω0為常數(shù)（ω0≥5）；i為虛數(shù)。

對于任意函數(shù)f（X），其小波變換可以定義為：

牧馬河流域西鄉(xiāng)站年徑流周期分析計算與出圖過程，通過Matlab工具調(diào)用編輯的西鄉(xiāng)站年徑流周期分析程序來實現(xiàn)小波分析。

3 結(jié)果與分析

3.1 趨勢分析

圖1中圓點折現(xiàn)圖為年徑流量圖，三角折線圖為五年滑動平均圖，方程式為年徑流量圖趨勢線方程式。由圖可知，年徑流趨勢線k為-0.1505，即年徑流量趨勢為減小；利用皮爾遜相關(guān)性分析，在顯著性水平為0.01以下，皮爾遜相關(guān)性系數(shù)為-0.459，兩個變量為負相關(guān)且顯著相關(guān)；五年平均滑動圖趨勢線k為-0.1346，五年平均徑流量趨勢為減小。

圖1 年徑流及五年滑動平均圖

3.2 突變分析

圖2 UF小于0表明呈下降趨勢，當UF和UB兩條曲線出現(xiàn)交點，且交點在臨界直線之間，交點對應的時間即為突變開始時刻。從圖2可知，西鄉(xiāng)站年徑流量突變時間為1990年。

圖2 M-K突變檢驗圖

3.3 周期風險

在不同的時間尺度下，年徑流量的變化有周期性和時間域上豐、枯變化的特性。

圖3 牧馬河流域年徑流小波系數(shù)模方圖

圖3 為小波模方圖，由圖3來分析年徑流在小波變化域中波動能量大小的變化特性，進而反映哪些能量聚集中心主導年徑流在時間域上的波動變化的情形，分析可知，在5 a、10 a、23 a左右的時間尺度上能量聚集中心閉合較好并且分布密集，與圖5比較沒有30 a以上的，原因為時間序列尺度不夠長。

圖4 牧馬河流域年徑流小波變換實部圖

圖4 為小波系數(shù)實部的時頻分布圖，表示年徑流在小波變化域中以基準面上下起伏的波動情況，從而反映年徑流在時間尺度上豐、枯變化的特性，依圖可知西鄉(xiāng)站年徑流豐枯變換情況，豐水年在1980年、1989年、2011年左右，枯水年在1986年、2001年左右。

圖5 牧馬河流域年徑流小波變化方差圖

由圖5可知，小波方差在8 a、22 a、38 a左右的時間尺度上存在峰值，22 a左右的周期震蕩最強，是流域的第一主周期，第二、三周期依次約為8 a、38 a，由于西鄉(xiāng)水文站建站時間較短，年徑流數(shù)據(jù)序列較短，因此無法準確預測更長時間尺度的周期。

總體來看，西鄉(xiāng)水文站年徑流周期為22 a～23 a，三個豐水年之間的間隔為 9 a、22 a、31 a 與圖 5 分析的周期 8 a、22 a、38 a左右，契合度較高。但在長周期上有較大差別，分析主要原因是隨著時間序列的增加，小時間尺度上的突變點逐漸退化為非突變點，隨著時間尺度的降低，突變點增多，大時間尺度上的一些平常態(tài)，在低層次上轉(zhuǎn)化為突變態(tài)，而本文因時間序列較短，在長周期上判斷不夠準確，進一步說明年徑流周期隨時間尺度的不同而變化。

4 結(jié)論

采用采用線性回歸、滑動平均法、MK檢驗、皮爾遜相關(guān)分析和小波分析方法分析牧馬河西鄉(xiāng)水文站天然年徑流時間序列進行多時間尺度變化特性的分析。結(jié)果顯示:西鄉(xiāng)水文站的年徑流呈減小趨勢，年徑流突變時間為1990年左右，呈周期性變化，主周期為22 a～23 a左右，它主導著牧馬河年徑流長期變化的特性，同時也存在著8 a～10 a、38 a左右的較強周期變化。在1980年、1989年、2011年左右為豐水年，在1986年、2001年左右為枯水年。