水文序列跳躍變異點的滑動相關(guān)系數(shù)識別方法

吳子怡，謝 平，2，桑燕芳，雷 旭，袁樹堂，王 超

（1.武漢大學 水資源與水電工程科學國家重點實驗室，湖北 武漢 430072；

2.國家領(lǐng)土主權(quán)與海洋權(quán)益協(xié)同創(chuàng)新中心，湖北 武漢 430072；3.中國科學院 地理科學與資源研究所陸地水循環(huán)與地表過程重點實驗室，北京 100101；4.云南省水文水資源局，云南 昆明 650106）

1 研究背景

近幾十年來，氣候變化和人類活動對流域下墊面條件改變所引起的水文效應［1］，已成為全球水循環(huán)變化研究的焦點問題。變化環(huán)境下，認為水文序列不再來自同一總體，即失去了原有的一致性［2］。對于非一致性的水文序列，若直接采用傳統(tǒng)方法進行分析與計算，如水文頻率計算，其水文設計結(jié)果將存在較大誤差，必將增加工程水文設計的風險［3］，同時會對水資源評價、規(guī)劃、管理等一系列水文水資源問題的研究產(chǎn)生影響［4］。

處理非一致性水文序列時，首要問題是判斷其是否存在變異點。變異點主要是指在某個時間位置序列前后變化特性出現(xiàn)顯著的差異，存在轉(zhuǎn)折［5］、跳躍［4］等多種表現(xiàn)形式。其中，跳躍變異點反映出水文序列均值發(fā)生突變抬升或下降的現(xiàn)象，打破了序列變化的連續(xù)性，是非一致性水文序列研究中最關(guān)心的問題之一。國內(nèi)外許多學者圍繞跳躍變異點（突變點）的識別問題開展了研究。Kundzewicz等［6］對跳躍變異點識別方法進行了總結(jié)和討論，并指出了各方法的適用條件與優(yōu)缺點。其中，Mann-Kendall秩次相關(guān)法［7-8］與Pettitt檢驗法［9-10］是兩種應用較為廣泛的非參數(shù)方法。Bayesian變點分析法也被用于水文氣象時間序列的跳躍變異點識別［11］。基于Bayesian理論，Lee和Heghinian提出Lee-Heghinian法［12］，主要解決均值跳躍變異點識別問題。Brown-Forsythe檢驗［13］可減少跳躍變異點分析對數(shù)據(jù)過多假設的要求［14］，也被用于水文序列的跳躍變異點識別。Li等［15］提出了用可變模糊集方法對灤河流域降雨、徑流序列進行跳躍變異點識別。此外，還有大量的相關(guān)研究成果與理論方法［16］。然而，從目前研究可以看出，跳躍點識別方法多是從統(tǒng)計學角度進行描述，未能很好地直觀反映時間序列跳躍變異的現(xiàn)象和本質(zhì)，且未能合理考慮跳躍變異點前后均值跳躍幅度、離散程度等因素的影響。因此，實際中需要從序列產(chǎn)生跳躍變異的實質(zhì)出發(fā)，以尋求更加合理有效、易于理解和進行物理成因分析的時間序列跳躍變異點識別方法。

為此，本文借助非一致性序列確定性成分與隨機性成分的概念［17］，并以易于求解的相關(guān)系數(shù)為基礎，通過求解跳躍成分與原始序列的相關(guān)系數(shù)作為跳躍變異點識別指標，提出一種約束條件較少的跳躍變異點識別方法，通過設計統(tǒng)計實驗對該方法在不同因素影響下的效率進行討論，并與常用的Pettitt檢驗法、Brown-Forsythe檢驗法進行對比。最后將該方法應用于瀾滄江流域允景洪水文站不同時間尺度的徑流序列進行跳躍變異點識別與檢驗，并結(jié)合成因分析對所提方法進一步驗證。

2 原理與方法

對于任一長度為n的水文序列假設其存在某一跳躍變異點，則整個序列可由變異點分割為長度n1和n2的兩段子序列，其中n2段序列存在均值跳躍成分b：

原始序列X與跳躍成分序列Y具有一定的相關(guān)性，其相關(guān)系數(shù)r計算式為：

由式（3）可知，對于某一水文序列，相關(guān)系數(shù)可能為正或負，當其為正時，原始序列X與跳躍成分序列Y成正相關(guān)，表明序列具有向上的跳躍變異；反之，相關(guān)系數(shù)為負時，表明序列具有向下的跳躍變異。當序列長度和統(tǒng)計特性固定時，其相關(guān)系數(shù)的絕對值|r|與跳躍變異點前后子序列均值之差呈正相關(guān)，相關(guān)系數(shù)絕對值最大時所對應的點最有可能為跳躍變異點。此外，式（3）還表明，在跳躍變異點的滑動過程中，相關(guān)系數(shù)還受到序列長度（n），跳躍變異點位置（n1，n2）以及離散程度（σx）的影響。當離散程度較大時，相關(guān)系數(shù)值較小，而序列長度與變異點位置對相關(guān)系數(shù)值的影響不能直觀從式（3）確定，需結(jié)合統(tǒng)計實驗等方式加以明確。

基于上述分析，提出利用相關(guān)系數(shù)對水文序列跳躍變異點進行滑動識別與顯著性檢驗的方法，本文稱為“跳躍變異點滑動相關(guān)系數(shù)識別方法”，具體思路如下：（1）設定合理的起訖點，并逐步向序列尾端滑動，利用式（3）依次求取各點對應的跳躍成分序列Y與原始序列X的相關(guān)系數(shù)；（2）取所有結(jié)果中絕對值最大的相關(guān)系數(shù)|rmax|，對其進行假設檢驗，其中相關(guān)系數(shù)臨界值可由服從F分布的F統(tǒng)計量臨界值間接確定［18］；（3）根據(jù)所需選擇合適的顯著水平α，若|rmax|＜rα，即未通過假設檢驗，認為該序列無跳躍變異點；若|rmax|＞rα，即通過假設檢驗，認為|rmax|所對應的點位為該水文序列的跳躍變異點。

實際中，考慮到滑動點位過于靠近首端或尾端時會使得子序列樣本容量太小，導致抽樣誤差太大和結(jié)果不可靠，因此利用上述方法進行跳躍變異點診斷時通常不從序列起點開始。本方法中取滑動起始點為10，終點為n-10。

3 不同因素對序列跳躍變異點識別的影響分析

設計多組統(tǒng)計實驗，檢驗所提方法對序列跳躍變異點識別的性能，并分別選取應用較為廣泛的非參數(shù)檢驗法——Pettitt檢驗法［19］與同樣具有代表性的參數(shù)檢驗法——Brown-Forsythe檢驗法［14］進行對比研究。考慮到序列長度、離散程度等因素對跳躍變異點識別均有影響，因此在設計統(tǒng)計實驗時分別對各主要因素影響下三種方法的性能進行討論。生成的模擬序列均假設服從P-III分布，考慮的影響因素包括：序列長度n，跳躍變異點位置系數(shù)a，跳躍變異點前后均值變化系數(shù)Kave，Cv變化系數(shù)KCv，Cs變化系數(shù)KCs。跳躍變異點前（后）子序列長度為n1（n2），則表示變異前序列長度占原序列長度的比值，即變異點在序列中出現(xiàn)的位置，取值為（0，1）；跳躍變異點前（后）子序列均值為則離散系數(shù)為Cvs（Cvs′），則偏態(tài)系數(shù)為Css（Css′），則

由于統(tǒng)計實驗存在抽樣誤差，實際中需要考慮相應的誤差允許度δ允許（此處取δ允許=0.01）。設各次統(tǒng)計實驗結(jié)果為對應的已知跳躍變異點位置為各次統(tǒng)計實驗序列長度為nk={ }nk|k=1,2,3,…,m，則誤差表達式如下：

若δk≤δ允許，認為此次檢驗有效。將每種方法的有效檢出次數(shù)與該組實驗總次數(shù)之比記為該方法在此次實驗中的識別效率η（0≤η≤100%）。

3.1 n對跳躍變異點識別的影響設跳躍變異點位于序列中部（a=0.5），且取變異點前序列均值sˉ=500，Cvs=0.4，Css=0.8，取Kave=2，KCv=1，KCs=1?？紤]到水文分析與計算時通常要求樣本長度大于20，因此取20為序列長度n的下限。令n值分別在20～200，201～400，401～600三組中隨機產(chǎn)生，每組生成100個滿足上述所有條件的序列，并分別用3種方法進行跳躍變異點檢驗，結(jié)果見表1。

對比表1中不同方法的結(jié)果可以看出，隨著序列長度增長，各方法的檢驗效率均有明顯提高?；瑒覴識別法在序列長度不同的3組統(tǒng)計實驗中，對跳躍變異點的檢驗效率與Pettitt檢驗法、B-F檢驗法基本相同，部分略有偏高。特別是當n＞200時，無論序列長度n如何延長，3種方法的識別效率都已達到90%以上，可認為此情況下，n的變化將不再對檢驗結(jié)果產(chǎn)生影響。但實際水文分析與計算時，大多數(shù)實測水文序列的樣本長度在100年以內(nèi)，因此n=20～200時的識別效率結(jié)果顯示，提出的滑動R識別法與其他兩種常用方法性能相近，均適用于對實際水文序列的跳躍變異點進行識別。

表1 不同序列長度n下3種方法識別效率

3.2 a對跳躍變異點識別的影響取序列長度n=100，變異點前序列均值sˉ=500，Cvs=0.4，Css=0.8，取Kave=2，KCv=1，KCs=1。選取序列9個變異點位置，即令a分別取序列長度的0.1、0.2、0.3、0.4、0.5、0.6、0.7、0.8、0.9，每組生成100個滿足上述所有條件的序列，并分別得到3種方法對應的識別效率（圖1（a））。

從圖1（a）可以看出，當跳躍變異點位置處于序列較中間位置（0.3＜a＜0.7）時，3種方法的識別效率基本相同，均在70%左右。變異點處于較前段（a＜0.3）時，Pettitt方法的識別效率最低，滑動R識別法的效率明顯高于Pettitt方法，且較為穩(wěn)定。變異點處于后段（a＞0.7）時，滑動R識別法的效率升至最高，且高于其他兩種方法。整體來看，隨著a的變化，滑動R識別法的效率較為穩(wěn)定且保持在70%～80%左右，并對位于序列尾部的變異點有更好的識別效果，基本不受變異點位置變化的影響。因此，認為對于不同位置的跳躍變異點，提出的滑動R識別法均有較好的識別效率。

圖1 3種方法識別效率對比

3.3 Kave對跳躍變異點識別的影響設序列長度n=100、a=1/2，取跳躍變異點前序列均值sˉ=500，Cvs=0.4，Css=0.8，取KCv=1，KCs=1。令Kave分別取1.2、1.4、…、3.0，即-s′=600、700、…、1500，每組生成100個滿足上述所有條件的序列，并分別得到3種方法對應的識別效率（圖1（b））。從圖1（b）可以看出，隨著Kave的增大，Pettitt檢驗法、B-F檢驗法與滑動R識別法的識別效率均顯著增加，且變化幅度區(qū)別不大。該結(jié)果表明，當序列變異點的跳躍程度越明顯時，各方法對跳躍變異點識別的性能均更優(yōu)，跳躍變異點更容易被識別出。

3.4 KCv對跳躍變異點識別的影響設序列長度n=100，跳躍變異點位于序列中部，取跳躍變異點前序列均值pˉ=500、Cvs=0.4、Css=0.8，取Kave=2、KCs=1。令KCv分別取0.5、1.0、…、3.0，即Cvs′=0.2、0.4、…、1.2，每組生成100個滿足上述所有條件的序列，并分別得到3種方法對應的識別效率（圖1（c））。從圖1（c）可以看出，隨著KCv的增大，序列的離散程度增大，Pettitt檢驗法、B-F檢驗法與滑動R識別法的識別效率均有明顯下降，且下降趨勢相近。表明當序列的離散程度增大時，均值跳躍程度會減弱，導致無法更準確地對其進行識別與分離。

3.5 KCs對跳躍變異點識別的影響設定序列長度n=100，跳躍變異點位于序列中部，取跳躍變異點前序列均值sˉ=500、Cvs=0.4、Css=0.8，取Kave=2、KCv=1。令KCs分別取0.5、1.0、1.5、2.0，即Css=0.4、0.8、1.2、1.6，每組生成100個滿足上述所有條件的序列，并分別得到3種方法對應的識別效率（圖1（d））。從圖1（d）可以看出，隨著KCs的增大，Pettitt檢驗法、B-F檢驗法與滑動R識別法的效率雖略有起伏，但波動不大，基本維持在80%左右，說明該因素對序列變異點識別的影響不大。

結(jié)合上述統(tǒng)計實驗結(jié)果可知，在各因素影響下，滑動相關(guān)系數(shù)識別法對跳躍變異點的識別效率與常用的Pettitt檢驗法、B-F檢驗法基本相同，甚至在有些因素影響下識別效率還略高于另兩種方法，且其性能相對穩(wěn)定。因此，認為本文提出的滑動相關(guān)系數(shù)識別法可以有效地檢驗出水文序列的跳躍變異點。此外，Pettitt檢驗法與B-F檢驗法均未能如滑動相關(guān)系數(shù)識別方法直觀地反映出水文序列跳躍變異點前后均值跳躍變異幅度、離散程度等的影響。相比較而言，滑動相關(guān)系數(shù)識別法不僅易于求取，也更便于對跳躍變異點識別結(jié)果進行直觀判斷，還可由相關(guān)系數(shù)定量分析水文序列中各個參數(shù)對跳躍變異識別的影響，進而揭示某一跳躍變異點的出現(xiàn)主要是由哪個水文序列參數(shù)的變異所引起的，可為水文變異歸因分析提供可靠的依據(jù)。

4 實例分析

4.1 研究區(qū)概況瀾滄江發(fā)源于青海省南部的唐古拉山脈，流經(jīng)青、藏、滇三省（區(qū)），于云南省西雙版納州勐臘縣出境，出境后稱湄公河。允景洪水文站位于西雙版納州景洪市，是瀾滄江下游主要控制站（圖2），為國家重要水文站（一類精度水文站）。

瀾滄江流域水能資源豐富。在云南省境內(nèi)，允景洪水文站上游的瀾滄江干流水利工程建設從1986年開始［20］，至2014年建成投產(chǎn)的大型水利工程有6座。將時間順序劃分為多組開發(fā)運行時期（圖3）。其中，漫灣水電站具有季調(diào)節(jié)能力，小灣水電站與糯扎渡水電站是多年調(diào)節(jié)水庫。這些水利工程的建設對瀾滄江流域的徑流情勢造成了較大影響［21］，如大壩截流蓄水對干流年際徑流分配的影響，雨季水庫蓄水改變洪量帶來的徑流年內(nèi)分配變化［22］，特別是徑流年內(nèi)短期變化受干流梯級水電站影響明顯［21］。

圖2 允景洪水文站位置及流域示意

圖3 云南省境內(nèi)瀾滄江干流大型水電站開發(fā)運行時間軸

4.2 結(jié)果與分析選取允景洪站1956—2014年不同時間尺度的徑流序列，采用本文提出的滑動R識別法對各序列進行跳躍變異點檢驗，同時仍與Pettitt檢驗法、B-F檢驗法的結(jié)果進行對比。若所提方法與其余兩種方法所得變異點一致，則信任該結(jié)果；若3種方法所得跳躍變異點存在不一致，則選擇識別結(jié)果相同的兩種方法所給出的跳躍變異點，并結(jié)合各結(jié)果對應的跳躍成分變化圖進行合理性判斷。最終得到的各序列跳躍變異點的綜合診斷結(jié)果見表2。

表2 允景洪站多時間尺度徑流序列跳躍變異點識別檢驗結(jié)果

由表2結(jié)果可統(tǒng)計出每種方法的有效檢出率，即該方法識別結(jié)果與綜合判斷一致的結(jié)果數(shù)占全部識別結(jié)果數(shù)的比值。可以得到滑動R識別法對變異點的有效檢出率為86.7%，Pettitt檢驗法為66.7% ，B-F檢驗法為86.7%。其中，Pettitt檢驗法在a、e、g、j、k等多個序列中均存在識別誤差，而B-F檢驗法在h、i序列中存在誤差，滑動R識別法在c、d序列中存在誤差?？梢?，滑動R識別法在實測水文序列跳躍變異點識別與檢驗中有很好的效果。再給出綜合結(jié)果所對應的各序列模比系數(shù)（即原始數(shù)據(jù)與序列均值之比）跳躍成分圖（圖4），直觀反映出跳躍變異點識別的準確性。

為了進一步說明本文所提方法識別跳躍變異點的可靠性，從降水與水電站修建運行兩方面著手，分析跳躍變異點出現(xiàn)的物理成因。取允景洪水文站周邊5個氣象站1957—2014年的年-月尺度降水序列及區(qū)域面均降水序列，與處理徑流序列相同，也采用上述3種方法對各個降水序列進行跳躍變異點檢驗，結(jié)果見表3。

圖4 各序列模比系數(shù)跳躍成分圖

表3 允景洪站周邊氣象站多時間尺度降雨序列跳躍變異點識別檢驗綜合結(jié)果

分析表3可知，除景洪站外，大部分站點及所處區(qū)域降水均未出現(xiàn)跳躍變異點。不少學者也對此進行了研究分析［23-24］，認為降水過程變異不是該流域徑流變異的主要物理成因。

結(jié)合圖4與圖3對比分析，2004年允景洪站上游糯扎渡水電站與漫灣水電站二期工程開工建設，小灣進行大江截流，其勢必對水電站下游的年徑流量帶來擾動，盡管大朝山水電站已于2003年開始運行，但其水庫不具有多年調(diào)節(jié)性能，因此年徑流量仍主要受施工影響，在2004年出現(xiàn)向下的跳躍變異；年最大1、3、5月徑流也于2004年出現(xiàn)減少，最小月徑流卻未減少甚至有所增大（最小5月），說明汛期與枯期的徑流量主要受大朝山水電站年內(nèi)調(diào)節(jié)作用的影響，且削峰作用更明顯。年最大與最小1、3、5、7日徑流可反映汛期與枯期的極值情況。允景洪站年最大1、3、5、7日徑流序列在2002年左右出現(xiàn)明顯減少，說明2002年后，年最大n日洪量被明顯削弱，除可能受到小灣水電工程開工的影響外，還反映出漫灣水電站的汛期調(diào)蓄效果明顯。年最小1、3日徑流也在1992年左右出現(xiàn)小幅度跳躍式減小，這應與1992年大朝山水電工程開工以來，一系列水電工程的依次建設有關(guān)。

此外，圖4還直觀地反映出一個普遍現(xiàn)象：在各時間尺度徑流序列的尾端（2010年之后），非常明顯地顯示出汛期徑流量減少、枯期徑流量增加的特征。這充分說明瀾滄江干流多級水電站投入運行后，起到了顯著且良好的“削豐補枯”作用。但由于這些變化出現(xiàn)在序列末尾，超出變異點檢驗方法的應用范圍，因此未在檢驗結(jié)果中有所體現(xiàn)。這也說明，無論跳躍變異點識別與檢驗方法是否有效，都需結(jié)合水文序列時序圖與物理成因進行分析，才能獲得更準確且全面的水文過程變異信息。

5 結(jié)論

本文主要提出了水文序列跳躍變異點的滑動R識別方法，在數(shù)學公式嚴格推導基礎上，設計一系列統(tǒng)計實驗，對影響跳躍變異點識別的多種因素分別進行分析，并與常用的Pettitt檢驗法、B-F檢驗法進行識別效率的比較。運用本文所提方法，對瀾滄江流域內(nèi)允景洪水文站的不同時間尺度徑流序列進行跳躍變異點識別及檢驗。結(jié)合降水序列變異檢測與大型水電工程施工運行情況，對徑流過程的突變現(xiàn)象進行成因分析?？偨Y(jié)得出以下結(jié)論：（1）本文所提滑動R識別法可描述和直觀反映水文序列跳躍變異點前后均值跳躍變異幅度、離散程度等因素的影響；采用該方法時，相關(guān)系數(shù)容易求取，在受到序列長度、變異點位置、均值變化系數(shù)、Cv變化系數(shù)與Cs變化系數(shù)等多種因素影響時，其跳躍變異點識別效率與常用的Pettitt檢驗法、B-F檢驗法基本相同，甚至略高，且性能較為穩(wěn)定。（2）利用提出的滑動相關(guān)系數(shù)識別法對瀾滄江下游允景洪水文站多時間尺度徑流序列進行跳躍變異點識別與檢驗，同時綜合Pettitt檢驗法及B-F檢驗法的識別結(jié)果，以及各序列模比系數(shù)跳躍成分變化圖得到各序列的跳躍變異點結(jié)果。在這些綜合檢測結(jié)果中，滑動R識別法的效率達到86.7%，與B-F檢驗法相近，但遠高于Pettitt檢驗法的識別效率。（3）結(jié)合物理成因分析，認為允景洪站徑流變異受降水影響很小，主要與上游一系列梯級水電站的施工與運行有關(guān)，且水電站的調(diào)蓄作用，特別是年內(nèi)“削豐補枯”對徑流有明顯的影響；同時也印證了本文所提方法識別跳躍變異點的可靠性。

本文提出的滑動R識別法可有效識別并檢驗出水文、氣候等時間序列的跳躍變異點，其主要是側(cè)重于對均值變異的識別與檢驗，是否能有效識別出離散系數(shù)Cv和偏態(tài)系數(shù)Cs的跳躍變異點還需后續(xù)作進一步的研究。此外，對于存在多個跳躍變異點的水文序列，可在一次識別后對相應的跳躍性成分進行剔除，再進行二次識別，以此類推，逐一識別出序列的多個跳躍變異點。但其準確率及受序列方差等因素影響的問題還需進一步進行驗證。同時還需注意到，不同跳躍變異點識別方法都可能存在一定的誤差，因此需要進行物理成因分析，同時將多種方法通過賦權(quán)等方式結(jié)合使用，盡可能減少和避免單一方法造成的識別誤差。

［1］宋曉猛，張建云，占車生，等.氣候變化和人類活動對水文循環(huán)影響研究進展［J］.水利學報，2013，44（7）：779-790.

［2］宋松柏，李揚，蔡明科.具有跳躍變異的非一致分布水文序列頻率計算方法［J］.水利學報，2012，43（6）：734-739.

［3］梁忠民，胡義明，王軍.非一致性水文頻率分析的研究進展［J］.水科學進展，2011，22（6）：864-871.

［4］謝平等.變化環(huán)境下地表水資源評價方法［M］.北京：科學出版社，2009.

［5］陳遠中，陸寶宏，張育德，等.改進的有序聚類分析法提取時間序列轉(zhuǎn)折點［J］.水文，2011，31（1）：41-44.

［6］KUNDZEWICZ Z W，ROBSON A J.Change detection in hydrological records—a review of the methodology［J］.Hydrological Sciences Journal，2004，49（1）：7-19.

［7］MANN H B.Nonparametric tests against trend［J］.Econometrica，1945，13（3）：245-259.

［8］BISAI D，CHATTERJEE S，KHAN A，et al.Application of sequential mann-kendall test for detection of approximate significant change point in surface air temperature for Kolkata Weather Observatory，West Bengal，India［J］.International Journal of Current Research，2014，6（2）：5319-5324.

［9］PETTITT A N.A non-parametric approach to the Change-Point problem［J］.Journal of the Royal Statistical Society，1979，28（2）：126-135.

［10］MA Z，KANG S，ZHANG L，et al.Analysis of impacts of climate variability and human activity on streamflow for a river basin in arid region of northwest China［J］.Journal of Hydrology，2008，352（3/4）：239-249.

［11］PERREAULT L，BERNIER J，BOBéE B，et al.Bayesian change-point analysis in hydrometeorological time series.Part 1.The normal model revisited［J］.Journal of Hydrology，2000，235（3/4）：221-241.

［12］HEGHINIAN S M.A shift of the mean level in a sequence of independent normal random variables：a bayesian approach［J］.Technometrics，1977，19（4）：503-506.

［13］BROWN M B，F(xiàn)ORSYTHE A B.Robust tests for the equality of variances［J］.Journal of the American Statistical Association，1974，69（346）：364-367.

［14］ZHANG Y C，ZHOU C H，LI B L.Brown-Forsythe based method for detecting change points in hydrological time series［J］.Geographical Research，2005，24（5）：741-748.

［15］LI J，TAN S，WEI Z，et al.A New method of change point detection using variable fuzzy sets under environmental change［J］.Water Resources Management，2014，28（14）：5125-5138.

［16］DEAN W E，ANDERSON R Y.Application of some correlation coefficient techniques to time-series analysis［J］.Mathematical Geosciences，1974，6（4）：363-372.

［17］丁晶，劉權(quán)授.隨機水文學［M］.北京：中國水利水電出版社，1997.

［18］邰淑彩.應用數(shù)理統(tǒng)計［M］.武漢：武漢大學出版社，2005.

［19］HE T，LU Y，CUI Y，et al.Detecting gradual and abrupt changes in water quality time series in response to regional payment programs for watershed services in an agricultural area［J］.Journal of Hydrology，2015，525：457-471.

［20］王億春，胡云萍.淺析瀾滄江上游水利工程對允景洪水文站水文要素的影響［C］//云南省水利學會2014年度學術(shù)交流會論文集.2014.

［21］HE D M，F(xiàn)ENG Y，GAN S，et al.Transboundary hydrological effects of hydropower dam construction on the Lancang River［J］.Chinese Science Bulletin，2006，51（B11）：16-24.

［22］鐘華平，王建生.瀾滄江干流水電開發(fā)對徑流的影響分析［J］.水利水電技術(shù)，2010，41（12）：72-74.

［23］李斌，李麗娟，李海濱，等.1960-2005年瀾滄江流域極端降水變化特征［J］.地理科學進展，2011，30（3）：290-298.

［24］尤衛(wèi)紅，何大明，索渺清.瀾滄江的跨境徑流量變化及其對云南降水量場變化的響應［J］.自然資源學報，2005（3）：361-369.