基于啟發(fā)式分割和近似熵法的徑流序列變異診斷*

黃生志，黃 強(qiáng)，王義民，陳昱潼

(西安理工大學(xué)西北旱區(qū)生態(tài)水利工程國家重點實驗室培育基地，陜西 西安 710048)

隨著全球氣候的變化以及日益加劇的人類活動的影響，流域水文要素以及水文情勢諸如河川徑流、湖泊蓄水以及流域降雨等無論是在時間上還是空間上均有可能發(fā)生變異，致使水文資料的一致性遭到破壞[1-2]。水文資料一致性的破壞將會給水文分析、水文模擬、水文計算以及流域的水資源規(guī)劃等帶來巨大的影響[3-5]，進(jìn)而導(dǎo)致變異前的成果不能再繼續(xù)使用，否則將產(chǎn)生極大的偏差。對流域水文變異進(jìn)行診斷，不僅有助于掌握流域水文要素的演變規(guī)律，還能避免因使用一致性被破壞的水文資料進(jìn)行設(shè)計或分析而產(chǎn)生重大失誤或者損失。

基于此，國內(nèi)外不少學(xué)者對水文變異進(jìn)行了研究并取得了大量的研究成果。諸如Lee等[6]運用了貝葉斯推斷理論對水文變點進(jìn)行了研究。然而，該方法對資料的要求過高且對樣本分布的假設(shè)過多，因此在實際的應(yīng)用過程中操作性不強(qiáng)[7]；陳廣才等[8]針對潮白河水資源分區(qū)的45 a年徑流量序列，采用基于平穩(wěn)與線性的檢驗假定的滑動F檢驗法、滑動T檢驗法以及秩和檢驗法等統(tǒng)計的方法，進(jìn)行了徑流序列的變異診斷。事實上，水文序列受到氣候變化以及人類活動的影響，往往呈現(xiàn)出高度的非線性和非平穩(wěn)性，使用傳統(tǒng)的統(tǒng)計診斷方法對水文序列的變異點進(jìn)行檢驗則往往產(chǎn)生一定的偏差[9]。

由Bernaola-Galván[10]提出的啟發(fā)式分割法主要用于處理非線性及非平穩(wěn)時間序列的均值變異點。同傳統(tǒng)的變異檢驗方法相比，啟發(fā)式分割法能夠?qū)⒁粋€非平穩(wěn)的時間序列分割成多個平穩(wěn)的子序列，各子序列的均值互不相同并表征不同的物理背景，能較好地彌補(bǔ)傳統(tǒng)檢驗方法的不足。由于該方法在分割時采取了一分為二的迭代算法，因此，極大地減少了計算量，具有較好的實用性[9]。此外，由Pincus于20世紀(jì)90年代提出的近似熵(approximate entropy)，因其具有較好的抗干擾與處理非穩(wěn)態(tài)序列能力且所需數(shù)據(jù)量較少等優(yōu)點，故在機(jī)械設(shè)備故障的診斷及醫(yī)學(xué)等領(lǐng)域得到廣泛地應(yīng)用[11]。鑒于此，本文基于啟發(fā)式分割算法與近似熵理論，對渭河流域的年徑流量進(jìn)行變異診斷，并對其變異進(jìn)行歸因分析。

1 研究區(qū)域與數(shù)據(jù)

1.1 研究區(qū)域概況

渭河作為黃河流域最大的一級支流，發(fā)源于甘肅省內(nèi)渭源縣的烏鼠山，自西向東分別流經(jīng)甘肅的隴西、天水，陜西的寶雞、西安、渭南，最終于潼關(guān)匯入黃河[12]。渭河的干流全長約為818 km, 總流域面積約為13.5萬km2(圖1)。渭河流域位于34°-38°N，104°-110°E，屬于大陸性季風(fēng)氣候，春季溫暖少雨，夏季雨熱同期且有伏旱，秋季涼爽濕潤，冬季寒冷少雨。流域的最冷月平均氣溫大約是-3～-1℃，最熱月平均氣溫大約為23～26℃[12]。流域降雨空間分布不均，由東南向西北遞減，秦嶺南麓雨水充沛，最大年降水量超過1 000 mm，而平原河套地區(qū)年降水量僅為500 mm左右[12]。

圖1 渭河流域及相關(guān)水文站位置Fig.1 The location of the Wei River Basin and the related hydrological stations

流域降水年內(nèi)分配不均且年際變化大，汛期降水量約為年降水量的65%。自從國家實行西部大開發(fā)這一國策以及成立了關(guān)中-天水國家級經(jīng)濟(jì)開發(fā)區(qū)后，渭河流域的經(jīng)濟(jì)社會高速發(fā)展，工業(yè)用水、農(nóng)業(yè)用水量劇增。因此，人類大量開采與利用地下水及地表水，致使渭河流域的河川徑流顯著減少。此外受全球氣候變化的影響，近年來渭河流域的降雨呈現(xiàn)下降的趨勢，如此一來渭河流域水文情勢將發(fā)生巨大的變化，日益減少的水資源將難以滿足渭河流域經(jīng)濟(jì)社會的可持續(xù)發(fā)展。鑒于渭河流域用水安全的重要性，為了進(jìn)一步掌握渭河流域河川徑流的演變規(guī)律，本文對該流域的年徑流量進(jìn)行變異診斷，并對變異結(jié)果進(jìn)行歸因分析。本文將研究區(qū)域劃分為三個部分，分別為林家村水文站以上流域，該區(qū)域代表渭河流域的中上游；張家山水文站以上流域，該區(qū)域代表渭河最大的支流涇河流域；華縣水文站以上流域，該區(qū)域代表整個渭河流域。

1.2 數(shù)據(jù)

本研究所用數(shù)據(jù)為渭河流域21個氣象站的1960-2005年的日氣象數(shù)據(jù)，通過求和運算得到各個氣象站的降雨、氣溫、風(fēng)速、日照變化率等年氣象序列以及林家村、張家山以及華縣水文站的1960-2005年的日徑流序列，通過求和而得到各水文站的年徑流量，渭河流域及各站點所在的位置如圖1所示。以上氣象數(shù)據(jù)均來源于中國國家氣象中心，水文數(shù)據(jù)摘自水文手冊。所有數(shù)據(jù)在使用之前均進(jìn)行了三性審查，審查結(jié)果表明本文所采用的水文氣象數(shù)據(jù)具有較好的可靠性、一致性以及代表性，可以用作相應(yīng)分析。潛在蒸發(fā)量是通過Penman-Monteith 公式計算而成。

2 研究方法

2.1 啟發(fā)式分割算法

某一時間序列x(t)由N個點組成，其中一個分割點i從序列的左邊沿著該序列依次向右邊滑動，分割點左邊和右邊部分的平均值分別為μ1(i)和μ2(i)，其標(biāo)準(zhǔn)差分別是s1(i)和s2(i)，i點的合并偏差SD(i)可表示為[9,13-14]：

(1)

式中，N1、N2分別表示i點左邊與右邊部分的點數(shù)。其中，i點左右兩邊兩子序列均值的差異能夠使用t檢驗統(tǒng)計值T(i)來進(jìn)行度量[9]：

(2)

T值越大，則說明該點左右兩邊兩子序列的差異越明顯。計算T(t)中最大值Tmax所對應(yīng)的統(tǒng)計顯著性P(Tmax)，其計算公式如下[9]：

(3)

通過蒙特卡洛模擬可得：η=4.19lnN-11.54且δ=0.40。其中，N表示該序列的長度，v=N-2,Ix(a,b) 是不完全β函數(shù)。預(yù)先設(shè)定一個臨界值P0，當(dāng)P(Tmax)≥P0，則在該點處將此序列分割成左右兩個均值差異較大的子序列，否則不進(jìn)行分割[9]。對得到的新序列不斷進(jìn)行迭代并重復(fù)以上操作，直到子序列的長度小于0(0是最小分割尺度)時便停止對其分割。由此，便將原序列分割成幾個均值不同的子序列，而分割點則為該序列的均值變異點。

通常情況下，P0取值范圍為0.5～0.95，而0的取值則不應(yīng)小于25[15]。通過不斷調(diào)整的大小，能夠?qū)崿F(xiàn)對該時間序列不同尺度上的變異檢測[9]。

2.2 近似熵方法

近似熵方法實質(zhì)上是一種基于邊緣概率分布統(tǒng)計的定量描述某一時間序列的復(fù)雜程度的方法，該算法的具體步驟簡要介紹如下[11]。

假定某一長度為N的時間序列u(1),u(2), …,u(N)，在該時間序列中構(gòu)造一組m維數(shù)的向量X(1),X(2), …,X(N+m-1),其中向量為[11]：

X(i)={u(i),u(i+1),…,u(i+m-1)}

i=1,2,...,N-m+1

(4)

定義兩個向量X(i)與X(j)之間的距離d[X(i),X(j)]是兩個向量中對應(yīng)的元素的最大值[11]：

d[X(i),X(j)]=max[|u(i+k)-u(j+k)|]

k=0,1,...,m-1

(5)

對于每一個i，定義兩個向量X(i)與X(j)之間的關(guān)聯(lián)程度：

(j)]≤r的個數(shù)}/

(N-m+1)

(6)

(7)

(8)

然而，實際上N不可能為無窮大。當(dāng)N取有限值時，便可得到近似熵的估計值，表示為：

ApEn(m,r)=φm(r)-φm+1(r)

(9)

3 結(jié)果與討論

3.1 基于啟發(fā)式分割算法的年徑流變異結(jié)果分析

根據(jù)上述2.1的具體步驟，將P0設(shè)定為0.95，0的取值則設(shè)為25，分別得到林家村以上流域、張家山以上流域以及整個渭河流域的年徑流變異分析結(jié)果，如圖2所示。

由圖2(a)中的第一次分割可知，1971年對應(yīng)的T值最大，并且其相應(yīng)的P(Tmax)等于1要大于臨界值0.95。因此，1971年是林家村以上流域年徑流序列的第一個突變點；從第二次分割可知，1994年所對應(yīng)的T值最大，且其相應(yīng)的P(Tmax)等于0.998 5，也大于臨界值0.95，故1994年為該流域年徑流序列的第二個突變點。

由圖2(b)可得，張家山以上流域的年徑流序列并無突變點，因其最大T值的P(Tmax)小于0.95。對于整個渭河流域的年徑流序列，由圖2(c)中第一次分割可知，1969年的T值為最大值，其對應(yīng)的P(Tmax)等于0.992 9，大于臨界值0.95，因此1969年為第一個變異點；由其第二次分割可得，1993年的T值為最大值，其對應(yīng)的P(Tmax)等于0.971 8，大于臨界值0.95，故1993年為第2個變異。

3.2 基于近似熵的年徑流變異分析結(jié)果

根據(jù)上述2.2求近似熵的具體步驟，采用滑動窗口的方法，將窗口寬度設(shè)為5，滑動步長設(shè)為1，m和r分別設(shè)為2和0.15σ。其中，σ為原始序列的標(biāo)準(zhǔn)差，分別得到林家村以上流域、張家山以上流域以及整個渭河流域的基于近似熵的年徑流變異分析結(jié)果于圖3(a)-(c)中。

這時零售商的和供應(yīng)商的期望收益均是供應(yīng)鏈整體期望收益函數(shù)的仿射函數(shù)，同時各企業(yè)的最優(yōu)定貨量q*相同，這表明在無突發(fā)事件下基準(zhǔn)契約能實現(xiàn)供應(yīng)鏈協(xié)調(diào).

由圖3(a)可知，林家村以上流域年徑流1960-1971年間的近似熵較大，而1971-1994年間的近似熵較小，1971年明顯是一個變異點，而1994-2005年間的近似熵較大，易知1994年為另一變異點。由圖3(b)可知，張家山以上流域年徑流序列的近似熵沒有顯著的變異情況。因此，該流域的年徑流序列無變異點。而由圖3(c)中整個渭河流域年徑流的近似熵明顯可知，渭河流域年徑流序列明顯存在兩個變異點，分別為1969與1993年。該結(jié)果與基于啟發(fā)式分割算法的結(jié)果一致，進(jìn)一步論證了以上3個流域年徑流序列變異點的正確性，也充分展示了基于啟發(fā)式分割算法及近似熵理論在水文變異診斷上的優(yōu)越性。

3.3 流域年徑流變異點的印證

3.3.1 累積年徑流與年份關(guān)系圖 雖然通過啟發(fā)式分割算法及近似熵方法已經(jīng)發(fā)現(xiàn)林家村以上流域年徑流序列存在2個變異點(1971和1994年)，張家山以上流域年徑流序列無變異點存在，整個渭河流域年徑流序列也存在兩個變異點(1969和1993年)。

圖2 年徑流啟發(fā)式分割檢驗圖Fig.2 The test chart of annual runoff based on heuristic segmentation algorithm

為了進(jìn)一步證明這兩個算法的在變異診斷上的準(zhǔn)確性，繪制3個流域的累積年徑流量與年份的關(guān)系圖，其結(jié)果如圖4所示。

由圖4實心圓部分可以明顯看出林家村以上流域累積年徑流被2個變異點分成了3部分，且每部分的線性擬合度都很高。其中，1960-1971與1972-1994年2個部分的線性擬合相關(guān)系數(shù)高達(dá)0.99，而1995-2005的線性擬合系數(shù)為0.96。而這兩個變異點恰好與上述3.2節(jié)中基于啟發(fā)式分割法及近似熵理論所找到的關(guān)于林家村以上流域年徑流序列的變異點相一致。圖4空心五角星部分所表示的是張家山以上流域累積年徑流量與年份之間的關(guān)系，在整個研究階段1960-2005年間，所有的點均擬合于一條直線，且擬合精度高，其相關(guān)系數(shù)高達(dá)0.99。因此，在1960-2005年間，張家山以上流域年徑流序列無明顯變異點存在。由圖4的空心圓部分可知，整個渭河流域的累積年徑流量被1969與1993年這2個變異點分成了3個部分，每個部分點的線性擬合程度都很高，3個部分的相關(guān)系數(shù)均高達(dá)0.99。因此，通過累積年徑流量與年份的關(guān)系圖可知，整個渭河流域的年徑流序列存在1969和1993年這2個變異點，而這又與上述利用啟發(fā)式分割法及近似熵理論所發(fā)現(xiàn)的整個渭河流域年徑流量的變異的十分吻合。由此進(jìn)一步論證了基于啟發(fā)式分割法及近似熵理論在水文變異診斷中的準(zhǔn)確性；與此同時，進(jìn)一步印證了渭河流域年徑流序列存在變異的事實。其中，林家村以上流域與整個渭河流域的年徑流序列的確存在兩個變異點，而張家山以上流域年徑流序列并無變異點存在。此外，對于林家村以上流域和整個渭河流域被兩個變異點分割成的3個部分而言，其擬合直線的傾斜度均隨著時間推移而不斷減小，由此可知這兩個流域的徑流量有著明顯的下降趨勢，而實際上也確實如此。該研究結(jié)果與拜有存[16]的基于差異信息理論的渭河徑流序列變異診斷結(jié)果基本一致，他發(fā)現(xiàn)渭河流域的年徑流存在這兩個變異點。

3.3.2 變異前后各參數(shù)對比 分別計算存在變異點的林家村以上流域以及整個渭河流域變異前后徑流序列的均值、方差以及變差系數(shù)等參數(shù)，進(jìn)一步對比分析變異點前后序列的差異性，具體的參數(shù)值如表1所示。由表1可知：① 林家村以上流域由兩個變異點所分割成的3部分徑流序列無論是均值還是方差均有較大差別，就均值而言，1972-1994年間的年徑流均值較1960-1971年減少了42.2%，而年徑流均值1995-2005年比1960-1972年減少了77.3%，可以說3個部分的均值差別比較明顯。② 而就方差而言，1972-1994年間的年徑流序列的方差較1960-1971年減少了71.4%，而1995-2005年的年徑流序列的方差同1960-1972年相比減少了86.9%，由此可知，這3部分的方差差異性較大。但是，對于變差系數(shù)而言，1972-1994年間的變差系數(shù)值較1960-1972年變化不大，而1995-2005的年徑流序列的變差系數(shù)與1960-1972年相比減少了40.9%。③ 整個渭河流域的年徑流序列變異前后參數(shù)的變化規(guī)律與林家村以上流域類似，3個部分的均值與方差的差異性較大，1970-1994年間的變差系數(shù)值同1960-1970相比變化不大，而1994-2005的年徑流序列的變差系數(shù)較1960-1969年減少了25%。

圖3 年徑流的近似熵檢驗圖Fig.3 The test chart of annual runoff based on approximate entropy method

圖4 三個流域累積年徑流量與年份的關(guān)系圖Fig.4 The relationship between accumulative runoff and year in the three basin

總的來說，被兩變異點分割成的3個年徑流序列在均值與方差上差異明顯，并且均值、方差及變差系數(shù)有明顯的減小趨勢。由此可側(cè)面證明變異點的正確性。

表1 兩流域年徑流序列變異前后參數(shù)對比表Table 1 The comparison of parameters of the annual runoff before and after the change points in the two basin

3.4 徑流變異的歸因分析

3.4.1 氣候變化的影響 采用Mann-Kendall法求算變點前后降雨和潛在蒸發(fā)的變化趨勢，其結(jié)果如表2所示。

由表2可知，林家村以上流域與渭河流域的降雨和潛在蒸發(fā)變點前后的變化趨勢相似。降雨先減小后變大，但在整個時間序列表現(xiàn)為不顯著地下降趨勢。潛在蒸發(fā)量先增大后減小再增大，在整個時間序列表現(xiàn)為不顯著的上升趨勢。從常規(guī)分析上可知渭河流域與徑流聯(lián)系密切的降雨與潛在蒸發(fā)等氣象因子均未發(fā)生顯著變化，人類活動對徑流變異的影響或許更強(qiáng)。然而，渭河流域地處中國大陸性季風(fēng)區(qū)，受西太平洋副熱帶高壓系統(tǒng)的影響。1971年時，西太平洋副熱帶高壓系統(tǒng)有所減弱并發(fā)生一個El Nio事件致使中國大部分地區(qū)的大氣環(huán)流及氣候發(fā)生變異，從而導(dǎo)致渭河流域當(dāng)年降雨量劇減[17]，而此時恰好是林家村以上流域年徑流序列的第一個變異點所對應(yīng)的時間。而1994年時，另外一個ENSO事件發(fā)生，此時西太平洋副熱帶高壓系統(tǒng)異常強(qiáng)大，導(dǎo)致渭河流域當(dāng)年的降雨量較年平均降雨量減少了14.9%[18]，而此時恰好對應(yīng)于林家村以上流域及整個渭河流域年徑流序列的第2個變異點。

3.4.2 人類活動的影響 20世紀(jì)70年代初，渭河流域修建了不少水利工程，如1970年修建的羊毛灣水庫和馮村水庫，1972年修建的大峪水庫和石門水庫。這些水庫的修建，極大地改變河川徑流的時空分布規(guī)律，破壞天然的徑流規(guī)律，從而導(dǎo)致徑流序列變異。除了興建水利工程外，流域下墊面的變化也改變河川徑流的演變規(guī)律。50年代末60年代初，由于經(jīng)濟(jì)困難與政治運動的緣故，大量的樹木遭到破壞，尤其在1958年的“大煉鋼鐵”時期，集中砍伐了大量的林木，從而造成了嚴(yán)重的水土流失，由此也使得1960-1970年間的徑流量較大[19]。但是，隨著水土保持措施的開展，尤其是70-80年代，渭河流域?qū)λ亮魇ч_展了一系列有組織、有計劃、有步驟的防治工作，并取得了很好的成效，從而在一定程度上減小了渭河流域的水土流失量。因此，在1973-1994年間河川徑流較1960-1970年有一定的減少。此外，隨著人口的迅速增加以及經(jīng)濟(jì)的高速發(fā)展，渭河流域的工業(yè)用水和農(nóng)業(yè)用水增長迅速。其中，1990-2000年間，國民經(jīng)濟(jì)各行業(yè)用水增加了17.1億m3，年均增長率為3.2%，遠(yuǎn)高于全國平均1.1%和黃淮海流域的0.3%。90年代后的渭河流域國民經(jīng)濟(jì)平均用水量為42.62億m3，比90年代前增加了52.6%。大量地開采與開發(fā)地表水和地下水，是渭河流域徑流序列變異的又一重要因素。

表2 變點前后降雨和潛在蒸發(fā)的變化趨勢Table 2 The trend of precipitation and potential evaporation before and after change points

4 結(jié) 論

為了探究渭河流域河川徑流變異情況，本文采用1960-2005年林家村、張家山以及華縣3個水文站的日徑流數(shù)據(jù)，將研究區(qū)域劃分為林家村以上流域、張家山以上流域以及整個渭河流域3個部分?；趩l(fā)式分割算法與近似熵方法對各流域的年徑流序列進(jìn)行了變異診斷，并利用累積年徑流量與年份的關(guān)系圖對變異結(jié)果進(jìn)行判定分析，得到主要結(jié)論如下：

1)啟發(fā)式分割算法與近似熵方法均具有較好的檢驗均值突變點的能力，具有較強(qiáng)的抗噪抗干擾的能力，十分適應(yīng)于非線性非穩(wěn)態(tài)的徑流序列的變異診斷，能夠準(zhǔn)確地找到其變異點。

2)林家村以上流域存在2個變異點，分別是1971和1994年；張家山以上流域無變異點存在；整個渭河流域存在兩個變異點，分別是1969和1993年。

3)對于林家村以上流域和整個渭河流域而言，被兩個變異點分割成3個子序列的線性趨勢線的斜率隨著年份的增大而不斷減小，由此說明在這兩個流域，其徑流量有明顯的減小趨勢。

4)造成渭河流域河川徑流變異的主要原因是氣候變化與日益增長的人類活動。其中，氣候變化主要是1971年的El Nio事件和1994年的ENSO事件的發(fā)生致使當(dāng)時渭河流域降雨量劇減，從而引起徑流量的顯著減少。人類活動則具體包括水利工程建設(shè)、流域下墊面的改變以及由于人口和經(jīng)濟(jì)的增長而大量開采與開發(fā)地下水和地表水，由此而使得河川徑流發(fā)生變異。

[1] 鄧建偉，宋松柏，盧書超. 石羊河流域年徑流序列的變異診斷[J]. 西北農(nóng)林科技大學(xué)學(xué)報:自然科學(xué)版，2006，34(4) ：121-124.

[2] 李艷，陳曉宏，張鵬飛. 北江流域徑流序列年內(nèi)分配特征及其趨勢分析[J]. 中山大學(xué)學(xué)報:自然科學(xué)版， 2007， 46(5) ：113-116.

[3] 丁晶，劉權(quán)授. 隨機(jī)水文學(xué)[M]. 北京：中國水利水電出版社，1997.

[4] 王孝禮，胡寶清，夏軍. 水文序列趨勢與變異點的R/S分析法[J]. 武漢大學(xué)學(xué)報：工學(xué)版，2002，35(2)：10-12.

[5] 詹道江，葉守澤. 工程水文學(xué)[M]. 北京：中國水利水電出版社，2000.

[6] AUSTIN F S L, SYLVA M H. A shift of the mean level in a sequence of independent normal random variable-A Bayesian Approach [J]. Technometrics, 1997, 19(4): 503-506.

[7] 肖宜，夏軍，申明亮，等. 差異信息理論在水文時間序列變異點診斷中的應(yīng)用[J]. 中國農(nóng)村水利水電，2001(11)：28-30.

[8] 陳廣才，謝平. 水文變異的滑動F識別與檢驗方法[J]. 水文，2006，26(2)：57-60.

[9] 陳廣才，謝平. 基于啟發(fā)式分割算法的水文變異分析研究[J]. 中山大學(xué)學(xué)報:自然科學(xué)版，2008，47(5)：122-125.

[10] BERNAOLA-GALVN P, IVANOV P C, AMARAL L A N,et al. Scale invariance in the nonstationarity of human heart rate[J].Physical Review Letters, 2001, 87:168105.

[11] 王啟光，張增平. 近似熵檢測氣候突變的研究[J]. 物理學(xué)報，2008，57(3)：1976-1982.

[12] 張宏利，陳豫，任廣鑫，等. 近50年來渭河流域降水變化特征分析[ J].干旱地區(qū)農(nóng)業(yè)研究，2008，26(4)：236- 241.

[13] 李海彬，張小峰，胡春宏,等. 基于B-G分割算法的河川年輸沙量突變分析[J]. 水利學(xué)報，2010，41(12)：1387-1391.

[14] 封國林，龔志強(qiáng)，董文杰. 基于啟發(fā)式分割算法的氣候突變檢測研究[J]. 物理學(xué)報，2005，54(11)：5494-5499.

[15] 龔志強(qiáng). 基于非線性時間序列分析方法的氣候突變檢測研究[D]. 揚(yáng)州：揚(yáng)州大學(xué)，2006.

[16] 拜有存. 基于差異信息理論的渭河年徑流變異點診斷[J]. 人民黃河，2009，31(7)：29-30.

[17] BAI H，MU X M，WANG F，et al. Analysis on evolution law of meteorological and hydrological drought and wetting[J]. Agricultural Research in the Arid Areas，2012，30:237-241.

[18] ZHANG H L，CHEN Y， REN G X. The characteristics of precipitation variation of Weihe River Basin in Shaanxi Province during recent 50 years[J]. Agricultural Research in the Arid Areas，2008，26:236-240.

[19] 張淑蘭，王彥輝，于澎濤,等.定量區(qū)分人類活動和降水量變化對涇河上游徑流變化的影響[J]. 水土保持學(xué)報，2010，24(4)：53-58.