關(guān)河水庫流域水文序列突變?cè)\斷與非一致性研究

陳 鑫，劉艷麗，刁艷芳，關(guān)鐵生，劉翠善，王國慶，金君良，鮑振鑫，賀瑞敏

(1.南京水利科學(xué)研究院水文水資源與水利工程科學(xué)國家重點(diǎn)實(shí)驗(yàn)室，南京 210098； 2.水利部應(yīng)對(duì)氣候變化研究中心，南京 210029；3．山東農(nóng)業(yè)大學(xué)水利土木工程學(xué)院，山東 泰安 271018)

0 引 言

水文序列是開展水文學(xué)研究與相應(yīng)設(shè)計(jì)規(guī)劃的基礎(chǔ)。近年來，受氣候變化以及人類活動(dòng)的影響，水文循環(huán)的物理機(jī)制以及循環(huán)過程發(fā)生了顯著的變化，水文序列也呈現(xiàn)出非一致性的特點(diǎn)[1]，即在不同時(shí)間觀測(cè)下的水文序列會(huì)呈現(xiàn)出不同的水文特征。由于水文序列的結(jié)果直接影響到水利設(shè)施的運(yùn)行調(diào)度與風(fēng)險(xiǎn)防護(hù)，因此，研究水文序列非一致性的診斷方法對(duì)于研究水文序列在變化環(huán)境下的特征以及相應(yīng)水文模型的構(gòu)建具有較為深刻的意義。

對(duì)于水文序列的非一致性，目前的研究主要集中在水文序列的趨勢(shì)性，變異性以及周期性上[2]，其中又以水文序列的變異性研究為主。對(duì)于水文序列的趨勢(shì)性檢驗(yàn)，目前已存在多種方法，包括相關(guān)系數(shù)法，Mann-Kendall秩次相關(guān)檢驗(yàn)，Hurst系數(shù)法等；對(duì)于突變點(diǎn)的檢驗(yàn)，一般分為參數(shù)以及非參數(shù)方法[3]，其中參數(shù)方法有滑動(dòng)t檢驗(yàn)[4]，滑動(dòng)F檢驗(yàn)等[5]，非參數(shù)方法有R/S分析法[6]，Mann-Kendall檢驗(yàn)法[7,8]，有序聚類分析法[9]，Yamamoto法[10]等；對(duì)于周期性檢驗(yàn)，目前常采用小波分析法[11]，方差分析法[12]以及功率譜分析法等。

同時(shí)，由于氣候變化以及人類活動(dòng)的影響，水文序列開始呈現(xiàn)出一定程度的自相關(guān)性。如果序列中存在正向的自相關(guān)性，那么序列的趨勢(shì)顯著性也就會(huì)被放大，使得原本趨勢(shì)不顯著的序列被識(shí)別為趨勢(shì)變化顯著的序列，從而影響最終的檢驗(yàn)結(jié)果[13]。因此，應(yīng)對(duì)序列進(jìn)行預(yù)制白處理，以消除序列中的其相關(guān)性成分，使得檢驗(yàn)結(jié)果更為準(zhǔn)確。目前常用的方法包括前置移除法以及參數(shù)校正法[14]，而以去趨勢(shì)預(yù)置白法(TFPW)為代表的前置移除法在降水、蒸發(fā)以及徑流分析方面已得到了廣泛的應(yīng)用[15-17]，且取得了較好的研究成果，基于去趨勢(shì)預(yù)置白方法的突變檢驗(yàn)方法也因其對(duì)于自相關(guān)性的處理較好，漸漸被學(xué)者所接受。

為了系統(tǒng)了解氣候變化以及人類活動(dòng)影響下關(guān)河水庫流域水文時(shí)間序列的非一致性變化特點(diǎn)，本文采用了Mann-Kendall秩次相關(guān)檢驗(yàn)法檢驗(yàn)了各站降水以及徑流序列的趨勢(shì)性，并采用與TFPW結(jié)合的3種突變檢驗(yàn)方法來確定其突變點(diǎn)，同時(shí)利用Morlet小波對(duì)序列的周期性進(jìn)行分析，以達(dá)到對(duì)流域內(nèi)水文序列變化的系統(tǒng)性分析。

1 水文時(shí)間序列非一致性研究方法

1.1 去趨勢(shì)預(yù)置白方法

本文將采用去趨勢(shì)預(yù)置白方法(TPFW)，以消除序列自相關(guān)性對(duì)結(jié)果的影響，從而準(zhǔn)確地對(duì)關(guān)河水庫的水文序列進(jìn)行分析。其原理如下：

設(shè)某一待檢序列為Xt(t=1，2，…，n；n為序列長度)，則有:

(1)

Yt=Xt-β·t

(2)

(3)

(4)

1.2 水文序列突變點(diǎn)檢驗(yàn)方法

(1)Mann-Kendall趨勢(shì)檢驗(yàn)法。Mann-Kendall趨勢(shì)檢驗(yàn)法是目前水文領(lǐng)域應(yīng)用較為廣泛的非參數(shù)檢驗(yàn)方法，其通過序列的對(duì)偶數(shù)xi

(5)

(6)

按照時(shí)間序列的逆序列重復(fù)上述過程，同時(shí)使得UBk=-UFk,UB1=0，并分別畫出兩條曲線。若UFk>0，則序列存在上升趨勢(shì)，反之，序列存在下降趨勢(shì)。如果UFk和UBk兩條曲線出現(xiàn)交點(diǎn), 且交點(diǎn)在臨界線之間, 那么交點(diǎn)對(duì)應(yīng)的時(shí)刻便是突變開始的時(shí)間。

(2)Pettitt檢驗(yàn)法。Pettitt檢驗(yàn)法是由Pettitt提出的基于Mann-Whitney統(tǒng)計(jì)量的變異點(diǎn)檢驗(yàn)方法，其基本公式如下：

(7)

Kt,N=max|Ut,N|

(8)

(9)

式中：p為統(tǒng)計(jì)量；sgn為符號(hào)函數(shù)，當(dāng)xt-xi>0時(shí)，sgn取1，當(dāng)xt-xi=0時(shí)，sgn取0，當(dāng)xt-xi<0時(shí)，sgn取-1。若存在t時(shí)刻使得Kt,N=max|Ut,N|，且p≥0.95,則認(rèn)為該點(diǎn)為序列的顯著變異點(diǎn)。值得注意的是，該方法雖在單一突變點(diǎn)的檢驗(yàn)上具有較好的準(zhǔn)確性，但在多突變點(diǎn)檢驗(yàn)的表現(xiàn)較差。

(3)有序聚類法。有序聚類法通過尋求最優(yōu)分割點(diǎn)以推求突變點(diǎn)，其公式如下：

(10)

(11)

Ψ0(t)=π-1/4eiw0te-t2/2

S=minS0(γ)=min(Vγ+Vn-γ),(2≤γ≤n-1)

(12)

式中：min表示取極小值，當(dāng)S取極小值時(shí)，對(duì)應(yīng)的γ為最優(yōu)二分割點(diǎn)，以此作為突變點(diǎn)。

1.3 Morlet小波

由于Morlet小波在時(shí)間與頻率的局部化之間存在較好的平衡[18]，因此本文采用Morlet小波函數(shù)作為小波變化的基函數(shù)，研究水文序列的周期性。Morlet小波公式為：

Ψ0(t)=π-1/4eiw0te-t2/2

(13)

式中：t為時(shí)間；w0為無量綱頻率。

2 關(guān)河水庫流域水文序列診斷

2.1 關(guān)河水庫概況

關(guān)河水庫位于漳河流域的濁漳河北源，控制流域面積1 783 km2，其控制流域如圖1所示。其地理位置在長治市北側(cè)的五項(xiàng)縣城東五華里的關(guān)河峽口處，屬于濁漳河的控制性工程。由于受到了黃土高原季風(fēng)的影響，關(guān)河水庫流域水域?qū)儆诖箨懶詺夂颍瑫r(shí)由于流域內(nèi)的植被覆蓋及水土保持工程較差，流域內(nèi)的水土流失問題也較為嚴(yán)重。氣候方面，關(guān)河水庫所處流域因受黃土高原季風(fēng)的影響，多年平均降雨量在550 mm左右， 7月平均氣溫達(dá)到29.3 ℃，1月平均氣溫為-13.5 ℃，溫差較大。

為系統(tǒng)分析關(guān)河水庫流域降雨以及徑流的變化及其成因，本文選擇關(guān)河水庫以上流域中的5個(gè)雨量站和2個(gè)水文站的實(shí)測(cè)資料，分別利用3中突變分析方法對(duì)經(jīng)TFPW處理后的站點(diǎn)數(shù)據(jù)進(jìn)行分析，最終確定流域內(nèi)降雨以及徑流的突變點(diǎn)位置?？紤]到關(guān)河水庫站以及石棧道站控制流域下墊面條件不同的情況，本文以石棧道站作為受氣候變化條件單一影響下的時(shí)間序列，以關(guān)河水庫站為受氣候變化以及人類活動(dòng)綜合影響的時(shí)間序列，分別對(duì)兩種條件干預(yù)下的徑流序列變化進(jìn)行對(duì)照分析。

圖1 關(guān)河水庫控制流域Fig.1 Watershed map of study area

2.2 序列的趨勢(shì)以及突變?cè)\斷

為了解流域各站實(shí)測(cè)時(shí)間序列的變化，對(duì)原始序列的趨勢(shì)進(jìn)行MK趨勢(shì)檢驗(yàn)，同時(shí)為降低因單一突變?cè)\斷方法以及不同站點(diǎn)影響因素的不確定性造成的影響，分別采用Mann-Kendall，Pettitt以及有序聚類法3種突變檢驗(yàn)方法，在無TFPW前處理與有TFPW前處理的條件下對(duì)流域內(nèi)各站降水以及徑流序列進(jìn)行檢驗(yàn)，以最終確定其突變點(diǎn)，同時(shí)利用Morlet小波對(duì)序列的周期性進(jìn)行分析，以達(dá)到對(duì)流域內(nèi)水文序列變化的系統(tǒng)性分析。原始序列MK趨勢(shì)檢驗(yàn)結(jié)果如表1所示，未經(jīng)TFPW處理后的突變檢驗(yàn)結(jié)果如表2所示，降水和徑流序列在3種突變檢驗(yàn)方法下顯示出的突變檢驗(yàn)結(jié)果存在一定的差異，且存在多突變點(diǎn)的情況?？紤]到不同的突變檢驗(yàn)方法的各自局限性，這些突變點(diǎn)可能是序列的自相關(guān)性“加強(qiáng)”的偽突變點(diǎn)或者是真的存在多個(gè)突變點(diǎn)的情況，需要進(jìn)一步借助于能去除序列自相關(guān)性的TFPW方法進(jìn)行分析。經(jīng)TFPW方法前處理后的序列，在3種突變檢驗(yàn)方法下的突變檢驗(yàn)結(jié)果如表3所示，結(jié)果表明經(jīng)TFPW去除序列的自相關(guān)性以后，各序列的突變點(diǎn)識(shí)別表現(xiàn)得較為一致。最終，關(guān)河水庫流域降水突變年份確定為1978年，兩站點(diǎn)實(shí)測(cè)流量的突變年份均為1979年。

表1 原始序列MK趨勢(shì)檢驗(yàn)結(jié)果Tab.1 The result of MK trend test

表2 原始突變方法檢驗(yàn)結(jié)果Tab.2 The result of original mutation method test

表3 經(jīng)TFPW預(yù)處理后突變點(diǎn)檢驗(yàn)結(jié)果Tab.3 The result of mutation method test with TFPW

2.3 降水徑流的周期分析

小波分析是水文時(shí)間序列分析的一項(xiàng)重要工具，其主要應(yīng)用于水文序列的變化趨勢(shì)以及周期分析等[19]。對(duì)于多時(shí)間尺度分析，小波可將序列進(jìn)行分解研究，并提取出水文序列的變化周期，進(jìn)而反映出水文時(shí)間序列的變化趨勢(shì)。其中Morlet小波由于在時(shí)間與頻率的局部化之間有著較好的平衡，因此本文選用Morlet小波對(duì)流域內(nèi)兩站點(diǎn)流量以及流域面雨量進(jìn)行多尺度的時(shí)間分析，其結(jié)果如圖2，圖3所示。

圖2 關(guān)河流域兩徑流站點(diǎn)Morlet小波變換實(shí)部時(shí)頻分析圖Fig.2 The result of Morlet wavelet in two runoff station of Guanhe basin

圖3 關(guān)河流域降水量Morlet小波變換實(shí)部時(shí)頻分析圖Fig.3 The result of Morlet wavelet in precipitation station of Guanhe basin

小波變換實(shí)部時(shí)頻圖可將流域內(nèi)降水以及徑流序列的豐枯以及突變點(diǎn)變化情況反映出來。由圖2(a)可以看出，關(guān)河站流量時(shí)間序列在1979年突變點(diǎn)之前存在著幾處明顯的變化周期：一是以4～6 a為特征時(shí)間尺度的變化，其平均變化周期為5 a，在1966-1979年經(jīng)歷了4個(gè)增減變化；二是以7～9 a為特征時(shí)間尺度的變化，平均變化周期為8 a，經(jīng)歷了2個(gè)增減變化；三是以22～24 a為特征時(shí)間尺度的變化，平均變化周期為23 a，經(jīng)歷了1個(gè)增減變化。但在1979年突變點(diǎn)之后，該站徑流時(shí)間序列的周期特性發(fā)生了改變：在經(jīng)歷了5 a左右的不穩(wěn)定時(shí)期后，突變點(diǎn)前3種特征時(shí)間尺度變化逐漸消失，僅存在一種特征時(shí)間尺度，即以11～14 a為特征時(shí)間尺度的變化，平均周期為12 a，且從實(shí)部時(shí)頻分析圖可以看出，其在一個(gè)周期內(nèi)的變化幅度有了明顯的下降。從2012年以后未封閉的等值線圖來看，關(guān)河水庫站的流量呈現(xiàn)出一定程度的上升趨勢(shì)，但需要未來新的時(shí)間序列進(jìn)行驗(yàn)證?？偟膩砜?，關(guān)河水庫站的流量隨著時(shí)間的變化，其增減變化程度下降，且周期性現(xiàn)象產(chǎn)生了明顯的變化。

與關(guān)河水庫站不同，石棧道站的流量時(shí)間序列分析結(jié)果呈現(xiàn)出另一種特點(diǎn)。突變點(diǎn)之前的徑流變化周期主要有兩種：一是以4～6 a為特征時(shí)間尺度的變化，其在突變點(diǎn)前經(jīng)歷了5個(gè)周期的增減變化；二是以17～19 a為特征時(shí)間尺度的變化，其在突變點(diǎn)之前經(jīng)歷了2個(gè)周期的增減變化。1979年之后，以4～6 a為特征時(shí)間尺度的變化仍然存在，但其震蕩幅度逐漸降低，主要變化周期變?yōu)橐?2～14 a為特征時(shí)間尺度的變化?？傮w來看，突變點(diǎn)之后的周期變化有所改變，震蕩幅度有一定程度的減小，但其變化程度不如關(guān)河水庫站明顯。

相對(duì)流量變化，降水的變化則更具有規(guī)律性。降水序列的震蕩變化主要以4～6，9～12 a以及22～27 a為主要震蕩周期，平均周期分別為5，10以及25 a。2012年之后的檢驗(yàn)結(jié)果同關(guān)河水庫以及石棧道的結(jié)果相同，均呈現(xiàn)出上升趨勢(shì)。相較于徑流變化周期的不規(guī)則性，降水序列的變化幅度較小，且變化周期較為穩(wěn)定。

從Morlet小波圖像結(jié)果來看，徑流序列的主要周期變化點(diǎn)與突變點(diǎn)檢驗(yàn)的結(jié)果基本接近，徑流的周期性檢驗(yàn)結(jié)果與突變點(diǎn)檢驗(yàn)的結(jié)果相符。

3 水文序列變化的特點(diǎn)及驅(qū)動(dòng)機(jī)制分析

通過對(duì)流域各雨量站以及流域平均雨量的分析可知，流域降水量時(shí)間序列的突變點(diǎn)位于1978年，且流域內(nèi)兩個(gè)水文站點(diǎn)流量時(shí)間序列的突變點(diǎn)位于1979年，兩者突變發(fā)生的時(shí)間較為接近。同時(shí)對(duì)突變點(diǎn)之前的降水以及徑流的時(shí)間序列進(jìn)行Pearson相關(guān)性檢驗(yàn)，結(jié)果如表4所示。關(guān)河站流量與流域降雨量的相關(guān)系數(shù)為0.73，石棧道站流量與流域降雨量的相關(guān)系數(shù)為0.72，且從時(shí)間序列的圖像本身來看，三者均在1979年左右出現(xiàn)了顯著的拐點(diǎn)，因此可以認(rèn)定，流域降水量的變化是引起兩站點(diǎn)流量時(shí)間序列產(chǎn)生突變的原因。但突變點(diǎn)之后，兩站點(diǎn)的流量序列卻呈現(xiàn)出與降水時(shí)間序列不同的特點(diǎn)。從相關(guān)性分析結(jié)果來看，關(guān)河水庫站的流量時(shí)間序列與降水量的相關(guān)性下降至0.6，而石棧道站流量與流域降水量的相關(guān)性為0.71。

表4 Pearson相關(guān)性檢驗(yàn)結(jié)果Tab.4 The result of Pearson correlation test

從突變點(diǎn)前后的變化來看，突變點(diǎn)前多年平均降雨為567 mm，極值比為2.4，突變點(diǎn)后為490 mm，極值比為2.3。突變點(diǎn)前后的多年平均降水有一定程度的減少，但年際變化基本保持不變。作為受自然條件單一因素影響的石棧道站，突變點(diǎn)之前的平均流量為2.51 m3/s，流量極值比為5.97；突變點(diǎn)后平均流量為1.26 m3/s，平均流量下降幅度為49.8%。在氣候變化單一因素的影響下，該站流量呈現(xiàn)出均值減少，年際變化增加的特點(diǎn)。而作為受到氣候因素以及人類活動(dòng)影響雙驅(qū)動(dòng)下的關(guān)河站，其在突變點(diǎn)之前的平均流量為7.59 m3/s，極值比為7.82，突變點(diǎn)后平均流量為2.08 m3/s，極值比為7，平均流量下降幅度為72%?？梢钥闯觯谑艿綒夂蜃兓约叭祟惢顒?dòng)雙重因素的影響下，該站流量呈現(xiàn)出均值大幅減小，年際變化減弱的特點(diǎn)。

對(duì)比石棧道站以及關(guān)河站測(cè)得的實(shí)測(cè)流量可以發(fā)現(xiàn)，關(guān)河站流量下降的幅度更大，其序列變化的趨勢(shì)與特點(diǎn)也與石棧道站有著明顯的差別，因此，關(guān)河站流量突變點(diǎn)后的特性并不是氣候因素導(dǎo)致的。從水利工程的角度來看，關(guān)河水庫自建成以來，主要承擔(dān)了流域內(nèi)防洪以及灌溉的主要任務(wù)，水量維持較好，多年來水庫壩上水位始終維持在980 m左右。因此，相較于受水利工程影響較小的石棧道站，關(guān)河水庫站測(cè)得的水文序列更為穩(wěn)定，這也說明了關(guān)河水庫對(duì)徑流的調(diào)節(jié)作用影響較大。而石棧道站由于位于流域中段，且周圍并無大型的水利設(shè)施，因此其流量變化主要受到江水變化的影響，流量時(shí)間序列與降水量的相關(guān)性依然保持不變，其周期性變化與降水量的周期變化相對(duì)于關(guān)河水庫站也較為一致。

因此，從整體來看，降水量的變化是引起流域徑流產(chǎn)生突變點(diǎn)的原因，但受氣候變化單一因素以及氣候與人類活動(dòng)雙因素影響下的徑流出現(xiàn)了不同的特點(diǎn)?？偟膩碚f，受氣候變化單一影響下的石棧道站徑流呈現(xiàn)出“總量減小，峰值降低，年際變化增大，周期略有變化”的特點(diǎn)，而受到氣候與人類活動(dòng)雙因素影響下的關(guān)河站徑流整體呈現(xiàn)出“總量大幅下降，峰值明顯減少，年際變化減小，周期不穩(wěn)定且延長”的特點(diǎn)。

4 結(jié) 語

(1)相較于原始數(shù)據(jù)，經(jīng)過TFPW方法處理后的徑流數(shù)據(jù)在不同方法突變點(diǎn)識(shí)別的結(jié)果上更為統(tǒng)一，且檢驗(yàn)結(jié)果符合Morlet周期圖的分析結(jié)果，因此TFPW方法在徑流數(shù)據(jù)突變點(diǎn)的識(shí)別上具有一定的校正作用。

(2)根據(jù)流域內(nèi)多個(gè)測(cè)站的突變分析結(jié)果，關(guān)河水庫流域降雨突變點(diǎn)為1978年，徑流的突變點(diǎn)為1979年，且突變前兩者存在著較強(qiáng)的相關(guān)性，因此，降水量的變化是引起關(guān)河水庫流域兩測(cè)站徑流數(shù)據(jù)產(chǎn)生突變點(diǎn)的主要原因。

(3)相較于降水序列的變化，徑流序列在突變點(diǎn)前后的變化更為顯著，其在突變點(diǎn)前后的非一致性也更為明顯。

(4)徑流量突變點(diǎn)后兩站流量均呈現(xiàn)出明顯的下降趨勢(shì)，徑流變化的周期性也發(fā)生了改變，但關(guān)河站相較于石棧道站變化更為明顯。

(5)受氣候變化單一影響下的石棧道站主要體現(xiàn)在徑流量與極值變化上，其徑流呈現(xiàn)出“總量減小，峰值降低，年際變化增大，周期略有變化”的特點(diǎn)；而受到氣候與人類活動(dòng)雙因素影響下的關(guān)河站徑流主要體現(xiàn)在徑流總量，周期性與序列震蕩幅度上，整體呈現(xiàn)出“總量大幅下降，峰值明顯減少，年際變化減小，周期不穩(wěn)定且延長”的特點(diǎn)。

□