周研來,梅亞東,張代青

(武漢大學(xué)水資源與水電工程科學(xué)國家重點實驗室,湖北 武漢 430072)

水文隨機模擬[1-2],是指根據(jù)水文系統(tǒng)觀測資料的統(tǒng)計和隨機變化規(guī)律,建立能預(yù)估系統(tǒng)未來水文情勢的隨機模型,由隨機模型通過統(tǒng)計試驗獲得大量的模擬序列,再進行水文系統(tǒng)分析計算,解決系統(tǒng)的規(guī)劃、設(shè)計、運行與管理問題的方法。水文隨機模擬技術(shù)的關(guān)鍵是建立合理可靠的隨機水文模型。目前,國內(nèi)外研究的隨機水文模型主要包括回歸類模型[3]、解集類模型[4]、非參數(shù)類模型[5-6]和非線性類模型[7-9]等。Sharma等[5]用核密度估計方法估計序列的條件分布,在核函數(shù)選為高斯函數(shù)時導(dǎo)出了基于非參數(shù)方法進行隨機模擬的NP(1)模型。王文圣等[6]又將NP(1)模型推廣到單變量和多變量的多階情形,即NP(p)模型,并研究了模型在徑流模擬中的可行性。肖義等[10]和張濤等[11]就洪水特征變量間的相關(guān)性基于Copula函數(shù)提出了一種洪水過程隨機模擬方法。閆寶偉等[9]建立了基于Copula函數(shù)的一階季節(jié)性隨機模型,嘗試將其應(yīng)用于徑流的隨機模擬,所建模型能較好地保持原序列的偏態(tài)特性、非線性相關(guān)性和概率密度特征等統(tǒng)計特性。一階自回歸模型[3]和典型解集模型[4]由于構(gòu)造簡單,在水文模擬中備受青睞,為避開一階自回歸模型和典型解集模型中對水文序列的正態(tài)轉(zhuǎn)化(該過程難免會出現(xiàn)部分信息失真),筆者嘗試用Copula函數(shù)構(gòu)造一種新的徑流過程解集隨機模型,即先用Copula函數(shù)構(gòu)造年徑流量與最大月徑流量之間的聯(lián)合分布,通過聯(lián)合分布達到隨機模擬年徑流量的目的,再優(yōu)選實測徑流過程,將年徑流量分解為各分量的徑流過程。

1 基于Copula函數(shù)的解集隨機模型

1.1 Copula函數(shù)的構(gòu)造

Copula函數(shù)是定義域為[0,1]均勻分布的多維聯(lián)合分布函數(shù),由Sklar定理可以將多個隨機變量的分布函數(shù)連接起來構(gòu)造聯(lián)合分布,對于二維情況,可以表述為[8,12]

式中:H(x,y)為隨機變量X和Y的聯(lián)合分布函數(shù);Cθ為參數(shù)為θ的Copula函數(shù);FX(x),FY(y)分別為隨機變量X和Y的邊緣分布函數(shù)。

Copula函數(shù)在多變量水文分析計算中的應(yīng)用及研究進展,文獻[12]已有詳細(xì)介紹。

徑流過程的2個特征量(最大月徑流量和年徑流量)具有一定的相關(guān)性。兩變量耿貝爾聯(lián)合(Gumbel-Hougaard Copula)函數(shù)可用于描述具有上尾相關(guān)特性多變量之間的相關(guān)關(guān)系[13]。例如它可以很好地描述當(dāng)一個變量增幅較大時,另一個與之存在相關(guān)性的變量隨之明顯增大的多變量相關(guān)關(guān)系[10]。基于Gumbel-Hougaard Copula函數(shù)的數(shù)學(xué)特性,采用該函數(shù)來描述最大月徑流量和年徑流量,其數(shù)學(xué)表達式為[8,10,13]

式中:FX(x)和FY(y)均為PⅢ型分布函數(shù);θ越大,說明u,v的相關(guān)性越強,θ的估計見文獻[9]。

1.2 模型算法

許多水文要素是累加而成的,如年水量是由月水量累加而得,15 d洪水總量是由15 d的日洪水量累計而得。一般而言,將各分量累加起來的總量叫做聚集;反之,從總量中可以分解出分量。解集隨機模型的實質(zhì)在于將總量隨機解集成各分量,其顯著的特點在于保持水量平衡,即各分量的水量相加嚴(yán)格等于總水量。Copula函數(shù)對最大月徑流量 X和年徑流量Y進行聯(lián)合描述,通過聯(lián)合分布的隨機抽樣方法對存在相關(guān)性的最大月徑流量和年徑流量成對地取樣,從而達到隨機模擬最大月徑流量和年徑流量的目的。參照文獻[12]的峰量隨機抽樣方法以及典型解集模型[14]得出基于Copula函數(shù)的解集隨機模型模擬徑流過程的步驟如下:

步驟1 根據(jù)最大月徑流量和年徑流量的聯(lián)合分布可以得到當(dāng)最大月徑流量為某指定值時年徑流量的條件分布:

步驟2 產(chǎn)生服從[0,1]均勻分布的2個獨立隨機數(shù)r1和r2。

步驟3 令r1為最大月徑流量發(fā)生的不超過x值的概率FX(x),即 u=FX(x)=r1,根據(jù) x=(u)得到x值。

步驟4 令r2為最大月徑流量X=x時年徑流量y的條件概率分布值,即r2=CX(y X =x),從而可根據(jù)y=(r2X=x)計算得到 y。

步驟5 由步驟3和步驟4即可抽取1對最大月徑流量和年徑流量的組合(x,y)。

步驟6 以模擬的最大月徑流量和年徑流量組合(x,y)為控制條件,選擇最大月徑流量和年徑流量的比值與x/y值接近的實測徑流過程,并采用解集[11]方法進行分解,得到模擬的徑流過程:

式中:Xi(t),xi(t)分別為第i年第t個截口的模擬徑流值和實測徑流值;Yi為第i年的模擬年徑流量;A為將年徑流量(億m3)轉(zhuǎn)換為總流量(m3/s)的單位換算系數(shù);Ki(t)為第i年第t個截口的分解系數(shù)。

步驟7 重復(fù)步驟2至步驟6共M次,即可模擬出M年的徑流過程。

2 模型應(yīng)用研究

利用黃河中上游河口鎮(zhèn)水文站1919—1988年共70年的實測月徑流量資料,檢驗基于Copula函數(shù)的解集隨機模型的可行性和有效性。由于所選資料序列較長,某種意義上可以代表整體的統(tǒng)計特性,以此來驗證模型,可信度較高。同時選取季節(jié)性一階自回歸模型和典型解集模型作為對照,以說明本文所建模型的實用性。

在實測徑流資料中,每年選擇1對最大月徑流量和年徑流量數(shù)據(jù),構(gòu)成樣本容量為70的聯(lián)合觀測值系列,并據(jù)此建立最大月徑流量X和年徑流量Y的兩變量聯(lián)合分布。X和Y的邊緣分布采用PⅢ型分布,分別為FX(x)和 FY(y),采用線性矩法估計FX(x)和FY(y)的參數(shù),得到實測年徑流量的統(tǒng)計特征值,見表1。將估計的參數(shù)假定為河口鎮(zhèn)站最大月徑流量和年徑流量總體分布的參數(shù)。依據(jù)Kendall秩相關(guān)系數(shù)和Copula函數(shù)的參數(shù) θ的解析關(guān)系[8,12],得到最大月徑流量和年徑流量的參數(shù)θ=4.25。圖1為聯(lián)合觀測值的經(jīng)驗分布和理論分布的相關(guān)關(guān)系。據(jù)圖1可知,θ=4.25時,最大月徑流量和年徑流量聯(lián)合觀測值的經(jīng)驗分布和理論分布擬合較好,θ=4.25可作為聯(lián)合總體分布的參數(shù)估計量,所建立的X和Y的聯(lián)合分布是合理的。由Sklar定理及式(2)得到X和Y的聯(lián)合分布F(x,y):

?

圖1 聯(lián)合觀測值經(jīng)驗分布和理論分布的比較

繪制最大月徑流量和年徑流量的兩變量聯(lián)合分布函數(shù)及概率密度函數(shù)三維圖,如圖2、圖3所示。

圖2 X和Y的聯(lián)合分布函數(shù)

圖3 X和Y的聯(lián)合概率函數(shù)

此外,從聯(lián)合分布隨機抽取7萬對最大月徑流量和年徑流量數(shù)據(jù),點繪得到它們的頻率散點圖,如圖4所示。由圖4可見,當(dāng)模擬的最大月徑流量為某一頻率時,可能發(fā)生多種頻率的年徑流量,最大月徑流量和年徑流量之間不存在嚴(yán)格的一一對應(yīng)關(guān)系;從回歸關(guān)系和相關(guān)關(guān)系來看,最大月徑流量和年徑流量頻率的回歸系數(shù)為0.92,相關(guān)系數(shù)為0.85,說明最大月徑流量和年徑流量組合既存在較好的同頻率特性又存在一定的隨機性;從頻率分布狀況來看,最大月徑流量和年徑流量偏小值的頻率居多,中等值的頻率次之,而較大值頻率分布較少,這與水文要素的枯、平、豐水年的頻率分布狀況十分吻合。

圖4 隨機模擬X和Y頻率的散點圖

根據(jù)圖4中任意一點對應(yīng)的最大月徑流量和年徑流量,可求得最大月徑流量與年徑流量的比值,依據(jù)該比值選擇具有相似最大月徑流量與年徑流量之比的實測徑流過程,采用典型解集方法進行分解,模擬多年月徑流過程(如7萬年),其年徑流量的統(tǒng)計參數(shù)見表1。

統(tǒng)計分析河口鎮(zhèn)水文站1919—1988年實測的月徑流過程資料,可知1967年為豐水年(最大月徑流量為3401m3/s,年徑流量為420億m3),見圖5。圖5還給出了以1967年實測月徑流過程為基礎(chǔ)進行縮放后得到的2條月徑流過程線。由圖5可見,模擬的徑流過程在形狀上與實測過程十分相似。因此,隨機模擬的徑流過程系列能夠反映實測徑流的形狀。

圖5 以1967年實測徑流過程模擬的徑流過程線

對基于Copula函數(shù)的解集隨機模型是否能很好地保持實測系列的年徑流量統(tǒng)計特征進行檢驗,并與水文領(lǐng)域常用的季節(jié)性一階自回歸模型和典型解集模型[14]進行比較。對于每個模型,在其模擬的徑流過程系列對應(yīng)取出7萬個年徑流量,分別組成1000個樣本容量為70的年徑流量樣本。對每個年徑流量樣本分別估計統(tǒng)計特征參數(shù)(均值,Cv和Cs),并分別計算統(tǒng)計特征參數(shù)的均值。采用綜合相對均方誤差RMSE[9-10]作為模型評價指標(biāo):

評價指標(biāo)RMSE的值越小,說明隨機模型越能夠保持實測系列的統(tǒng)計特征。表1給出了實測系列和模擬系列的年徑流量統(tǒng)計特征值及評價指標(biāo)值。由表1可知,季節(jié)性一階自回歸模型模擬徑流過程的年徑流量統(tǒng)計特征效果相對較差;基于Copula函數(shù)的解集隨機模型的RMSE均低于季節(jié)性一階自回歸模型和典型解集模型的RMSE。究其原因,基于Copula函數(shù)的解集隨機模型避開了季節(jié)性一階自回歸模型和典型解集模型中對水文序列的正態(tài)轉(zhuǎn)化(難免會出現(xiàn)部分信息失真),而是直接隨機模擬年徑流總量,然后選用實測徑流過程分解為模擬的徑流過程,因此能夠很好地保持實測系列的徑流過程形狀及年徑流量的統(tǒng)計特征。

3 結(jié) 語

構(gòu)造了最大月徑流量和年徑流量的聯(lián)合分布,采用基于聯(lián)合分布的隨機抽樣方法模擬了大樣本的徑流過程,從而提出了一種新的徑流過程隨機模擬方法,建立了基于Copula函數(shù)的解集隨機模型。該模型避開了季節(jié)性一階自回歸模型和典型解集模型中對水文序列的正態(tài)轉(zhuǎn)化和相依形式的主觀假定,很好地考慮了最大月徑流量和年徑流量之間的相關(guān)性及兩者頻率組合的隨機性,選用實測徑流過程作為隨機模擬徑流過程的時程分配形式,能很好地保持實測系列的徑流過程形狀及年徑流量的統(tǒng)計特征。實例計算結(jié)果表明,基于Copula函數(shù)的解集隨機模型的模擬效果至少可以與常用的典型解集模型相當(dāng),適用于水文時間序列的隨機模擬,從而為水文水資源的隨機模擬提供了一種新途徑。

[1]芮孝芳,劉方貴,邢貞相.水文學(xué)的發(fā)展及其所面臨的若干前沿科學(xué)問題[J].水利水電科技進展,2007,27(1):75-79.

[2]胡慶芳,王銀堂,劉克琳,等.基于改進的兩參數(shù)月水量平衡模型的月徑流模擬[J].河海大學(xué)學(xué)報:自然科學(xué)版,2007,35(6):638-642.

[3]葛朝霞,王會容.黃河上游月徑流量與海溫的關(guān)系研究[J].河海大學(xué)學(xué)報:自然科學(xué)版,2003,31(6):605-608.

[4]KOUTSOYIANNIS D,MANETAS A.Simple disaggregation accurate adjusting procedures[J].Water Resources Research,1996,32(7):2105-2117.

[5]SHARMA A,TARBOTON D G,LALL U.Streamflow simulation:a nonparametric approach[J].Water Resources Research,1997,33(2):192-308.

[6]王文圣,丁晶.基于核密度估計的多變量非參數(shù)隨機模型初步研究[J].水利學(xué)報,2003,34(2):9-14.

[7]PRAIRIE J R,RAJAGOPALAN B,FULP T J,et al.Modified K-NN modelfor stochastic streamflow simulation[J].Journal of Hydrologic Engineering,2006,11(4):371-378.

[8]SALVADORI G,de MICHLELE C.Frequency analysis via copulas theoretical aspects and applications to hydrological events[J].Water Resources Research,2004,40(12):1-17.

[9]閆寶偉,郭生練,劉攀,等.基于Copula函數(shù)的徑流隨機模擬[J].四川大學(xué)學(xué)報:工程科學(xué)版,2010,42(1):5-9.

[10]肖義,郭生練,熊立華,等.一種新的洪水過程隨機模擬方法研究[J].四川大學(xué)學(xué)報:工程科學(xué)版,2007,39(2):55-60.

[11]張濤,趙春偉,雒文生.基于Copula函數(shù)的洪水過程隨機模擬[J].武漢大學(xué)學(xué)報:工學(xué)版,2008,41(4):1-4.

[12]郭生練,閆寶偉,肖義,等.Copula函數(shù)在多變量水文分析計算中的應(yīng)用及研究進展[J].水文,2008,28(3):1-7.

[13]NELSON R B.An introduction to Copula[M].New York:Springer,1999.

[14]中華人民共和國能源部,中華人民共和國水利部.水利水電工程設(shè)計洪水計算手冊[M].北京:中國水利水電出版社,2001.