1960－2017年珠江流域下游徑流年際與年內變化特征

楊遠東，王永紅，蔡斯龍，劉 鋒

（1.海底科學與探測技術教育部重點實驗室，中國海洋大學 海洋地球科學學院，山東 青島266100；2.廣東省水文局，廣州 廣東510150；3.中山大學 海洋科學學院，廣州 廣東510275）

河川徑流年際變化與年內變化特征對于河道整治、農業灌溉、工農業用水、水庫調水等方面意義重大。分析徑流年際、年內特征的方法較多，例如白麗等［1］使用徑流變差系數（Cv）法對1957—2015年贛江流域徑流特征及變化規律分析，戴明龍等［2］對長江流域1952—2006徑流時空分布及變化規律進行了研究。楊敏等［3］使用M-K 趨勢檢驗法對1951—2015年洞庭湖水沙演變進行了評估，李瑞等［4］使用此方法對渭干河流域1960—2006年徑流量進行了研究，除此之外，常用的方法還有R／S 分析法［5-7］，豐枯等級劃分法［8-9］，滑動平均法［10-11］等。徑流年內分析主要方法有集中度（期）法［12-13］，年內分配不均勻系數法［14］等。

這些方法在研究珠江流域徑流變化時也經常使用，例如吳創收等［15］對近100年來珠江流域西江年徑流量進行了研究，發現其存在著2～3 a，12 a左右以及27 a左右的振蕩周期；陳立華等［16］對珠江流域西江1980—2010年徑流進行了分析，發現西江下游水文站徑流存在不顯著的下降趨勢；李艷等［17］對珠江流域北江1956—2000年徑流進行了分析，發現年徑流量系列呈不顯著緩慢上升趨勢，并存在2年左右的變化周期；劉金鳳等［18］對珠江流域東江1956—2015年徑流進行研究，發現近60 a東江流域年徑流變化趨勢總體平穩，未有顯著變化。這些研究往往只集中在某一流域，研究的內容也僅僅對徑流年內變化或年際變化單獨研究，未能將其綜合對比。本文采用徑流年際變差系數（Cv）、年際極值比（W）、M-K 趨勢分析、滑動平均法、Morlet小波分析法、集中度（期）等方法對珠江流域西江（高要）、北江（石角）、東江（博羅）水文站1960—2017年徑流進行分析，以期更好理解近幾十年內其年際與年內變化特征。

1 材料與方法

1.1 研究區概況

珠江發源于云南省馬雄山，由西江、北江與東江3大支流構成，全長2 400 km，流域面積45×104km2，流經云南、貴州、廣西、廣東、湖南、江西6 個省（區）和越南的北部，是中國第二大河，境內第三長河。珠江流域地處亞熱帶，北回歸線橫貫流域的中部，氣候溫和多雨，多年平均溫度在14 ℃～22 ℃之間，多年平均降雨量1 200～2 200 mm，多年平均徑流量3.10×1011m3多（1954—2016年），多年平均輸沙率2.80×107t多（1954—2016年），徑流水沙季節分配很不均勻，常年4—9月份洪水期多年平均徑流量約占全年的80%，輸沙率占91～95%。珠江流域西北兩江分別在廣東省三水市思賢滘東西兩側、東江在廣東省東莞市石龍鎮匯入珠江三角洲。高要（112°27′13″E，23°2′43″N）、石角（112°56′59″E，23°33′39″N）、博羅（114°17′42″E，23°9′50″N）分別為進入三角洲之前西、北、東江最主要的水文控制站。

1.2 數據來源

本文收集了珠江流域進入珠江三角洲之前西、北、東江的主要水文控制站高要、石角、博羅1960—2017年月度徑流量數據，長時間序列徑流數據來源于廣東省水文局，中山大學與《中國河流泥沙公報（http：∥www.mwr.gov.cn／）。珠江流域主要大壩建設數據來源于中國大壩工程學會網站（http：∥www.chincold.org.cn／）。高要、石角、博羅水文站歷史位置未發生改變，同時對3站徑流年際累計值進行分析，分析結果如圖1所示。3站徑流年際累計值基本為一條直線，其相關系數R 值接近1，徑流分析結果表明3站徑流數據一致性較強。

圖1 1960－2017年高要、石角、博羅水文站際累計徑流量

2 結果與分析

2.1 徑流年際變化特征

2.1.1 徑流年際變異性 河川徑流量的年際變化主要取決于降水量的年際變化，同時還受流域水源補給類型及流域內地貌、地質等條件的影響［19］，徑流年際變異性常用變差系數（Cv）與年際極值比（w）來表示，變差系數（Cv）是水文統計中的一個重要參數，用來說明水文變量長期變化的穩定程度［20］，其值為樣本序列的標準差與算術平均值的比值，Cv值越大，年徑流量的年際變化越大，是反映徑流年際變化不均勻性的重要指標［21-22］。變差系數（Cv）計算公式如（1）［23］：

式中：ki——模比系數，即第i年的年徑流量與平均年徑流量的比值；n——系列長度。

年際極值比（w）是指徑流量變化的絕對值比例，其大小用多年最大、最小年徑流量比值表示［24］，年際極值比計算公式如（2）［23］：

珠江流域西、北、東江1960—2017年徑流年際變化特征值計算結果見表1。由表1可知，西、北、東江徑流變差系數（Cv）值均較小，其值介于0.196～0.273之間，說明珠江流域徑流年際變率較小，徑流量穩定性較高，穩定性的大小順序為西江大于北江大于東江。西、北、東江徑流變差系數（Cv）與徑流年際極值比（w）變化趨勢相似，徑流年際極值比（w）值介于2.398～4.636之間，西江w 值較小，北江、東江w 值較為接近且較大。徑流變差系數（Cv）值與徑流年際極值比（w）均反映1960—2017年珠江流域徑流量年際變異程度東江高于北江高于西江。

表1 1960－2017年珠江流域西、北、東江徑流年際變化特征

2.1.2 徑流年際趨勢性 用于診斷連續水文序列趨勢性方法較多，其中應用最廣范、最為簡便直觀的方法為M-K 趨勢分析法與滑動平均法，本文分別用這兩種方法對珠江流域西、北、東江1960—2017年徑流趨勢性進行分析。M-K 趨勢分析法（Mann-Kendall法）是一種非參數統計檢驗方法，其優點是樣本不需要遵循一定的分布，也不受少數異常值的干擾，更適用于類型變量與順序變量，計算也比較簡便，最初由曼和肯德爾提出原理并發展了這一方法。現已成為世界氣象組織推薦并已廣泛使用的非參數檢驗方法，適用于水文、氣象等非正態分布的數據。其計算公式如公式（3）—（6），設徑流序列為x1，x2，…，xn，Sk表示第i個樣本xi＞xj（1≤j≤i）的累計數，

其中：

在時間序列隨機獨立的假設下，Sk的均值與方差分別為：

將Sk標準化為：

式中：UF1＝0，給定顯著性水平α，若｜UF｜＞Uα，表明序列存在明顯的趨勢變化。將此方法引用到反序列，序列順序變為xn，xn－1，…，x1，反序列標準化由UBk表示，UBk和UFk可組成2條曲線，在這里設定α 值為0.05，臨界值U0.05為±1.96。

珠江流域西、北、東江1960—2017 年徑流M-K趨勢分析結果如圖2所示。西江高要水文站1960—2017年UF與UB 統計曲線在信度區間內發生了5次交叉，根據統計學檢驗意義，1964，1986，1991，1999，2016年是高要水文站徑流可能的突變年份，其中，UF值并未突破0.05置信區間，說明徑流有上升、下降趨勢但趨勢不明顯；北江石角站1960—2017 年UF與UB 統計曲線在信度區間內發生了6次交叉，可能 的 突 變 年 份 分 別 為1961，1973，1984，1990，1999，2016年，但UF 值未突破0.05置信區間；東江博羅水文站1960—2017年UF與UB統計曲線在信度區間內發生交叉次數較多，根據圖像性質推測可能突變的年份分別為1963，1990，1999，2007，2017 年，UF值同樣未突破0.05置信區間。綜合分析，1960－2017年珠江流域西、北、東江支流徑流雖存在各自的上升下降突變點，但UF值均未突破0.05置信區間，徑流有輕微變化但總體變化趨勢并不顯著。

圖2 1960－2017年高要、石角、博羅水文站徑流M－K 檢驗

滑動平均法因其簡單直觀性，在水文學研究中應用較為廣泛。其主要思路是將時間序列x1，x2，…，xn的幾個前期值和后期值取平均，求出新的序列yt。其計算公式如公式（7）［25］：

式中：yt——新序列；k 為滑動平均尺度；xt＋i——參與此次生成新序列的舊序列值；k＝1時yt為3 a滑動平均值，k＝2時yt為5 a滑動平均值；若Xn有趨勢性，選擇合適的k，即可通過 將其趨勢清晰地表示出來。

本文計算了珠江流域高要、石角、博羅1960—2017年徑流3 a滑動平均值與5 a滑動平均值，計算結果如圖3所示。

3個水文站3 a滑動平均值與5 a滑動平均值相差不大，均可說明徑流的變化趨勢。總體看來，徑流都呈現上升—下降—上升—下降—上升5 個變化過程。高要水文站存在3 個上升期，分別為：1960—1970，1990—1997，2010—2017 年，且以后兩個階段上升速率較快，其余時段徑流呈下降趨勢，主要的變異年份為1964，1970，1988，1997，2011年；石角水文站同樣存在3個上升期，分別為：1965—1975，1990—1997，2010—2017年，3個階段上升速率差異不大，其余時段徑流呈下降趨勢，主要的變異年份為1964，1975，1983，1990，1998，2010 年；博羅水文站上升—下降階段較為復雜，以6～8 a為時間間隔上升—下降相互交替出現，主要變異年份為1964，1975，1984，1990，1999，2003，2008年。

圖3 1960－2017年高要、石角、博羅水文站年徑流及3 a，5 a滑動平均

分析1960年以來西、北、東江突變點出現的原因，推測主要是由大洪水的出現導致，西江在1968，1986，1994，1998，2005，2016 年爆發過大洪水，北江在1968，1982，1994 年爆發過大洪水，東江在1984，2007，2016年爆發過大洪水，一般在洪水發生前后2 a內出現突變點。此外東江由于徑流小而水庫庫容較大，人類對于水庫的儲放水活動也容易產生突變點。

對M-K 趨勢分析法與滑動平均法檢測出的突變點進行對比，結果如表2所示，兩種方法檢測出的突變時刻相差不大，滑動平均法檢測出的突變點往往更多，但在兩個時間端點（1965 年之前與2010 年以后），滑動平均法由于兩端點數據缺失（3 a滑動平均只能輸出1961—2016年數據，5 a滑動平均只能輸出1962—2015年數據）檢測性能往往較差，而M-K 趨勢分析法不受端點限制，檢驗出3個水文站在時間軸右端點2016年左右均發生突變。

表2 高要、石角、博羅水文站M-K 檢驗與滑動平均法突變點對比

此外，對高要、石角、博羅3個水文站1960—2017年徑流數據作線性回歸分析，研究其58 a間總體變化特征，線性回歸方程如圖3所示，高要、石角、博羅線性分析斜率k值分別為－0.16，0.95，0.35，相關系數r值分別為0.002，0.047，0.03，根據線性回歸計算結果分析可知，西江流域徑流總體呈下降趨勢，北、東江流域徑流總體呈上升趨勢。當α＝0.05，n＝58時，查表得相關系數r＝0.259，3 個流域｜r｜值均小于r0.05，由此證明3個支流雖存在上升或下降趨勢，但變化趨勢并不明顯。

2.1.3 徑流年際枯豐變化 為找出徑流序列豐、枯水年階段性變化規律，根據水利部信息中心編制的水文預報規范，將徑流量距平百分率Ki劃分為5個級別，5個級別分別是：Ki＜－20%為枯水，－20%≤Ki＜－10%為偏枯，－10%≤Ki≤10%為平水，10%＜Ki≤20%為偏豐，Ki＞20%為豐水。其中Ki計算公式如公式（8）：

珠江流域西、北、東江控制站高要、石角、博羅3個水文站年徑流距平百分率Ki計算結果如圖4 所示，年徑流分級統計結果見表3。3個水文站徑流豐枯差異較大，高要、博羅站平水年所占比例最大，1960—2017年平水年所占比例分別為37.9%，44.8%，石角站豐水年所占比例最大，1960—2017年豐水年所占比例為48.3%，接近調查年份總年數的一半。1960—2017年高要水文站豐水年與偏豐年總計18 a，而枯水年與偏枯年總計也為18 a，其平水年比例占比最大，說明西江流域旱澇災害發生較均勻且發生率較低，水情較為穩定；石角水文站豐水年與偏豐年總計31 a，枯水年與偏枯年總計9 a，豐水年比例遠遠高于枯水年，且豐水年比例也高于平水年，說明北江流域易發生洪澇災害；博羅水文站豐水年與偏豐年總計16 a，而枯水年與偏枯年總計也為16 a，且其平水年所占比例最大，說明東江流域旱澇災害發生較均勻且發生率較低，水情比較穩定。按此方法分析西、北、東江旱澇災害發生率最小的為西江、其次為東江，最易發生旱澇災害的為北江。

圖4 1960－2017年高要、石角、博羅站徑流距平百分率K i 值

表3 1960－2017年珠江流域高要、石角、博羅水文站年徑流分級統計

高要、石角、博羅3個水文站徑流量距平百分率Ki在不同年代之間也存在差異，高要水文站豐水年主要發生在1965—1975，1994—1998年期間，枯水年主要發生在1987—1992年期間；石角水文站主要發生在1965—1977，1990—2016 年期間，其中1990—2016期間也是枯水年多發年份；博羅水文站豐水年、枯水年多發的年份都為1990—2016年，豐枯年份基本交替出現。

2.1.4 徑流年際周期性變化 本文采用Morlet小波分析對珠江流域西江、北江、東江支流1956—2013年徑流量進行周期性分析。Morlet小波能從時域和頻域兩方面進行信號分辨，能夠準確地揭示隱藏在時間序列中的多種變化周期，充分反映出系統在不同時間尺度中的變化趨勢。Morlet小波分析計算公式如公式（9）—（10）。對于給定的能量有限信號f（t）∈L（R），其連續小波變換為：

式中：Wf（ɑ，b）——小波變換系數；f（t）——一個信號或平方可積函數；ɑ，b∈R；ɑ≠0，ɑ 為伸縮尺度；b——平移參數。將小波系數的平方值在b 域上積分，得小波方差如下式：

根據公式計算得出珠江流域流域西、北、東江小波系數實部等值線見圖5。小波系數實部等值線圖能反映徑流序列不同時間尺度的周期變化及其在時間域中的分布。

根據圖5可知，珠江流域西江徑流存在25～32 a，15～20 a，3～8 a這3類尺度的周期變化特征；北江徑流存在25～32 a，10～20 a，4～10 a這3類尺度的周期變化特征；東江徑流存在25～32 a，10～20 a，4～7 a這3類尺度的周期變化特征。3個支流均存在25～32 a的周期變化，58 a間西、北、東江分別出現6，5，5次周期震蕩，該時段的周期變化較為穩定，具有全域性；第二震蕩周期基本為10～20 a，西、北、東江分別出現9，11，12次周期震蕩；第三震蕩周期基本為4～10 a，西、北、東江震蕩次數較多，平均每2～3 a就會出現一次豐季或枯季。

2.2 徑流年內變化特征

本文借鑒徑流量年內分配的向量法定義單站徑流時間分配特征的參數徑流集中度和集中期。徑流集中度指各月徑流量按月以向量方式累加，其各分量之和的合成量占年徑流量的百分數，其意義是反映徑流量在年內的集中程度。集中期是指徑流向量合成后的方位，反映全年徑流量集中的重心所出現的月份，以12個月分量和的比值正切角度表示，具體計算過程見劉賢趙等［12］文獻。計算公式如公式（11）—（13）：

式中：RCDyear——年徑流集中度；RCPyear——年徑流集中期；Ryear——年徑流量；Rx，Ry——12個月的分量之和所構成的水平、垂直分量；ri——第i月的徑流量；θi——第i月徑流的矢量角度。RCDyear越大表示徑流量越分散，越小徑流量越集中，RCPyear表示一年中最大徑流量出現的時間，表4為各個月份包含及代表角度值。

圖5 1960－2017年高要、石角、博羅水文站徑流小波分析實部等值線

珠江流域高要、石角、博羅水文站1960—2017年間10 a尺度的徑流集中度與集中期計算結果如表5所示。從徑流集中度的大小關系可以看出西江年徑流量的分散程度略高于北江，高于東江。根據年徑流集中期對應的角度，發現西江最大徑流量主要出現于7月份，北江最大徑流量主要出現在5，6 月份，東江最大徑流量主要出現6，7月份。

RCDyear與RCPyear在1960—2017 年也表現出一定的變化趨勢，1960—2009 年間年徑流集中度較為相似，2010—2017徑流集中度明顯減小，尤其以西江高要水文站減小最為劇烈，反映了最近10 a間徑流年內分配變得更加分散。而對于年徑流集中期，各個水文站變化具有不一致性，高要水文站徑流集中期對應角度雖都位于165°～195°之間，但徑流集中期對應的角度呈減小趨勢，反映了西江最大徑流月份不變，但在最大徑流月內最大徑流日期向前移動；石角水文站徑流集中期對應角度無明顯增大或減小趨勢，上下震蕩，1980—1989年徑流集中期為5月份，其余時間徑流集中期為6月份；同石角水文站相似，博羅水文站流集中期對應角度也屬于上下震蕩類，1980—1989，2010—2017年徑流集中期為6月份，其余時間徑流集中期為7月份。

同時，對1960—2017 年高要、石角、博羅水文站月度徑流作等值線圖（圖6），從圖6中可以清楚地看出徑流年內的分散程度與最大、最小徑流出現的時間。對比徑流等值線圖與徑流集中度與集中期表格，發現其演變階段具有一致性，徑流等值線圖年內徑流較分散年份，徑流集中度往往也較低，反之亦然，而徑流集中期角度的波動與等值線圖最大值的上下移動也有較好的對應性。

圖6 1960—2017高要、石角、博羅水文站徑流等值線

表4 RCPyear計算各月份包含及代表角度值 （°）

3 討論

3.1 徑流年際變化

本文從變異性、趨勢性、周期性與豐枯變化特征等方面說明了珠江流域下游主要水文站西江高要站、北江石角站、東江博羅站徑流年際變化特征。其中變異性、豐枯變化特征反映了徑流的年際波動特征，即反映徑流的穩定性，綜合變異性與豐枯變化特征，珠江流域下游徑流變異性與豐枯變化最小的都為西江，其徑流穩定性最高，而北江、東江在變異性與豐枯變化兩個特征方面排序不同，在變異性方面，東江雖高于北江但兩者徑流年際變差系數（Cv）、年際極值比（w）數值都較為相近，而在豐枯變化特征方面，北江變化程度遠遠高于東江，綜合兩方面，認為東江穩定性高于北江，穩定性由大到小為西江、東江、北江，且三江都基本存在25～32 a，10～20 a，4～10 a的徑流周期變化。而在徑流趨勢性方面，3個水文站雖存在眾多變異點，但由于M-K 分析均為突破0.05 置信區間，且線性回歸｜r｜值均小于r0.05，證明徑流無明顯趨勢性。

表5 高要、石角、博羅水文站年徑流集中度與集中期

分析徑流年際變化的原因，流域降水的多少對于徑流影響巨大［26］，特別是在人類活動影響較小的地區，其作用更為顯著，降水通過影響流域水分在下墊面的垂向和橫向的再分配，進而作用于流域徑流的變化［27］。除此之外，流域植被覆蓋率也是影響流域徑流的因素之一，但在珠江流域，其對徑流變化的貢獻并不顯著［28］。

受距離海洋的遠近及海拔的影響，珠江流域西北部長期降水分布較平均，降水極值情況發生的幾率相對較低；東南部長期降水較集中，降水極值情況發生的幾率相對較高［29］。受降水空間分布特征的影響，西江徑流較為穩定，東江、北江相對較差。除自然因素外，人類對于徑流年際變化也有一定的影響，但相對于自然作用其影響程度較小，人類主要通過流域修建大壩、影響流域植被覆蓋率等方式影響流域徑流，在珠江流域，大壩的修建是人類的主要作用方式，例如由于東江徑流量較小，且其水庫庫容較大（見表6），在新豐江水庫等的調節作用下，徑流穩定性得到改善，較少發生大洪水或極端干旱，而北江徑流量大于東江，水庫庫容也較東江較小，人類干預作用較小，穩定性較差。

表6 珠江流域較大規模水庫建設情況

3.2 徑流年內變化影響因素

上文通過集中度與集中期分析了西江高要站、北江石角站、東江博羅站徑流年內分配情況，研究結果表明，徑流集中度總體呈下降趨勢，西江、東江下降速度較快，其下降時期分別為2010—2017，1970—1979年，北江下降速度較慢，主要下降期為2010—2017年，徑流集中度的下降與流域水庫修建等人類活動關系密切，人類活動一般通過建設各類水利工程、變更土地利用類型等改變下墊面條件來間接影響流域徑流［30］，且其對于徑流年內的影響程度遠高于年際，例如通過水庫、淤地壩、水窖等蓄水工程，攔蓄地表徑流、消減洪峰流量、調節徑流分配。目前，在對中國眾多流域的研究中都發現，人類活動對于流域徑流量變化的貢獻率都遠遠超過50%［31-33］。

水庫的修建延長了匯流時間，消減了洪峰，使年內流量過程線變得平緩。20世紀60年代，大規模的水庫尚未建成，此時徑流集中度可近似代表原始狀態下的年內徑流集中程度，集中度主要受各個流域降水的影響，其由高向低依次為西江、北江、東江，2010—2017年處于人類的劇烈活動下，眾多水庫已建成或正在修建，這一階段徑流集中度由高到低依次為北江、西江、東江，2個階段西江、北江、東江徑流集中度的變化比例分別為23%，11%，28%，而計算目前三江水庫庫容所占徑流量比例，分別為25%，8%，72%，兩者大小關系較一致，說明了大壩修建在調節徑流年內分配方面有重要作用。

與徑流集中度相似，徑流集中期同樣受人類活動影響，20世紀60年代珠江流域受人類活動影響較小，徑流主要受降水控制［34］，由于各個流域間氣候存在差異，造成北江最大徑流出現時間早于東江早于西江，受人類活動的影響，1960—2017年間徑流集中期有向前移動的趨勢，但目前各個支流年內最大徑流出現的時間仍為北江早于東江早于西江。

4 結論

本文基于珠江流域西、北、東江主要控制站高要、石角、博羅水文站1960—2017年月度徑流數據，利用徑流年際變差系數（Cv）、年際極值比（w）、M-K 趨勢分析、滑動平均、集中度（期）等方法對三站長序列徑流年際變化與年內變化特征進行研究，主要結論如下：

（1）珠江流域西、北、東江徑流變差系數（Cv）值均較小，其值介于0.196～0.273之間，說明珠江流域徑流年際變率較小，徑流量穩定性較高，流變差系數（Cv）與徑流年際極值比（w）均反映珠江三條支流年徑流穩定性順序西江＞北江＞東江。

（2）M-K 趨勢分析法與滑動平均法均檢測出珠江流域徑流存在若干突變點，徑流有輕微變化但總體變化趨勢并不顯著。

（3）徑流豐枯變化表明西、東江流域旱澇災害發生較均勻且發生率較低，北江流域易發生洪澇災害。

（4）Morlet小波分析表明，近幾十年珠江流域西、北、東江徑流普遍存在25～32 a，10～20 a，4～10 a的周期變化。

（5）珠江流域徑流年內分配極不均勻，年徑流集中度（RCDyear）西江略高于北江，高于東江，且近年來集中度都有減小趨勢；年徑流集中期（RCPyear）表明西江最大徑流量主要出現于7月份，北江最大徑流在5，6月份交替，東江最大徑流量在6，7月份交替出現，近年來集中期都有向前（更早月份）移動趨勢。