基于線性矩法的區(qū)域枯水頻率分析

基于線性矩法的區(qū)域枯水頻率分析

1概率權(quán)重矩的線性組合線性矩陣來自綠化巖等提出的概率權(quán)重矩陣（vm）。對于隨機變量x，r級概率權(quán)重矩陣設(shè)置如下。βr=∫10x(F)FrdF,r=1,2,3??(1)式中F(x)為隨機變量X的分布函數(shù);x為概率為F時的反函數(shù)值.隨機變量X的第r+1階線性矩定義為λr+1=∫10x(F)Ρ*r(F)dF,r=0?1??(2)式中Ρ*r(F)=r∑k=0(-1)r-k(r+k)!(k!)2(r-k)!Fk,利用該多項式和式(1),式(2)可化為λr+1=r∑k=0(-1)r-k(r+k)!(k!)2(r-k)!βk,r=0,1,?(3)從式(3)可以看出,隨機變量X的線性矩是概率權(quán)重矩的線性組合.一般在實際應(yīng)用中,需利用有限的樣本數(shù)來計算樣本線性矩.假設(shè)隨機變量X有n個樣本,按從小到大排序,即x1∶n≤x2∶n≤…≤xn∶n,xi∶n(i=1,2,…,n)為在n個值中第i小的值.與隨機變量X的線性矩類似,前r+1階樣本線性矩定義為lr+1=r∑k=0(-1)r-k(r+k)!(k!)2(r-k)!bk,r=0,1,\:,n-1(4)式中bk為概率權(quán)重矩βk的無偏估計值,bk=n-1n∑i=k+1(i-1)(i-2)\:(i-k)(n-1)(n-2)\:(n-k)xi∶n,k=0,1,\:,n-1。樣本L-CV定義為t=l2/l1(5)樣本線性矩系數(shù)定義為tr=lr/l2,r=3,4,5,??n(6)2區(qū)域頻率分析2.1枯水徑流資料及長度東江流域發(fā)源于江西省尋鄔縣的椏髻體,地處亞熱帶,氣候溫和,多年平均降雨量為1750mm,年徑流量較為豐沛,但年內(nèi)變化較大.東江流域為珠江流域的主要水系,地處經(jīng)濟發(fā)達的華南地區(qū),該流域除了向廣東省惠州、深圳等市供水外,還承擔著向香港特區(qū)供水的重任.近些年來華南地區(qū)枯季干旱缺水情況屢屢發(fā)生,枯水期遭受了嚴重的供水危機,再加上沿海地區(qū)咸潮的頂托影響,使得東江流域枯水徑流的研究具有很重要的現(xiàn)實意義.本文選取其中16個子流域作為研究對象,流域面積大約為37.2—2091km2,且基本不受人類活動的影響.這些子流域大多具有較長的徑流資料,日平均流量資料的長度為9—45年,各子流域的平均資料長度約為22年,具體情況見表1.從表1可以發(fā)現(xiàn),在這些流量站中,有5個站點的資料長度較長,大多接近或超過40年,一般可以滿足枯水頻率分析的需要;其余11個站點的資料長度較短,直接由這些站點來推求設(shè)計枯水流量,其精度一般難以得到保證.本文將16個站點組成統(tǒng)一的研究區(qū)域,采用區(qū)域頻率分析方法來探討該研究區(qū)域的頻率分布,利用資料長度較長站點的水文信息來提高資料長度較短站點的設(shè)計枯水流量的預(yù)測精度,進而達到“以空間換取時間”的目的.在枯水徑流的頻率分析中,對于日平均流量資料而言,可采用1、3、7、10、15、30、60、90等不同的時段平均值作為研究對象.由于最小七日平均流量可以更好地反映枯水徑流的特性,大多以其作為枯水徑流的統(tǒng)計特征量.本文結(jié)合東江流域的具體特點,利用滑動平均法得到每個水文年(每年的4月—次年的3月)內(nèi)的七日平均流量序列,以每個水文年的最小值組成年最小七日平均流量序列,對這些序列進行枯水徑流的頻率分析研究.2.2區(qū)域頻率分析區(qū)域頻率分析一般包括以下4個主要步驟:各個站點的一致性檢驗、水文相似區(qū)檢驗、區(qū)域分布函數(shù)的選擇和區(qū)域頻率計算等.(1)各站點的枯水流量核算,計算每軸點的密度值的計算如果某個區(qū)域有N個站點,令ui=[t(i)t(i)3t(i)4]T為第i個站點的線性矩系數(shù)矩陣.區(qū)域線性矩系數(shù)的算術(shù)平均值為ˉu,定義為ˉu=1ΝΝ∑i=1ui,則第i個站點的一致性檢驗系數(shù)定義為Di=13Ν(ui-uˉ)ΤA-1(ui-uˉ)(7)式中:A=∑i=1Ν(ui-uˉ)(ui-uˉ)Τ若計算得到的Di大于某個標準值,則該站點就與其它站點不一致.Hosking等人研究發(fā)現(xiàn),判定某個站點是否與其它站點不一致的標準并不是固定不變的,而是與區(qū)域站點數(shù)有關(guān).計算得到的研究區(qū)域16站點的樣本矩系數(shù)t、t3、t4以及Di值見表1.從中可看出,水背和下陂兩站的Di值皆大于16個站點時的判斷標準值3,故將這兩個站點從研究區(qū)域中剔除,將剩余14個站點作為一個研究區(qū)域進行分析.另外,從表1也可看出,紅花塔、黎嶂峰、水背和桃西4個站點的t3(相當于偏態(tài)系數(shù))和t4(相當于峰度系數(shù))有時為負值,可能是由于這些站點的枯水流量序列本身原因所造成,也有可能是序列長度較短所致.(2)水文相似區(qū)檢驗假定研究區(qū)域有N個站點,第i個站點的樣本長度為ni,則以樣本長度為權(quán)重的區(qū)域平均線性矩系數(shù)分別為tR=∑i=1Νnit(i)/∑i=1Νni,t3R=∑i=1Νnit3(i)/∑i=1Νni,t4R=∑i=1Νnit4(i)/∑i=1Νni(8)以樣本長度為權(quán)重的樣本線性矩系數(shù)t的均方差為V1={∑i=1Νni[t(i)-tR]2/∑i=1Νni}1/2(9)為了判別由這些站點所組成的區(qū)域是否為水文相似區(qū),采用MonteCarlo統(tǒng)計試驗方法得到模擬值V序列.以穩(wěn)健性較好的4參數(shù)Kappa分布作為總體分布,并令區(qū)域平均線性矩系數(shù)為其分布函數(shù)的線性矩.假設(shè)模擬區(qū)也有N個站點,每個站點都符合Kappa分布,且與實際站點具有同樣的樣本長度,由式(8)和(9)計算所生成的樣本模擬值V.重復(fù)模擬Nsim次得到Nsim個模擬值V,根據(jù)Nsim個樣本計算V的均值μV和均方差σV.則水文相似區(qū)檢驗系數(shù)H為Η=(V1-μV)σV(10)如果H值較大,則該區(qū)域為非水文相似區(qū).Hosking等人認為,如果H<1,則可以將該區(qū)域看成為水文相似區(qū);如果1≤H≤2,則可能為非相似區(qū);如果H≥2,則一般為非相似區(qū).一般情況下,Nsim等于500次就可以得到較好的精度.計算得到的區(qū)域線性矩為:LR=1,tR=0.2374?t3R=0.1178?t4R=0.1438用于模擬樣本V的Kappa分布參數(shù)為:ξ=0.8447,α=0.3343,k=0.0374,h=-0.1468.實際區(qū)域V1值為0.0432,500次模擬后得到的均值μV和均方差σV分別為0.0363和0.0075,則水文相似區(qū)檢驗系數(shù)H為0.93.由于得到的水文相似區(qū)檢驗系數(shù)H值小于1,可將由該14個子流域所組成的區(qū)域看作為水文相似區(qū),采用統(tǒng)一的分布函數(shù)進行區(qū)域頻率分析.(3)區(qū)域平均4階線性矩系數(shù)鑒別出水文相似區(qū)后,需要確定采用何種分布線型作為該水文相似區(qū)的統(tǒng)一分布.一般采用區(qū)域分布函數(shù)的檢驗系數(shù)Z加以判別,其定義為ΖDΙSΤ=(τ4DΙSΤ-t4R+B4)/σ4(11)其中B4=Νsim-1∑m=1Νsim(t4[m]-t4R),σ4=[(Νsim-1)-1{∑m=1Νsim(t4[m]-t4R)2-ΝsimB42}]1/2式中DIST為某種分布函數(shù);τ4DΙSΤ為由DIST分布模擬得到的區(qū)域平均值;t4R為由實際區(qū)域得到的區(qū)域平均值;假定檢驗區(qū)域分布線型為4參數(shù)的Kappa分布,令區(qū)域平均線性矩系數(shù)為其分布函數(shù)的線性矩,利用MonteCarlo試驗方法生成第m個樣本,t[m]4為其樣本的區(qū)域平均4階線性矩系數(shù);B4和σ4分別為從Kappa分布模擬得到的t4R的偏差和均方差.一般而言,如果ZDIST≤1.64,則可以采用該分布作為統(tǒng)一的區(qū)域分布函數(shù),如果滿足條件的分布函數(shù)不止一個,通常取ZDIST最小的分布函數(shù)作為最優(yōu)的分布.除此之外,很多學者也利用線性矩系數(shù)相關(guān)圖來確定分布函數(shù),通過比較區(qū)域樣本線性矩系數(shù)t3R和t4R分布函數(shù)的線性矩系數(shù)τ3和τ4的相關(guān)程度來判斷分布函數(shù)的類型,如果某個分布函數(shù)的線性矩系數(shù)和與區(qū)域線性矩系數(shù)相一致,則該分布函數(shù)即為最優(yōu)的分布函數(shù).本文考慮的分布函數(shù)包括5種常見的三參數(shù)分布函數(shù):GeneralizedPareto分布(GPA)、GeneralizedExtreme-Value分布(GEV)、GeneralizedLogistic分布(GLO)、LogNormal分布(LN3)和PearsontypeⅢ(PⅢ).5種分布函數(shù)的ZDIST計算結(jié)果見表2,依據(jù)判定準則,GEV、GLO、LN3和PⅢ皆可以作為東江流域枯水的區(qū)域頻率分布,但LN3分布最優(yōu).5種分布函數(shù)的τ3和τ4關(guān)系見圖1,區(qū)域平均線性矩系數(shù)分別為t3R=0.1178,tR4=0.1438,它們和GEV、LN3和PⅢ皆比較接近,但最接近于LN3分布,該結(jié)果與表2中的ZDIST值相吻合.另外,Pandey等人于2001年利用4階線性矩系數(shù)τ4來確定序列長度較短樣本的分布函數(shù)的線型,其判定原則是如果某個分布函數(shù)的4階線性矩系數(shù)和樣本的4階線性矩系數(shù)之差的絕對值最小,那么該分布函數(shù)即為最優(yōu)的分布函數(shù),該方法已經(jīng)被Pandey等人從理論和實際應(yīng)用中加以論證.表3也給出東江流域的計算結(jié)果,從中看出,LN3分布為最優(yōu)的分布函數(shù).結(jié)合分布函數(shù)的檢驗系數(shù)ZDIST、線性矩系數(shù)相關(guān)圖和4階線性矩系數(shù)τ4判定準則,選擇LN3分布作為該區(qū)域的統(tǒng)一分布線型.(4)預(yù)測值的計算確定出區(qū)域分布函數(shù)為LN3后,以LN3分布作為統(tǒng)一的分布線型計算各站在不同重現(xiàn)期下的枯水流量.利用指標洪水法推求的各站預(yù)測值為Qi(F)=Qˉiq(F)(12)式中Qi(F)為第i個站點頻率為F時的預(yù)測值;Qˉi為第i個站點的尺度量,此處為多年平均的年最小七日平均流量,即為表1中的一階線性矩l1;q(F)為頻率為F時的區(qū)域最小七日平均流量值.利用區(qū)域線性矩系數(shù)1、tR和t3R,得到LN3分布的參數(shù)為ξ=0.9496,α=0.4106,k=-0.2419,求得在各種重現(xiàn)期T下的區(qū)域枯水流量值q(F).將該區(qū)域化的計算結(jié)果乘上某站點的平均流量Qˉi就得到該站點某個頻率下的最小七日平均流量.2.3gorten公式將LN3分布作為區(qū)域化的分布函數(shù),得到各種重現(xiàn)期下的枯水流量,比較各個站點的實測最小七日平均流量與區(qū)域頻率模擬得到的最小七日平均流量.各個站點的經(jīng)驗頻率采用Gringorten公式來模擬Ρ(i)=i-0.44n+0.12(13)式中i為從小到大排列的序號;n為樣本總數(shù).利用式(13)計算出各個站點實測枯水流量序列所對應(yīng)的經(jīng)驗頻率,用所得到的經(jīng)驗頻率基于區(qū)域LN3分布求得各個站點的模擬值.比較同頻率下的實測值和模擬值,結(jié)果表明,所有站點的兩者之間的平均相對誤差為9.39%.圖2為研究區(qū)域內(nèi)14個站點的實測和模擬結(jié)果的散點圖,除了少數(shù)點據(jù)以外,實測和區(qū)域模擬結(jié)果基本吻合.因此,采用LN3分布來分析東江流域的最小七日平均流量的分布可以得到較高的精度.3研究區(qū)域的劃分線性矩法具有較好的穩(wěn)健性,為一類優(yōu)越的分布參數(shù)確定方法,現(xiàn)已在洪水和枯水頻率分析中得以廣泛應(yīng)用,本文利用線性矩法進行東江流域最小七日