王永紅, 彭 錦, 黃 暢, 吳宏旭, 蔡斯龍, 劉 鋒

(1.中國海洋大學 海洋地球科學學院, 海底科學與探測技術教育部重點實驗室, 山東 青島 266100;2.廣東省水文局, 廣州 510150; 3.中山大學 海洋工程與技術學院, 廣東 珠海 519082)

隨著社會經(jīng)濟的發(fā)展，人類活動對于河口三角洲水沙變化的影響日益增強。強人類活動導致了水沙情勢的異變，引起水道的洪潮水位、徑流量和輸沙率以及分流比等水文要素的變化，容易造成洪水災害以及河段沖淤逆轉等，從而直接影響河口三角洲社會經(jīng)濟的可持續(xù)發(fā)展。在珠江三角洲，強烈人類活動尤為顯著，例如中上游水庫建設，森林砍伐，水土保持，河網(wǎng)采砂，口門圍墾以及河道疏浚等，都會造成河口三角洲區(qū)域水沙的異變，影響大灣區(qū)的水安全。因此各個相關部門都極為重視珠江流域以及河口三角洲的水沙變化及其原因，特別是自然和人類活動貢獻的定量化，以及水沙發(fā)生異變的時間節(jié)點以及控制因素。

目前的研究手段基本可以定量區(qū)分自然因素和人類活動因素對于水沙變化的影響，認為目前人類活動對于水沙變化的貢獻率在80%以上。例如在對珠江流域1954—2018年的水文氣象數(shù)據(jù)以及高要、石角和博羅站水沙數(shù)據(jù)分析的基礎上，發(fā)現(xiàn)人類活動特別是大壩建設造成了96%的輸沙量變化，而4%是由于氣候變化[1]。在對1950—2000年的氣候變化和人類活動對珠江流域(西北東江)徑流和泥沙進行定量計算中，發(fā)現(xiàn)人類活動對三江的貢獻值不同，貢獻值從高到低依次為西東北江，其中西江的人類活動貢獻值最大為96%[2]。對1956—2013年珠江流域(西北東江)輸沙率進行研究得出2000年后人類活動對輸沙率減少作用更加顯著，貢獻率達到80%[3]。通過對高要(西江)、石角(北江)和博羅(東江)水文站1954—2011年的徑流和輸沙率的相關關系計算得出氣候變化和人類活動對輸沙率的貢獻值，認為人類活動的貢獻值從30%增加到80%[4]。

對于珠江流域以及河口三角洲水沙的階段性變化以及影響因素也有了一定的認識，例如一些學者在對西江流域(梧州站和高要站)的輸沙進行研究時發(fā)現(xiàn)，輸沙變化的主控因素是水土保持和水利建設等人類活動[5]，如西江流域(梧州站)輸沙量可以分為兩個階段，其發(fā)生階段性變化主要受流域大型水利樞紐工程的建設攔沙等控制因素影響[6]。另外，有研究認為珠江流域(西北東江)輸沙率分為3個階段，人類活動(大壩建設)是造成輸沙率變化的主要原因[7]。在對東江流域(博羅站)的輸沙量的研究中發(fā)現(xiàn)，其輸沙變化分為3個階段，其控制因素主要為水庫和采砂等人類活動[8]。三水、馬口站泥沙變化分為兩個階段，其變化的控制因素主要為上游大型水庫建設和森林覆蓋率的提高[9]。

對于珠江流域的水沙變化的研究還有很多，但大都集中在珠江流域和河口三角洲區(qū)域的水沙階段性變化趨勢研究。雖然也有一些研究指明水沙的異變，但是只是在某一個站位，例如只談了西江[10]，或者馬口三水等[11]。由于珠江流域和河網(wǎng)的復雜性，西江、北江和東江流域的氣候以及人類活動程度不同，目前還沒有給出每一個江的具體水沙異變時間節(jié)點以及明確的異變機制。主要的原因是研究者使用的數(shù)據(jù)時間長度序列不一樣，或者研究的水文站位比較單一，空間分布不足，分析方法也不相同，從而造成對于水沙數(shù)據(jù)分段以及異變研究結果的不一致性。本文收集珠江流域下游5個典型水文站的1960—2019年最新的水沙數(shù)據(jù)，利用水沙突變檢測和雙累積曲線的方法，對珠江流域和網(wǎng)河區(qū)的水沙情勢異變分別進行判斷，并且使用異變分級的方法，給出三江不同流域的水沙異變點和顯著變化點，并闡明異變的主控因子。

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

1.1 水沙數(shù)據(jù)來源

本文收集了高要、石角、博羅、馬口和三水5個典型水文站(圖1)1960—2019年年度流量與輸沙率數(shù)據(jù)，數(shù)據(jù)來源于廣東省水文局(歷史測量數(shù)據(jù))，中山大學(歷史保存數(shù)據(jù))、《中國河流泥沙公報》(2003—2019版) (http:∥www.mwr.gov.cn/)。統(tǒng)計年份見表1。珠江流域主要大壩建設數(shù)據(jù)來源于中國大壩工程學會網(wǎng)站(http:∥www.chincold.org.cn/)。

表1 五大水文站歷史水沙統(tǒng)計年份

1.2 分析方法

1.2.1 Mann-Kendall檢驗法 M-K趨勢分析法(Mann-Kendall檢驗法)是一種非參數(shù)統(tǒng)計檢驗方法，其優(yōu)點是樣本不需要遵循一定的分布，也不受少數(shù)異常值的干擾，更適用于類型變量與順序變量，計算也比較簡便，最初由曼和肯德爾提出原理并發(fā)展了這一方法。現(xiàn)已成為世界氣象組織推薦并已廣泛使用的非參數(shù)檢驗方法，適用于水文、氣象等非正態(tài)分布的數(shù)據(jù)。其計算公式如公式(1—4) ，設長時間序列(流量、輸沙率)為x1，x2，…，xn;Sk表示第i個樣本xi>xj(1≤j≤i)的累計數(shù)：

圖1 珠江流域下游5個典型水文站和上游水庫位置

(1)

其中：

(2)

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

(3)

將Sk標準化為：

(4)

式中：UF1=0,給定顯著性水平α,若|UF|>Uα,表明序列存在明顯的趨勢變化。將此方法引用到反序列,序列順序變?yōu)閤n,xn-1,…,x1,反序列標準化由UBk表示,UBk,UFk可組成兩條曲線,在這里臨界值U0.05為±1.96,設定α值為0.05。當UBk,UFk兩條曲線出現(xiàn)交點,則交點處可能是變化點[11]。

1.2.2 Pettitt檢驗法 Pettitt檢驗法依據(jù)長時間序列趨勢性的變化為基礎，確定序列突變點。假設序長度為T，變化點為t，把假設序列劃分為前后兩段，即兩樣本為x1,x2,…,xt和xt+1,xt+2,…,xT, Pettitt突變檢驗的本質是檢驗序列中的兩個樣本是否來自于同一個樣本的Mann-Whitney統(tǒng)計量。

統(tǒng)計量Ut的計算公式為：

(5)

if(Xi-Xj)>0 sgn(Xi-Xj)=1

if(Xi-Xj)=0 sgn(Xi-Xj)=0

if(Xi-Xj)<0 sgn(Xi-Xj)=-1

(6)

在檢驗序列時，根據(jù)順序統(tǒng)計量理論引用了K-S兩樣本檢驗，檢驗兩個樣本累積分布函數(shù)的最大差值而求得的累積概率為P，即突變點位置。

Kt=max1≤t≤T|Ut|

(7)

(8)

如果累積概率越接近于1，說明存在突變點的趨勢越顯著；令α為置信度，若p>Pα，則突變點趨勢顯著。當α<0.5時，表示存在有效突變點；當α<0.05時，表示存在顯著突變點；當α<0.01時，表示存在極顯著突變點[12-13]。

1.2.3 判別水沙異變的程度 本研究根據(jù)Mann-Kendall檢驗和Pettitt檢驗的結果，對水沙變化程度進行層級劃分。(1) 當Mann-Kendall檢驗的UF曲線超過置信度區(qū)間且在0.05顯著水平內(nèi)與UB曲線交點，此時Pettitt檢驗曲線超過0.01顯著水平線時定義為異變點; (2) 當Mann-Kendall檢驗的UF曲線在置信度區(qū)間且在0.05顯著水平內(nèi)與UB曲線交點，同時Pettitt檢驗曲線僅超過0.5顯著水平線時定義為顯著變化點; (3) 當Mann-Kendall檢驗的UF曲線在置信度區(qū)間且在0.05顯著水平內(nèi)與UB曲線交點，同時Pettitt檢驗曲線未超過0.5顯著水平線時只表示為一般趨勢性變化點。本文結合水沙比(流量與輸沙率的比值，m3/kg)和水沙雙累積曲線的變化趨勢，重點分析水沙的異變點和顯著變化點。

2 結果與分析

根據(jù)Mann-Kendall檢驗、Pettitt檢驗以及雙累積曲線的方法，對西江的高要站和馬口站、北江的石角和三水站，以及東江的博羅站進行了異變點和顯著變化點的確定(圖2，表2)。異變點表明此點前后水沙發(fā)生了異變，顯著變化點前后表明水沙發(fā)生了顯著的變化，但不是異變(表2)。

2.1 西江高要站和干流馬口站水沙異變的主控因子

根據(jù)西江高要站徑流量的Mann-Kendall檢驗和Pettitt檢驗，Mann-Kendall檢驗中UF曲線在置信度區(qū)間且在0.05顯著水平內(nèi)與UB曲線有交點，同時Pettitt檢驗曲線未超過0.5顯著水平線，因此可以認為高要站的流量沒有異變點(圖2)。多年的徑流量變化范圍為3 367～10 198 m3/s，平均值為6 912 m3/s(圖3，表3)。西江高要站輸沙率Mann-Kendall檢驗顯示UF曲線超過置信度區(qū)間，且在0.05顯著水平內(nèi)與UB曲線交點，表明其異變點為2001年。Pettitt高要站輸沙率檢驗曲線超過0.01顯著水平線，確定的異變點為2000年。結合輸沙率過程圖發(fā)現(xiàn)在2000年后輸沙率開始進入一個新階段，因此判斷2000年為西江水沙異變點(圖3，表3)。異變發(fā)生之前，1960—2000年的輸沙率變化范圍為521～4 167 kg/s，平均值為2 249 kg/s，水沙比為3.4 m3/kg(表3)。異變之后2000—2019年的輸沙率變化范圍為255～1 823 kg/s，平均值為764 kg/s，水沙比為11 m3/kg。異變之后輸沙率約為異變前的三分之一。說明輸沙率在量級上已經(jīng)發(fā)生了根本的變化。從雙累積曲線也可看出在2000年前后斜率發(fā)生明顯變化，斜率從0.33變?yōu)?.10(表3)。

高要站在2000年前后發(fā)生異變的主要原因是西江上游水庫的建設。1960—2000年，輸沙率變化不大，主要受到氣候的控制。在1960—1985年，氣候和人類活動貢獻值分別為83%和17%[2]，此時輸沙率由于森林大量砍伐導致流域水土流失加劇而增大。此時廣西的水土流失面積在20世紀50—60年代為12 000 km2，而80年代為30 600 km2，是50—60年代的近3倍。而同期廣東省水土流失面積則從7 444 km2增加到了 17 070 km2，是原來的約2.3倍[14-15]。由于西江流域在廣西境內(nèi)集水面積共計20.24萬km2，占全流域集水面積的85.7%，因此廣西壯族自治區(qū)和廣東省的水土流失面積增加，使得高要站泥沙含量呈現(xiàn)上升趨勢。

而在1985—2000年，人類活動開始增強，但輸沙率仍然以氣候控制為主。氣候和人類活動貢獻值分別為63%和37%[2]。此階段水土流失程度減弱，廣西壯族自治區(qū)20世紀90年代的水土流失面積為28 100 km2，比80年代下降了8.2%，變化不大[14-15]。相反此時水庫大壩的影響增強，1985—2000年，西江上游主要修建的水庫為魯布革水庫、巖灘水庫和天生橋水庫(表4，圖1)。1992年巖灘水庫建成后高要站輸沙率則由1981—1991年的7 658萬t/a降到了1992—2002年的6 310萬t/a，降低了18%[16-17]。巖灘水庫的泥沙沉積速率為35 Mt/a，高要站的泥沙量減少了10 Mt/a[18]。因此可以看出巖灘水庫對于泥沙減少起到重要的調節(jié)作用[17]。不過巖灘水庫庫容較少，因此對于泥沙的調節(jié)作用仍然有限，沒有在1992年前后引起輸沙率的異變。

2000—2019年，雖然前期巖灘水庫調節(jié)作用有限，但在1997年建成的天生橋水庫的庫容量是巖灘水庫的3倍多(表4)，建成后經(jīng)過一定時間的蓄水和攔沙過程，導致2000年后的輸沙率大幅下降，致使輸沙率發(fā)生異變，之后在水庫的調節(jié)作用下，輸沙率保持穩(wěn)定波動。其中天峨站(巖灘水庫上游150 km處的一個水文站)在天生橋水庫建成后的輸沙量減少了30 Mt/a，發(fā)生了較大的變化，也可以幫助說明天生橋水庫的巨大調節(jié)作用[17]。值得注意的是2006年在西江上游建設了更大型的龍灘水庫(表4，圖1)，其庫容量約為天生橋水庫的2.5倍。水庫建成后使得2006—2019年輸沙減小之后保持基本穩(wěn)定(圖3)，但是并沒有引起輸沙率的異變。這是因為天生橋水庫位于西江上游(圖1)，而龍灘水庫位于天生橋水庫下游，天生橋水庫已經(jīng)攔截大部分泥沙，從而導致在2000年輸沙率發(fā)生異變，而龍灘水庫庫容量大，但建成后對下游的輸沙影響不如天生橋水庫，沒有引起輸沙率發(fā)生異變。

圖2 1960-2019年珠江流域五站位年流量與年輸沙率MK檢驗和Pettitt檢驗

表2 珠江流域五站點異變的判定

輸沙率異變后森林覆蓋面積持續(xù)上升，例如廣西壯族自治區(qū)森林覆蓋率由2000年的39%增長到2005年的53%[11]，因此水土保持使得輸沙率持續(xù)小幅度降低。另外，1995—2014年高要站附近也發(fā)生了挖沙的行為[19]，也對輸沙率減少起到一定的作用。此階段人類活動的貢獻值大幅上升，氣候和人類活動貢獻值分別為4%和96%[2]。

據(jù)馬口站流量的Mann-Kendall檢驗和Pettitt檢驗，發(fā)現(xiàn)在1986年Mann-Kendall檢驗的UF曲線在置信度區(qū)間且在0.05顯著水平內(nèi)與UB曲線交點，同時Pettitt檢驗曲線僅超過0.5顯著水平線，因此可以認為馬口站的流量在1986年存在顯著變化點(圖2)。1960—1986年的流量變化范圍為3 818～9 978 m3/s，平均值為7 498 m3/s，水沙比為3.3 m3/kg，水沙變化的相關系數(shù)為0.52。雙累積曲線斜率為0.33。顯著變化之后1986—2019年的流量變化范圍為4 403～9 602 m3/s，平均值為6 777 m3/s，較異變前減少9.6%。水沙比為8.0 m3/kg，水沙變化的相關系數(shù)為0.08，雙累積曲線斜率為0.17(表3)。馬口在1986年發(fā)生的徑流顯著變化，1986年初期發(fā)生的徑流量變化和當時劇烈的氣候動蕩變化有關[4]。

西江馬口站輸沙率的Mann-Kendall檢驗線在0.05顯著水平之外，顯示的異變點為1999年；Pettitt檢驗線0.01顯著水平線之外，異變點為1999年。因此判斷1999年為西江干流馬口站水沙異變點(圖3)。異變發(fā)生之前，1960—1999年的輸沙率變化范圍為506～4 186 kg/s，平均值為2 306 kg/s，水沙比為3.5 m3/kg，水沙變化的相關系數(shù)為0.36(表3)。異變之后1999—2019年的輸沙率變化范圍為202～1 680 kg/s，平均值為745 kg/s，水沙比為10.8 m3/kg，水沙變化的相關系數(shù)為0.46，輸沙率約為異變前的三分之一。從雙累積曲線可看出在1999年前后斜率發(fā)生明顯變化，斜率為從0.34變?yōu)?.10(表3)。因此馬口站的水沙變化基本和高要站相似。其不同階段輸沙異變的原因也和高要站一致。因此馬口的輸沙率異變主要受到上游水庫的建設誘發(fā)。不過馬口站的輸沙率異變的時間略微早于高要站，可能和馬口站和三水站的水流分水比在1988—1997年期間的顯著變化有關[11]。馬口站的采砂活動主要集中在1995—2006年[19]，因此采砂會對輸沙率有一定的影響，但并不是主要因素。

2.2 北江的石角站和干流三水站水沙異變和主控因子

根據(jù)徑流和輸沙率異變和顯著變化的診斷標準，和高要以及馬口站的判別方式一樣，確定北江石角站徑流量沒有異變點。多年的徑流量變化范圍為515～2 280 m3/s，平均值為1 335 m3/s(圖3)。石角站輸沙率沒有異變點，但是存在顯著變化點，結合輸沙率過程圖發(fā)現(xiàn)在1999年后輸沙率開始進入一個新階段，因此判斷1999年為北江石角輸沙率顯著變化點(圖3，表3)。1960—1999年，輸沙率整體平穩(wěn)，呈現(xiàn)上升趨勢，輸沙率平均值為190 kg/s，水沙比為7.9 m3/kg，水沙變化的相關系數(shù)為0.68，雙累積曲線斜率為0.15(表3)。1965年之前，幾乎不受人類活動影響，氣候貢獻90%，人類活動貢獻10%[2]。該階段降雨量減少，造成徑流量和輸沙率的降低[2]。1965年后人類活動逐漸增強。例如1971年修建的南水水庫，庫容量約為12億m3(表4)，引起了隨后的輸沙率的急劇減少(圖2)。

注：虛線為顯著變化點，實線為異變點。

表3 珠江流域五站位流量與輸沙率之間的相關關系

表4 珠江流域較大規(guī)模水庫建設情況

雖然繼續(xù)有水庫的建設，但是水庫庫容量都小于2億m3，對泥沙的攔截作用有限。同時森林砍伐引起的水土流失作用超過了水庫攔沙的影響[11]，覆蓋北江流域的廣東省水土流失面積從20世紀 50—60年代的7 444 km2增加到了80年代的17 070 km2，增加了1.3倍[14-15]，比1950s初期增加了2倍多[9]，使得河口三角洲區(qū)域輸沙率增加。而水庫的建設所攔截的泥沙量小于水土流失的量，因此泥沙含量呈現(xiàn)上升趨勢。1985—1993年珠江的降雨量、徑流量和輸沙率處于劇烈的減小階段[4]，本階段輸沙率受到氣候影響的貢獻大致為88%，人類活動影響的貢獻為12%[2]。因而輸沙率在該階段前期呈下降趨勢，但是輸沙率在此階段后期(1989年)呈現(xiàn)一定上升，主要還是受到氣候和徑流的控制。

1999—2019年，在1999年建設的飛來峽大壩(庫容約為20億m3)，造成了石角站輸沙率的顯著變化，這一階段的平均輸沙率為134 kg/s，水沙比為12.3 m3/kg，水沙變化的相關系數(shù)為0.67，雙累積曲線斜率為0.11(表3)。因此1999年的水庫建設造成了石角輸沙率的顯著變化。本階段氣候貢獻為36%，人類活動貢獻為64%[2]。這一階段國家水土保持政策使得森林覆蓋面積有一定增長，20世紀80年代中期開始，在國家水土保持政策的引導下，北江森林覆蓋率有一定的提高，北江區(qū)域面積為4.8萬 km2，森林覆蓋率由20世紀80年代3.5萬 km2增加到20世紀90年代3.7萬 km2，森林覆蓋率由73%上升為77%[7]，上升速度并不快，因此水土流失得到一定控制，對輸沙率的減少有一定的作用，但不是主控因子。另外，石角站從2003—2014年有持續(xù)的挖沙活動，也對輸沙率的減小有一定貢獻。雖然水沙存在顯著性變化，但是不是異變，水沙相關性一直相對較好。

根據(jù)異變和顯著變化的診斷標準，北江干流三水站徑流量的異變點為1991年。異變發(fā)生之前，1960—1991年的流量變化范圍為298～1 910 m3/s，平均值為1 213 m3/s。異變之后1991—2019年的流量變化范圍為1 139～2 938 m3/s，平均值為1 952 m3/s，較異變前增加61%。水沙比為原來的2倍，水沙變化的相關系數(shù)變差，雙累積曲線斜率減小(表3)。

北江三水站異變點為2003年(圖2)。異變發(fā)生之前，1960—2003年輸沙率變化范圍為18～578 kg/s，平均值為292 kg/s(表3)。異變之后2003—2019年的輸沙率變化范圍為44～296 kg/s，平均值為155 kg/s，水沙比為原來的2.3倍，水沙變化的相關性變差。泥沙量約為異變前的二分之一，說明輸沙率在量級上已經(jīng)發(fā)生了根本的變化。從雙累積曲線可看出在1991年和2003年前后斜率發(fā)生顯著變化，斜率從0.21變?yōu)?.09(表3)。

雖然三水站主要的徑流來自北江，但是其異變現(xiàn)象并不和石角站一致，石角站徑流量并不存在異變點和顯著變化點。三水站1991年前徑流量基本保持穩(wěn)定，徑流量突然增加的原因是由于三水站在1989—2006年發(fā)生了大量的采砂活動[19]，而北江河網(wǎng)區(qū)的采砂量遠遠大于西江網(wǎng)河區(qū)，致使北江網(wǎng)河區(qū)的河床出現(xiàn)大量下切，導致三水站分流比升高，更多的徑流流入到三水站，引起了徑流量的突然增大，是以前的1.6倍，因而發(fā)生異變。而此時我們發(fā)現(xiàn)馬口站并沒有同時發(fā)生徑流異變，主要是因為馬口站的徑流量較大，例如20世紀80年代之后，馬口站年徑流量減少12%，而三水站年徑流量增加54%[11]，因此由于和三水站分水分沙引起的變化量并不足以使得馬口站產(chǎn)生徑流量的異變。

同樣，三水站的輸沙率并不和徑流量異變時間一致，而是在2003年發(fā)生異變，1960—2003年水沙變化的相關系數(shù)為0.56，異變后2003—2019年水沙相關性變差，相關系數(shù)為0.30(表3)。這是因為雖然在三水站附近的挖沙活動集中于1989—2006年[19]，加上1999年的飛來峽水庫，實際上已經(jīng)造成了三水站輸沙率的持續(xù)降低，但三水站的輸沙率在1999—2003年仍然保持一定的水平(圖3)。而2003年，三水站的輸沙率有了斷崖式下降，并且在2003年之后輸沙率都在較低水平震蕩，說明2003年后的異變主要是這段時間的上游來沙量減少所致。2000年左右，廣東省將水土保持面上資金納入省財政專項，每年固定投入1 500萬元。近幾年來，廣東省每年綜合治理水土流失面積約400 km2[20]。根據(jù)時間點變化，判斷上游的水土保持是三水站輸沙率發(fā)生異變的主要原因。另外2003—2004年遭遇極其干旱的天氣，降雨量只有1 423～1 315 mm，鄰近的2001年、2002年、2005年、2006年的降雨量都在1 766～2 110 mm[21]，加劇了上游來沙量的減少。

2.3 東江博羅站的水沙異變和主控因子

根據(jù)異變和顯著變化的診斷標準，東江博羅站徑流量不存在異變點。多年的流量變化范圍為284～1 313 m3/s，平均值為741 m3/s，水沙比為16.4 m3/kg。東江博羅利用輸沙率異變點為1988年(圖2)。異變發(fā)生之前，1960—1988年的輸沙率變化范圍為10～161 kg/s，平均值為87 kg/s，水沙比為9.7 m3/kg。異變之后1988—2019年的輸沙率變化范圍為7～128 kg/s，平均值為43 kg/s，水沙比為22.6 m3/kg。泥沙量約為異變前的二分之一。說明輸沙率在量級上已經(jīng)發(fā)生了根本的變化。從雙累積曲線可看出在1988年前后斜率發(fā)生顯著變化，斜率從0.11變?yōu)?.06(表3)。

輸沙率異變之前的1960—1988年，此階段輸沙率變化較為平穩(wěn)，輸沙率有略微下降。水沙變化的相關系數(shù)為0.60(表3)。此階段主要由人類活動主導，人類活動的貢獻為89%[2]。河流輸沙隨著流域水庫累計庫容的增加而減少。主要受到上游新豐江水庫(1962年)、楓樹壩水庫(1973年)和白盆珠水庫(1985年)的影響(表4)。3個水庫的共同作用下雖然使輸沙率減小，但是總的庫容量并不大，對于白盆珠水庫來說，其多年平均流量為37.5 m3/s，而博羅站本身徑流量為741 m3/s，遠遠大于其調節(jié)流量[22-23]。從輸沙率的變化過程圖來看(圖3)，輸沙率在2004—2006年較高，也可以證明水庫的調節(jié)作用有限，因此水庫建設對于博羅站輸沙率異變影響不大，不是異變的主要原因。而在1988年后，1988—2003年人類活動的貢獻為從89%降為75%[2]。在國家水土保持政策引導，大范圍的治山治水、開發(fā)“四荒”以及退耕還林，水土流失得到較好控制。自20世紀90年代以來，共治理東江中上游水土流失區(qū)域961.3 km2，植被覆蓋率從治理前的30%提高到65%，入河泥沙則從治理前800～1 000萬t/a降到350～400萬t/a[23]。2000年左右，廣東省每年綜合治理水土流失面積約400 km2[20]。此階段并沒有大型水庫建設，根據(jù)時間點變化，水土保持是造成博羅站的徑流量的異變的主要控制因素。由于在博羅站附近1998—2008年階段才有大規(guī)模持續(xù)的挖沙行為[19]，特別是2002—2005年間，新角水閘—譚公廟河段(博羅站上游)河床采砂量高達1 800萬 m3，河床平均下切5.7 m；受此影響博羅站同流量級水位下降嚴重，潮汐動力顯著增強[24]，但博羅站的輸沙率在1988年就已經(jīng)發(fā)生了異變，因此挖沙對博羅站輸沙率發(fā)生異變的影響很小。

3 討 論

從研究結果看，長序列的最新的水沙數(shù)據(jù)，以及對前人大量的工作成果總結，更加全面準確地認識珠江流域以及河口三角洲發(fā)生的水沙異變及其原因，基本解決了由于數(shù)據(jù)時間長度序列不一樣、站位不同以及分析方法不同造成的計算和判別結果不一致的問題。特別是西江和北江及其干流、以及東江的同步數(shù)據(jù)分析，更加清晰地理解這些地區(qū)的水沙異變的時空變化及其原因。

水沙異變主要是在徑流量和輸沙率持續(xù)變化后由于大型的事件導致兩者發(fā)生了數(shù)量上本質性的改變，除了三水的徑流量發(fā)生異變后徑流量大幅上升，異變后徑流量平均值是異變前的1.6倍，其他各站徑流量保持平穩(wěn)。各站的輸沙率發(fā)生異變或者顯著變化后普遍減少，例如高要站和馬口站輸沙率發(fā)生異變后，輸沙率只有異變之前的三分之一。三水和博羅站發(fā)生異變后，輸沙率只有異變之前的二分之一。

近70 a來，珠江流域和河口三角洲的徑流量除了個別站位，基本仍然受到氣候的控制，而泥沙則受到人類活動的控制，上游水庫建設、水土保持、以及河網(wǎng)區(qū)挖沙仍然是影響不同流域河流泥沙異變的主要原因。不過在珠江流域及河口三角洲區(qū)域，三江及其干流的水沙異變情勢并不相同。在西江高要、北江石角、東江博羅，以及西江干流馬口和北江干流三水，只有三水站的徑流在1991年發(fā)生異變，主要是由于人類挖沙導致三水河槽剖面迅速下降，經(jīng)過思賢滘的分水比增大引起了徑流量異變。這和前人根據(jù)1954—2011年徑流量和輸沙率的數(shù)據(jù)確定三水站徑流1990年發(fā)生突變相差不大[11]。但是前人認為馬口站在2002年前后徑流量減少12%，發(fā)生異變，但本文根據(jù)判別標準，認為這種量級的變化并不構成異變，因此馬口站徑流只有顯著變化，并沒有發(fā)生異變，主要是因為馬口站的徑流量大，減少的水量不至于引起其異變。

5個站位有4個站位的輸沙率量發(fā)生了異變，而每個站位的異變時間和原因并不相同。高要站和馬口站在2000年和1999年發(fā)生了輸沙率的異變，主控制因子為上游水庫和大壩的建設。這與前人的研究結果相差不大。例如在西江流域(梧州站)以1960—2015年年輸沙量為研究對象，采用R/S法和Mann-Kendall檢驗法分析確定2000年為異變點[6]。對西江高要站1957—2016年的輸沙數(shù)據(jù)通過滑動T檢驗、有序聚類法和累積距平法綜合判斷西江高要站輸沙量的異變點為1999年[25]。

北江干流的三水站和東江博羅站的輸沙率分別在2003年和1988年發(fā)生異變，主要是由于三水站附近挖沙以及東江上游水土保持造成的。只有石角站沒有發(fā)生異變，但是在1999年發(fā)生了顯著變化，1999年的顯著減小是由于上游大壩的建設引起的。前人根據(jù)1957—2000年輸沙量計算確立了石角站1982年為輸沙量的明顯變化點[25]，而1983年和1993年為兩個水沙變化節(jié)點[5]，認為其主要控制因素是采砂活動和水庫建設，和本文的研究結果有所差異，這些節(jié)點只能是一般趨勢性變化節(jié)點。前人研究表明博羅站輸沙量(1954—2011年)在20世紀80年代以來下降是一突變現(xiàn)象，并確定1984年前后輸沙量發(fā)生異變，異變后輸沙量減少46%[11]，這與本文結果相差不大。

另外，馬口三水的徑流量和輸沙率的異變和其分流分沙比的異變并不相同。前人研究表明由于北江下游河道下切幅度大于西江河道下切幅度，使馬口三水站分流比于1988年發(fā)生異變、分沙比1992年發(fā)生異變[11]。不過三水的徑流量異變發(fā)生在1991年，而輸沙率異變發(fā)生在2003年。因此徑流量發(fā)生異變是在馬口三水的分流比異變后發(fā)生，主要和挖沙導致的河道下切有關。但是輸沙率在2003年發(fā)生異變，卻和分沙比1992年的異變沒有關聯(lián)，因為三水站輸沙率的變化主要受到上游水土保持的劇烈變化有關，其分沙比引起的輸沙率變化還不足以引起異變。

4 結 論

(1) 根據(jù)Mann-Kendall檢驗法和Pettitt檢驗法兩種診斷方法，將水沙變化分為異變和顯著變化兩個層級。當Mann-Kendall檢驗的UF曲線超過置信度區(qū)間且在0.05顯著水平內(nèi)與UB曲線交點，此時Pettitt檢驗曲線超過0.01顯著水平線時定義為異變點。而當Mann-Kendall檢驗的UF曲線在置信度區(qū)間且在0.05顯著水平內(nèi)與UB曲線交點，同時Pettitt檢驗曲線僅超過0.5顯著水平線時定義為顯著變化點。

(2) 對于流量來說，馬口站流量在1986年發(fā)生顯著變化，主要受氣候控制；顯著變化前水沙比為3.3 m3/kg，顯著變化后水沙比為8.0 m3/kg，水沙相關系數(shù)從0.52變?yōu)?.08；雙累積曲線斜率也從0.33變?yōu)?.17。三水站的流量在1991年發(fā)生異變，異變后徑流量較異變前增加61%，主要是由于人類挖沙導致三水剖面迅速下降，經(jīng)過思賢滘的分水比增大引起了徑流量異變。異變前水沙比為5.2 m3/kg，異變后水沙比為10.5 m3/kg，水沙相關系數(shù)從0.61變?yōu)?.44；雙累積曲線斜率也從0.22變?yōu)?.12。其他站位沒有發(fā)生徑流異變以及顯著變化；

(3) 5個站位的輸沙量只有石角站沒有發(fā)生異變，但是在1999年發(fā)生了顯著變化；石角站發(fā)生輸沙率顯著變化后，較第一階段減少29.5%。水沙比從7.9 m3/kg變?yōu)?2.3 m3/kg；雙累積曲線斜率從0.15變?yōu)?.11。受到上游水庫的建設控制。

(4) 除石角外，其余4個站位的輸沙率量都發(fā)生了異變，而每個站位的異變時間和原因并不相同。高要站和馬口站在1999—2000年發(fā)生了輸沙率的異變，輸沙率只有異變之前的三分之一，異變前后兩站水沙比變化幅度相似，從3.4～3.5 m3/kg異變?yōu)?1.0～10.8 m3/kg；雙累積曲線斜率兩站變化幅度相同從0.34變?yōu)?.10，主要受到上游水庫和大壩的建設的控制。北江干流的三水站和東江博羅站的輸沙率分別在2003年和1988年發(fā)生異變，主要是由于上游水土保持造成的。異變后輸沙率只有異變之前的二分之一，異變前后兩站水沙比變化幅度相似，從5.6～9.7 m3/kg異變?yōu)?3.4～22.6 m3/kg；雙累積曲線斜率也從0.11～0.21變?yōu)?.06～0.09。