張繼鵬，彭 楊，時(shí)玉龍，趙曉東，丁夢(mèng)霞

(1.華北電力大學(xué) 可再生能源學(xué)院，北京 102206；2.華北水利水電大學(xué) 水利學(xué)院，鄭州 450045)

0 引 言

含沙量是一個(gè)主要的河流特征變量，綜合反映河流的徑流與輸沙特性。含沙量變化直接影響水庫(kù)泥沙淤積、河道沖淤、水質(zhì)及水生環(huán)境，因此含沙量隨機(jī)模擬對(duì)流域水沙調(diào)控、水質(zhì)與水環(huán)境管理等方面具有重要意義。

自回歸模型和解集模型是目前常用的水文隨機(jī)模型[1,2]，但這兩類模型的共同點(diǎn)是采用線性函數(shù)關(guān)系描述水文現(xiàn)象[3]，難以反映水文時(shí)間序列復(fù)雜的非線性特征。近年來(lái)非線性、非參數(shù)隨機(jī)模型成為水文隨機(jī)模擬研究的熱點(diǎn)，提出了非線性自回歸模型[4,5]、人工神經(jīng)網(wǎng)絡(luò)模型[6,7]、基于核密度估計(jì)的非參數(shù)隨機(jī)模型[3,8]以及非參數(shù)解集模型[9,10]等，為水文非線性時(shí)間序列模擬提供了多種途徑。水文非線性時(shí)間序列模擬的核心問(wèn)題是描述時(shí)間序列相鄰時(shí)段的非線性關(guān)系。Copula函數(shù)是一種描述變量相關(guān)關(guān)系的有效方法，由于其對(duì)變量的邊緣分布沒(méi)有限制，構(gòu)造聯(lián)合分布較為靈活，目前已廣泛應(yīng)用于水文多變量分析中，如洪水頻率[11,12]、降雨頻率[13,14]、干旱特征分析[15,16]、洪水遭遇[17-19]等。相比較而言，Copula函數(shù)在水文隨機(jī)模擬的應(yīng)用成果不多，且主要集中在降雨[20]、徑流[21,22]和洪水[23]等方面。迄今為止，還未見到Copula函數(shù)應(yīng)用于日含沙量時(shí)間序列的隨機(jī)模擬。因此，本文將采用Copula函數(shù)來(lái)描述相鄰兩日含沙量之間的相關(guān)結(jié)構(gòu)，在此基礎(chǔ)上建立了日含沙量隨機(jī)模擬模型，并將該模型應(yīng)用于長(zhǎng)江屏山水文站日含沙量過(guò)程的隨機(jī)模擬中。

1 基于Copula函數(shù)的日含沙量隨機(jī)模擬模型

1.1 Copula函數(shù)

Copula函數(shù)是定義域?yàn)閇0,1]的多維聯(lián)合分布函數(shù)。依據(jù)Sklar定理[24]，一個(gè)二維聯(lián)合分布函數(shù)可以由兩個(gè)獨(dú)立的邊緣分布函數(shù)和一個(gè)Copula連接函數(shù)組成，即

F(x,y)=Cθ[FX(x),FY(y)]=Cθ(u,v)

(1)

式中：Cθ為Copula函數(shù)；θ為Copula函數(shù)的參數(shù)，u=FX(x)，v=FY(y)分別為隨機(jī)變量X和Y的邊緣分布。

Copula函數(shù)分為橢圓型、二次型和Archimedean型三大類，其中Archimedean Copula函數(shù)結(jié)構(gòu)相對(duì)簡(jiǎn)單，計(jì)算方便，可以構(gòu)造出多種形式、適應(yīng)性強(qiáng)的組合變量聯(lián)合分布函數(shù)，在水文領(lǐng)域中應(yīng)用最多。常見的二維Archimedean Copula主要有[24,25]：

(1)Gumbel-Hougaard (GH)Copula函數(shù)。

(2)

(2)Clayton Copula函數(shù)。

C(u,v)=(u-θ+v-θ-1)-1/θ,θ∈(0,∞)

(3)

(3)Frank Copula函數(shù)。

(4)

(4)Ali-Mikhail-Haq (AMH)Copula函數(shù)。

(5)

1.2 相鄰兩日含沙量聯(lián)合分布

根據(jù)式(1)，相鄰兩日含沙量的聯(lián)合分布可表示為：

FS(st,st-1)=Cθt[Fst(st),FSt-1(st-1)]=Cθt(ut,ut-1)

(6)

式中：st-1和st分別為第t-1日和第t日的含沙量；FS(st,st-1)為st-1和st的聯(lián)合分布；ut-1=FSt-1(st-1)和ut=FSt(st)分別為st-1和st的邊緣分布；θt為第t日的Copula函數(shù)的參數(shù)，根據(jù)相關(guān)性指標(biāo)法[25]，可利用其與Kendall秩相關(guān)系數(shù)τ之間的關(guān)系求出；第t日的Kendall秩相關(guān)系數(shù)τt采用下式計(jì)算：

(7)

1.3 基于條件分布的隨機(jī)模擬模型

根據(jù)式(1)，第t日的含沙量st的條件分布可表示為：

(8)

式中：Gθt(ut,ut-1)為條件Copula函數(shù)。由此可得到基于Copula函數(shù)的日含沙量隨機(jī)模擬模型，即：

(9)

根據(jù)式(9)就可以實(shí)現(xiàn)對(duì) 的隨機(jī)模擬，具體步驟如下：

(1)采用正態(tài)分位數(shù)變換法對(duì)實(shí)測(cè)含沙量數(shù)據(jù)進(jìn)行標(biāo)準(zhǔn)正態(tài)化；

2 實(shí)例應(yīng)用

選取金沙江下游屏山站1955-2010年共56年的實(shí)測(cè)日含沙量資料作為建模的基本資料，采用所建模型對(duì)屏山站日含沙量進(jìn)行隨機(jī)模擬。屏山站是金沙江干流出口控制站，流域控制面積458 592 km2，多年平均徑流量1 460 億m3，多年平均輸沙量2.39 億t，多年平均含沙量1.66 kg/m3。為簡(jiǎn)化計(jì)算，模擬時(shí)每年2月份均取28 d。

2.1 日含沙量聯(lián)合分布建立

圖1 各日與的秩相關(guān)系數(shù)Fig.1 Daily rank correlation coefficients of

(10)

式中：Pei，Pi分別為經(jīng)驗(yàn)聯(lián)合頻率和理論聯(lián)合頻率；i和n分別為樣本序號(hào)和樣本容量。

各日的OLS值與對(duì)應(yīng)的最優(yōu)Copula函數(shù)分別如圖2和圖3所示。由圖2可見，在一年內(nèi)，OLS值的最大值為0.025 1，最小值為0.008 1，平均值為0.013 7，說(shuō)明所選定的二維Copula聯(lián)合分布計(jì)算出的理論聯(lián)合分布概率與經(jīng)驗(yàn)頻率擬合較好。圖4給出了屏山站日OLS值接近平均值時(shí)聯(lián)合分布的經(jīng)驗(yàn)概率與理論概率的擬合情況，從圖4中可以看出，二者擬合情況較為理想。

2.2 日含沙量的模擬及檢驗(yàn)

采用上述模型對(duì)屏山站日含沙量過(guò)程進(jìn)行了長(zhǎng)序列模擬，模擬了100組與實(shí)測(cè)樣本容量(56 a)相等的日含沙量序列，并對(duì)模擬結(jié)果進(jìn)行檢驗(yàn)。

圖2 各日OLS值Fig.2 Daily OLS values

圖3 各日最優(yōu)Copula函數(shù)Fig.3 Daily best-fit copula functions

圖4 經(jīng)驗(yàn)與理論聯(lián)合概率對(duì)比Fig.4 Comparison of empirical and theoretical joint probabilities

圖5 屏山站模擬與實(shí)測(cè)序列各日統(tǒng)計(jì)參數(shù)Fig.5 Simulated and measured statistics of daily sediment concentrations in the Pingshan station

圖5給出了屏山站實(shí)測(cè)和模擬含沙量序列各日統(tǒng)計(jì)參數(shù)，即均值、標(biāo)準(zhǔn)差Sd、變差系數(shù)Cv、偏態(tài)系數(shù)Cs以及滯時(shí)為1的秩相關(guān)系數(shù)τ。表1統(tǒng)計(jì)了屏山站日含沙量模擬序列各日統(tǒng)計(jì)參數(shù)絕對(duì)誤差和相對(duì)誤差的平均值以及相對(duì)誤差≤10%的通過(guò)率。由圖5和表1可見，模擬序列各日的統(tǒng)計(jì)參數(shù)與實(shí)測(cè)序列吻合較好，各統(tǒng)計(jì)參數(shù)的通過(guò)率均在98%以上，最大平均相對(duì)誤差4.22%，其中均值和秩相關(guān)系數(shù)τ的通過(guò)率均達(dá)到100%，最大相對(duì)誤差分別為0.83%和0.49%。因此，基于二維Copula函數(shù)的屏山站日含沙量過(guò)程隨機(jī)模擬結(jié)果能夠反映實(shí)測(cè)日含沙量序列的統(tǒng)計(jì)特征，模擬精度也較高。

表1 各統(tǒng)計(jì)參數(shù)通過(guò)率及平均相對(duì)誤差Tab.1 Passing rate and average relative error of statistical parameters

2.3 與其他模型比較

為進(jìn)一步說(shuō)明本文所建模型的模擬效果，本文還將模擬結(jié)果與分期平穩(wěn)AR模型的模擬結(jié)果進(jìn)行比較。表2給出這兩種模型的日含沙量模擬序列各統(tǒng)計(jì)參數(shù)與實(shí)測(cè)序列的相對(duì)均方根誤差RMSE[27,28]，其中RMSE采用下式計(jì)算：

(11)

式中：δo,j為實(shí)測(cè)序列第j日的統(tǒng)計(jì)參數(shù)；δs,j為模擬序列第j日的統(tǒng)計(jì)參數(shù)；m為總?cè)諗?shù)。式(11)表明，RMSE值越小，模擬序列與實(shí)測(cè)序列之間的擬合程度越好。

由表2可見，基于Copula函數(shù)的日含沙量模擬的均值和Cv值對(duì)應(yīng)的RMSE值與分期平穩(wěn)AR模型的基本一致，但Sd、Cs和τ值對(duì)應(yīng)的RMSE值明顯小于分期平穩(wěn)AR模型，表示基于Copula函數(shù)的日含沙量模擬序列在偏態(tài)性和非線性相關(guān)關(guān)系方面明顯優(yōu)于分期平穩(wěn)AR模型。因此，與分期平穩(wěn)AR模型相比，基于Copula函數(shù)的日含沙量模擬序列的擬合結(jié)果基本一致，但在偏態(tài)性和非線性相關(guān)關(guān)系方面具有明顯的優(yōu)勢(shì)。

表2 兩種模型統(tǒng)計(jì)參數(shù)RMSE值比較Tab.2 Comparison of RMSE values of statistical parameters in two models

3 結(jié) 論

本文考慮相鄰兩日含沙量之間的相關(guān)性，建立了基于Copula函數(shù)的日含沙量隨機(jī)模擬模型。針對(duì)日含沙量模擬時(shí)段較多，日含沙量邊緣分布難以統(tǒng)一確定的特點(diǎn)，采用正態(tài)分位數(shù)變換對(duì)實(shí)測(cè)含沙量資料進(jìn)行了標(biāo)準(zhǔn)正態(tài)化處理。將該模型應(yīng)用到金沙江下游屏山站日含沙量隨機(jī)模擬，并將計(jì)算結(jié)果與分期平穩(wěn)AR模型計(jì)算結(jié)果進(jìn)行比較，結(jié)果表明，本文所建模型能較好地保持實(shí)測(cè)日含沙量資料的統(tǒng)計(jì)特性，模擬精度較高，模擬序列的均值，標(biāo)準(zhǔn)差，變差系數(shù)，偏態(tài)系數(shù)以及秩相關(guān)系數(shù)τ的通過(guò)率均在98%以上，平均相對(duì)誤差最大不超過(guò)4.5%，且在偏態(tài)性和非線性相關(guān)關(guān)系方面優(yōu)于分期平穩(wěn)AR模型，說(shuō)明本文所建模型可用于隨機(jī)模擬長(zhǎng)序列日含沙量過(guò)程，并能為其他水文資料的長(zhǎng)序列模擬方法提供參考。