梁冀雨，林炳章

(南京信息工程大學應用水文氣象研究院，南京 210044)

1 研究背景

相較于傳統(tǒng)的單站頻率分析，地區(qū)頻率分析法通過整合運用某一特定地區(qū)內(nèi)所有站點歷史資料序列的總體平均特征，優(yōu)化對該區(qū)域內(nèi)各站點的頻率分析結(jié)果，使之更為準確、穩(wěn)健[1]。在使用地區(qū)頻率分析法時，需滿足假設(shè)：研究區(qū)域內(nèi)所有站點除了有一個不同的尺度系數(shù)外，頻率分布的線型和參數(shù)應保持一定程度的一致[2]。滿足上述假設(shè)時，稱該研究區(qū)域為一致區(qū)。一致區(qū)的劃分是地區(qū)頻率分析的基礎(chǔ)，區(qū)域的一致性將直接影響頻率估計值的準確性[3]。因此,準確有效的區(qū)域一致性判別方法具有重要的實踐意義和理論價值。

在地區(qū)頻率分析法的早期應用中，Dalrymple[4]于1960年提出了一種基于樣本離差系數(shù)(Cv)或偏態(tài)系數(shù)(Cs)的區(qū)域一致性檢驗方法，被大量應用于早期的降雨以及徑流的地區(qū)頻率分析中。但是該法對經(jīng)驗閾值的依賴較大，導致在應用于某些特定類型的洪水資料序列時，對異質(zhì)站點的分辨能力較差。為了解決這一問題，Wiltshire[5]在1984年提出兩種改進方案，其一仍然基于樣本離差系數(shù)，其二基于區(qū)域的頻率分布曲線。前者因與分布無關(guān)，適用范圍更廣。之后，Lu等[6]于1992年利用歸一化的概化極值分布(GEV)的方差特性構(gòu)建了新的一致性檢驗方法。在這些一致性檢驗方法的基礎(chǔ)上，Hosking和Wallis[7]在1993年提出了至今應用最為廣泛的Hosking-Wallis一致性檢驗(后文簡稱HW檢驗)。HW一致性檢驗以樣本線性矩離差系數(shù)(L-Cv)和線性矩偏態(tài)系數(shù)(L-Cs)代替?zhèn)鹘y(tǒng)的離差系數(shù)和偏態(tài)系數(shù)，分析研究區(qū)域內(nèi)各站點資料序列線性矩系數(shù)的離散程度，能夠很好地與地區(qū)線性矩法(RLMA)相結(jié)合，對異質(zhì)性的鑒別能力較強。

雖然HW一致性檢驗至今仍被大量應用，但是卻存在理論上的缺陷。Hosking和Wallis在先前的研究中發(fā)現(xiàn)，站點資料間的相關(guān)性會降低地區(qū)頻率分析中一致性檢驗的效力[3],其后2人提出的HW檢驗也未針對該缺陷進行調(diào)整。Castellarin等[8]發(fā)現(xiàn)HW檢驗的效力也會隨著資料相關(guān)性的增大而減小，即當資料存在相關(guān)性時，HW檢驗可能將存在異質(zhì)站點的區(qū)域判定為一致。Castellarin等給出的解決方案是一個由站點平均相關(guān)系數(shù)確定的經(jīng)驗修正系數(shù)，但該經(jīng)驗修正系數(shù)在每次使用時需事先根據(jù)資料律定參數(shù)，無法直接應用。

本文根據(jù)HW檢驗的理論基礎(chǔ)，重點分析資料相關(guān)性影響其效力的原因，并使用蒙特卡羅模擬法揭示資料相關(guān)性對HW檢驗的影響。同時利用江西省年最大降雨量序列(AMS)研究實際資料的相關(guān)性特征，通過結(jié)合資料和理論，給出一種不受資料相關(guān)性影響的一致性檢驗。最后再次通過蒙特卡羅模擬對新的一致性檢驗方法進行評估。結(jié)果顯示，該法可以為地區(qū)頻率分析法，特別是地區(qū)線性矩法提供更有效的理論支持。

2 Hosking-Wallis一致性檢驗

為了將線性矩[9]應用于地區(qū)頻率分析，Hosking和Wallis在1993年提出了基于資料線性矩系數(shù)的HW檢驗。HW檢驗將研究區(qū)域內(nèi)各站點資料的各階線性矩系數(shù)的離散程度，與符合一致區(qū)定義的、與研究區(qū)域具有相同站點數(shù)以及資料長度的人工模擬資料相比較。2者的差別若在一定的范圍之內(nèi)，就認定研究區(qū)域為一致區(qū)。

(1)

于是,區(qū)域內(nèi)各站點的平均加權(quán)標準差V為：

(2)

最后進行蒙特卡羅模擬，生成500組模擬區(qū)域數(shù)據(jù)，根據(jù)式(2)計算得到各組模擬數(shù)據(jù)的平均加權(quán)標準差Vsim。利用Vsim的均值μV和標準差δV，對研究區(qū)域的平均加權(quán)標準差V進行類正態(tài)化處理，得到HW檢驗統(tǒng)計量H：

(3)

為了得到μV、σV，需進行蒙特卡洛模擬。首先根據(jù)式(1)，使用研究區(qū)域各站點資料計算得到的前4階地區(qū)平均線性矩系數(shù)，擬合4參數(shù)的kappa分布[10]。接著模擬生成500個人工區(qū)域，每個人工區(qū)域都擁有和研究區(qū)域相同的站點數(shù)、資料長度，同時符合先前得到的kappa分布。值得注意的是，由此得到的每個人工區(qū)域都符合一致區(qū)的假設(shè)，并且站點間不存在相關(guān)關(guān)系。最后，對每一個人工區(qū)域分別計算平均加權(quán)標準差Vsim，并得到500個Vsim的均值和標準差。

Hosking和Wallis認為[9]，當H小于1時，研究區(qū)域為可接受的一致區(qū)；當H大于1并小于2時，研究區(qū)為可能的異質(zhì)區(qū)域；而當H大于2時，研究區(qū)域為確定的異質(zhì)區(qū)域。

另外，在考慮線性偏態(tài)系數(shù)和線性峰度系數(shù)的情況下，V具有2種變式，即：

(4)

式中：V2、V3分別代表同時考慮線性離差系數(shù)、線性偏態(tài)系數(shù)，以及線性偏態(tài)系數(shù)、線性峰度系數(shù)的情況。

此時，為區(qū)別3者不同，式(2)中的V改稱V1。而由3者得到的H對應記作H1、H2和H3。事實上，3者結(jié)構(gòu)相似，計算得到的H值往往代表相同的結(jié)果，因此Hosking等[9]建議多數(shù)情況下可僅考慮H1的結(jié)果。其后，Viglione等[11]發(fā)現(xiàn)，H1相較于H2以及H3，具有更強的鑒別異質(zhì)性(heterogeneity)能力，且表現(xiàn)更為穩(wěn)定。綜上，本文僅對H1進行分析，相應的統(tǒng)計量，后文簡稱H以及V。

3 實測與生成雨量資料的相關(guān)性

相關(guān)系數(shù)是表征2個資料序列間相關(guān)程度的統(tǒng)計量。其中最為常用的是統(tǒng)計學家K. Pearson于1895年首次提出的皮爾森積矩相關(guān)系數(shù)，簡稱相關(guān)系數(shù)。對于2個成對的隨機變量序列X、Y，其皮爾森相關(guān)系數(shù)的表達式為：

(5)

式中：rX,Y、cov(X,Y)分別表示隨機序列X、Y的皮爾森相關(guān)系數(shù)和協(xié)方差。

由式(5)可知，皮爾森相關(guān)系數(shù)可以被看作2個隨機變量之間標準化的協(xié)方差，當隨機序列X、Y為標準正態(tài)分布時，2者的標準差都為1，式(5)右邊可簡化為cov(X,Y)。因此，選擇皮爾森相關(guān)系數(shù)作為相關(guān)性的度量能夠在進行后文的蒙特卡羅模擬實驗時，通過分解協(xié)方差矩陣的方法[12]，生成具有特定相關(guān)性的資料序列。

為了比較HW檢驗中生成的模擬年最大降雨量資料與實測資料相關(guān)性的區(qū)別，本文以與HW檢驗相同的設(shè)置[9]進行蒙特卡羅模擬實驗。從15 a到50 a，每一個序列長度生成1 000個人造區(qū)域，每個人造區(qū)域都包含100個站點。為做到與實測資料一致，令每個站點的模擬資料都服從4參數(shù)的kappa分布[10]，并擁有江西全省平均的24 h年最大雨量資料線性矩系數(shù)：t=0.18，t3=0.20，t4=0.15。同時因?qū)Y料進行歸一化處理，所以所有生成資料的均值為1。之后，求出每個人工區(qū)域的站點間平均相關(guān)系數(shù)，最后計算每個序列長度1 000個人工區(qū)域平均相關(guān)系數(shù)的均值，繪于圖1。

圖1 實測與模擬雨量序列平均相關(guān)系數(shù)隨序列長度的分布

圖1中，圓形標記為實測年最大雨量平均相關(guān)系數(shù)隨資料長度的分布，三角標記為HW檢驗中使用的模擬資料的情況。根據(jù)圖1，模擬相較實測資料折線圖下降趨勢更為明顯、平滑，且始終位于實測資料的下方，表明實測資料中普遍存在的相關(guān)性大于HW檢驗中生成的與之比較的模擬數(shù)據(jù)，且相差達到25%以上。

4 資料相關(guān)性對Hosking-Wallis一致性檢驗結(jié)果的影響

4.1 蒙特卡羅模擬中年最大雨量相關(guān)性對H統(tǒng)計量的影響

為了定量研究年最大雨量資料的相關(guān)性對H統(tǒng)計量的影

圖2 資料相關(guān)性對Hosking-Wallis一致性檢驗的影響

由圖2可知，相關(guān)性對不同的地區(qū)分布的影響程度是大致相同的，資料間的平均相關(guān)系數(shù)越大，H統(tǒng)計量曲線的位置越靠下，其值越小。當區(qū)域平均相關(guān)系數(shù)從0.2增加到0.8時，H值減小量均達到1.5以上，已經(jīng)超過了判定有效一致區(qū)的H值區(qū)間[0,1]的長度，這大大降低了HW檢驗對資料相關(guān)區(qū)域異質(zhì)性的辨別能力。

4.2 不同站點間相關(guān)性空間分布情況對H統(tǒng)計量的影響

(5)

式中：R為資料的相關(guān)系數(shù)矩陣，其形式?jīng)Q定了生成的資料間相關(guān)系數(shù)的均一的。

(6)

圖3展示了地區(qū)分布為GEV，區(qū)域站點數(shù)為66時，不同相關(guān)系數(shù)分布情況對H值的影響。圖3的橫縱坐標與圖2相同，改進后的模擬結(jié)果以圓點虛線分別繪于2張圖上。圖3(a)中黑色實線和三角形點虛線分別表示以式(5)的形式生成的平均相關(guān)系數(shù)為0.4、0.8的10 000個模擬區(qū)域的HW檢驗結(jié)果的均值，與圖2(e)相同；而圓點虛線代表的人造區(qū)域中，僅有1/2的站點間具有0.8的相關(guān)系數(shù)，根據(jù)式(6)，其區(qū)域平均相關(guān)系數(shù)為0.4。圖3(b)與圖3(a)類似，其圓點虛線代表的人造區(qū)域中 的站點間擁有0.6的相關(guān)系數(shù)，區(qū)域平均相關(guān)系數(shù)為0.2。圖3(a)、(b)中，圓點虛線幾乎與黑色實線重合，而與三角形點虛線相距較遠。由此不難發(fā)現(xiàn)，在區(qū)域平均相關(guān)系數(shù)相同時，相關(guān)性對H值的影響是近似相同的，而與站點間相關(guān)性空間分布無關(guān)。

基本指令是可編程控制編程中的重要組成部分，有 LD、LD NOT、AND、AND NOT、OR 等?；局噶疃嘤糜诤唵慰刂?，編程較為靈活；對于繁瑣的控制，若僅使用基本指令，編程較復雜，可用功能指令配合編程使其程序簡化，且所編的程序較易閱讀。

圖3 不同相關(guān)系數(shù)分布對H值的影響

5 對Hosking-Wallis一致性檢驗的修正

5.1 修正Hosking-Wallis檢驗的算法及步驟

由第4節(jié)知,HW檢驗的偏差主要來自研究區(qū)域資料間固有的相關(guān)性與檢驗中生成的人工資料序列間相關(guān)性的不同。而即使站點資料序列間相關(guān)系數(shù)分布情況不同，只要整個區(qū)域的平均相關(guān)系數(shù)相同，其對HW檢驗結(jié)果的影響就基本相同。基于以上結(jié)論，對HW檢驗進行修正，使檢驗中生成的人工資料擁有與研究區(qū)域相同的平均相關(guān)系數(shù)，改進后的HW檢驗步驟如下。

(7)

式中：ρij為研究區(qū)域相關(guān)系數(shù)矩陣第i行第j列的元素；N為區(qū)域站點數(shù)。

(2)由式(5)得到N維平均相關(guān)系數(shù)矩陣R。

(3)計算研究區(qū)域各站點的各階線性系數(shù)和區(qū)域的各階平均線性離差系數(shù)。

(4)根據(jù)式(2)計算統(tǒng)計量V。

(5)生成與研究區(qū)域相同站點數(shù)、資料長度的人工區(qū)域資料yik，i=1,2,…,N，k=1,2,…,ni。ni為第i個站點的資料長度。運用正態(tài)Copula函數(shù)，使人工區(qū)域中站點資料服從多元正態(tài)分布，擁有均值零，協(xié)方差矩陣R。

(6)由(3)中求得的地區(qū)平均線性矩系數(shù)擬合4參數(shù)kappa分布。將資料yik轉(zhuǎn)化為符合擬合的kappa分布的資料Qik：

Qik=KAP[φ(yik)]

(8)

式中：KAP為擬合得到的kappa函數(shù)的分位數(shù)函數(shù)；Φ為標準正態(tài)分布的累積分布函數(shù)。

(7)對人工區(qū)域資料序列Qik計算Vsim。

(8)重復1 000次步驟(5)～(7)，計算得到的1 000個Vsim的均值μV和標準差σV，并根據(jù)式(3)計算得到修正后的統(tǒng)計量H*。

5.2 資料相關(guān)性對修正前后Hosking-Wallis檢驗的影響

圖4 改進前后HW檢驗結(jié)果對比

6 結(jié) 語

(1)利用蒙特卡洛模擬實驗，得出HW一致性檢驗中生成的人工資料序列間的平均相關(guān)系數(shù)隨資料長度的分布，并將其與實測資料進行對比。結(jié)果表明人工資料間的相關(guān)性均大于實測資料25%以上。

(2)從HW檢驗的定義出發(fā)，定性分析資料相關(guān)性影響H統(tǒng)計量值的原因，并通過蒙特卡羅模擬實驗，定量分析資料相關(guān)性以及不同的相關(guān)性分布對H值大小的影響。定量分析發(fā)現(xiàn)當資料間存在較大相關(guān)性時，HW檢驗對異質(zhì)性的檢驗能力大幅下降，而當區(qū)域平均相關(guān)系數(shù)相同時，相關(guān)性導致H值的減小量相同。

(3)針對雨量資料間的平均相關(guān)系數(shù)，對HW檢驗進行修正。實踐證實，資料相關(guān)性對修正后HW檢驗結(jié)果的影響大大降低，同時修正后的HW檢驗依舊保留了其對區(qū)域一致性的鑒別能力。這種新的一致性檢驗方法適用條件等同HW檢驗，無需額外調(diào)整，計算速度快，可直接運用于地區(qū)線性矩頻率法，幫助劃分水文一致區(qū)。相較傳統(tǒng)HW檢驗更為符合實測雨量資料情況、檢驗結(jié)果更加可靠。

(4)修正后的一致性檢驗通過優(yōu)化水文一致區(qū)劃分，提升了針對各水文時間序列資料，特別是降雨、徑流序列頻率分析的效率及其結(jié)果的準確性。將改進后的降雨、徑流頻率分析結(jié)果與氣象降水量預報或流域水文模型徑流量預報相比較，能夠更為準確地劃定預報降水或洪水的量級、重現(xiàn)期，為洪水預警預報提供更加堅實的理論基礎(chǔ)。

□