吳昊昊，倪 晉

(1.安徽省·水利部淮河水利委員會水利科學研究院，安徽 蚌埠 233000；2.安徽省水利水資源重點實驗室，安徽 蚌埠 233000；)

1 研究背景

參數模型大多數基于模型驅動，并且對于水文序列的概率分布和相關關系做了一定程度的假設。為盡可能接近于歷史序列的真實分布，水文學者們提出大量非參數隨機模型，從數據驅動出發(fā)，避免了參數估計以及對研究對象概率分布和相依結構的假定，能夠較好地捕捉水文序列間的相依特性。

非參數核密度估計模型作為非參數模型的重要研究方向，近年來已取得較大發(fā)展。Sharma等[1]假定徑流為一個具有時間依賴性的馬爾柯夫過程，利用條件概率密度函數的核密度估計生成模擬徑流序列，提出非參數一階馬爾柯夫模型(NP(1))。結果表明這種基于數據驅動的非參數徑流合成方法，比隨機水文學中使用的傳統(tǒng)模型更靈活，并且能夠再現線性和非線性相依關系。王文圣等從數據驅動角度出發(fā)，建立了單變量多階核密度估計模型，在屏山站日徑流過程模擬中證明了其適用性[2]；之后利用該模型展開金沙江流域李莊-屏山區(qū)間年徑流過程的模擬研究，通過和自回歸模型對比驗證了非參數模型的優(yōu)越性[3]。Sharma等[4]基于上述模型無法在較大的(季節(jié)性到年際性)時間滯后中保留相依關系，尤其是對觀測到的分布特征如概率密度函數中強不對稱或多模態(tài)性等表現不佳，使用可變核和徑流聚集變量改進了NP模型，提出了一種月徑流序列年際相關性的非參數模型(NPL)。王文圣和丁晶[5]嘗試建立多變量非參數模型，并隨機生成屏山站和宜賓-屏山區(qū)間兩站日徑流模擬序列，實例應用表明該模型的模擬效果較好。陳大春[7]利用粒子群算法優(yōu)化基于最小二乘交叉法建立的帶寬系數目標函數和可變核帶寬方法，構建了非參數核密度估計模型，并應用于烏魯木齊河月徑流模擬。

考慮到NPL模型在國內的應用研究較少，且現有文獻幾乎未深入探討過NP和NPL在徑流模擬方面的適用性對比，本文以黃河流域蘭州、龍門和白馬寺3個水文站點的天然月徑流序列資料為研究對象，分別建立各站點月徑流的NP和NPL模型，比較分析兩模型在研究區(qū)的適用效果，以期為水利工程規(guī)劃設計管理提供依據。

2 模型原理

2.1 核密度估計理論

從概率密度函數f(x)未知總體中抽取獨立同分布樣本，當隨機變量X為d維時，總體密度f(x)的核密度估計定義為[9]：

(1)

式中，X=(X1，X2，…，Xd)T；Xi=(Xi1，Xi2，…，Xid)T(i=1，2，…，n)；d—隨機變量X的維數；S—向量X的d×d維對稱樣本協方差矩陣；K(u)—核函數；h—帶寬系數；n—樣本容量。

本文中選用標準高斯函數作為核函數，最小二乘交叉證實檢驗(LSCV)法尋求帶寬系數h[10]，公式為：

(2)

式中，Lij=(Xi-Xj)TS-1(Xi-Xj)；Xi，Xj數據不同；S—數據樣本協方差陣；p—模型階數；n—數據組數。

2.2 非參數p階馬爾科夫徑流模型(NP)

當xt(t=1，2，…，n)為單變量水文相依時間序列，且xt依賴于前p個值xt-1，xt-2，…，xt-p，令Vt=(xt-1，xt-2，…，xt-p)T，則xt的條件概率密度函數為[11]：

(3)

式中，fV(Vt)—p維邊緣概率密度函數；f(xt，Vt)—p+1維聯合分布密度函數。基于高斯核密度估計，得到：

(4)

(5)

式中，S—(xt，Vt)的(p+1)×(p+1)階對稱樣本協方差矩陣；Sx—xt的樣本方差；SxV—xt與Vt的1×(p+1)階樣本協方差矩陣；SV—Vt的p×p階對稱樣本協方差矩陣。Vt=(xt-1，xt-2，…，xt-p)T和xt來自于實測樣本數據，其中t=p+1，p+2，…，n。

將式(4)、(5)代入式(3)整理得到

(6)

式中，

(7)

(8)

(9)

(10)

應用式(6)能夠隨機模擬序列xt，模擬公式如式(11)：

(11)

式中，et—均值為0，方差為1的獨立高斯隨機變量。模型階數p根據AIC準則確定。

2.3 NPL模型

為精確表達模擬徑流月份之間和年際之間的相關性，Sharma[4]提出一種季節(jié)性徑流序列隨機模擬方法(NPL)，它引入了確保準確描述季節(jié)與年際相關性的聚集變量和可變核，試圖在模擬的徑流序列中再現這種長期相關特性。

考慮時間t處的徑流量為xt，例如月徑流前12個月的月徑流可以表示為x1，x2，…，x12，接下來的12個月可以表示為x13，x14，…，x24，以此類推。因此聚集徑流變量Zt定義為：

(12)

式中，m—聚集水平，代表模擬月份過去m個月的徑流量和。本文選用模型階數p=1、聚集水平m=12進行介紹，則用于模擬的條件概率密度如下：

(13)

式中，fm(·)—變量集的邊緣概率密度。用式(1)進行高斯核密度估計，變量集(Xt，Xt-1，Zt)的聯合概率密度可以表示為

(14)

式中，h—帶寬系數；n—樣本容量；S—變量集(Xt，Xt-1，Zt)的協方差矩陣，其中S11—Xt的樣本方差，S12—Xt與Xt-1的樣本協方差，S1z—Xt和Zt的樣本協方差，S22—Xt-1的樣本方差，Szz—Zt的樣本方差，S2z—Xt-1和Zt的樣本協方差。

式(14)相當于n個正態(tài)概率密度函數的比例和，根據多元正態(tài)分布的已知關系，可以進一步簡化為

(15)

其中，

(16)

(17)

式中，參數同式(14)意義一致。

類似地，(Xt-1，Zt)的聯合密度函數估計為

(18)

式中，參數同式(14)意義一致。將式(15)和式(18)代入式(13)，整理可得

(19)

其中，

(20)

式中，bi—與每個內核相關聯的條件平均數，如式(16)所示；ωi—構成條件概率密度每個核相關的權重；其他參數同式(14)意義一致。

應用式(21)能夠隨機模擬序列Xt，模擬公式為：

(21)

式中，Et—均值為0，標準差為1的高斯隨機變量。

由于高斯核函數自身具有對稱無界特點，可能會大量導致部分負徑流值的生成，對模擬序列造成邊界影響。NPL模型采用可變核修正方法對帶寬系數進行優(yōu)化，具體修正步驟見文獻[13]。

3 實例應用

3.1 資料的選取與審查

本文以黃河流域蘭州、龍門和白馬寺3個水文站點的天然月徑流序列資料為研究對象，開展NP和NPL模型在徑流模擬中的應用研究，水文數據均滿足可靠性、代表性和一致性要求，可用于模型的相關分析計算。

3.2 模型參數

各站點的模型階數根據AIC準則識別，經計算蘭州、龍門和白馬寺的模型階數取p=1，因此建立黃河流域3個水文站點的一階NP和NPL模型應用于月徑流隨機模擬。根據式(2)采用優(yōu)化方法最小化LSCV求得NP和NPL模型各站點各月優(yōu)化帶寬h，結果見表1。

表1 NP和NPL模型各站月徑流序列各月優(yōu)化帶寬h估計

3.3 模型檢驗和對比分析

本節(jié)采用短序列法檢驗模型適用性，針對3個水文站點分別生成300組樣本，每組樣本容量對應站點序列長度。為探討分析模擬效果，NP和NPL模型的適用性檢驗和對比分析利用均值、標準差、偏態(tài)系數等統(tǒng)計值，部分統(tǒng)計特征參數如圖1—5所示，見表2—4。

圖1 各站NP和NPL模型均值

表2 蘭州站NP和NPL模型月徑流統(tǒng)計特征值

表3 龍門站NP和NPL月徑流統(tǒng)計參數

表4 白馬寺站NP和NPL月徑流統(tǒng)計參數

由圖1可得：均值方面，蘭州站NPL模型模擬效果要優(yōu)于NP模型，而在龍門站和白馬寺站，NP模型模擬序列的各月均值要更接近于實測序列，在部分月份NPL模型模擬效果不佳。分析表2—4中的均值統(tǒng)計特征能夠得到：蘭州站NPL模型均值未控制在兩個均方差標準下的月份達50%，有33%的月份均值在一個均方差標準下；NP模型超出2個均方差標準下的月份比例要略低于NPL模型，表明對于蘭州站NP模型和NPL模型的模擬序列沒有很好地保持實測的均值統(tǒng)計特性。龍門站NPL模型有50%的月徑流模擬序列均值未能控制在兩個均方差下，控制在一個均方差標準下的月份達33%；NP模型模擬效果要優(yōu)于NPL模型，有50%的月徑流模擬序列控制在一個均方差下。白馬寺站NP模型除3月和4月外其余月份的均值均控制在一個均方差下，而NPL模型未能控制在兩個均方差標準下的月份占41%，總體上來看NPL在均值統(tǒng)計特性方面模擬效果不如NP模型。

如圖2所示，標準差方面，在蘭州站、龍門站和白馬寺站，非汛期NP和NPL模型模擬序列的標準差與實測序列基本能保持一致，汛期NP模型表現更好，實測序列的標準差處于模擬樣本中位數附近浮動。結合表2—4中的標準差統(tǒng)計特征能夠得到：蘭州站NP模型占75%的月份均值能控制在一個均方差標準下，相比NPL模型能夠控制在兩個均方差標準下的月份僅占25%，表明NP模型相比較更適合蘭州站用于徑流模擬。龍門站兩模型模擬序列標準差除個別月份基本都能較好的控制在兩個均方差標準下，與實測序列的統(tǒng)計特征值差距較小。白馬寺站NPL模型模擬序列均值控制在兩個均方差標準下的月份占41.7%，NP模型除3月份其余月份均控制在2個均方差標準下，說明白馬寺站NP模型相比NPL模型對歷史實測序列的標準差統(tǒng)計特征保持得更好。

圖2 各站NP和NPL模型標準差

由表2—4中的偏態(tài)系數統(tǒng)計特征能夠得到：蘭州站NPL模型無法控制在兩個均方差標準下的月份達50%，NP模型除8月外其余月份均能控制在兩個均方差標準下，模擬效果更加理想。龍門站NPL模型的偏態(tài)系數在大多數月份能夠控制在兩個均方差標準下，NP模型除1月、2月和3月外其余月份均良好地保持在一個均方差標準下。白馬寺站NPL模型超出兩個均方差標準的月份達41.7%，NP模型有75%的月徑流序列偏態(tài)系數控制在兩個均方差標準下。

由兩個模型的二階自相關系數圖(如圖3所示)和表2—4中的R1和R2統(tǒng)計特征能夠得到：R1方面，NP模型和NPL模型模擬序列的統(tǒng)計特征值與實測序列R1均有較大偏差；根據計算NPL模型各站點的模擬序列一階自相關系數控制在兩個均方差標準下的月份占58.3%及以上，NP模型在蘭州站和白馬寺站基本保持R1統(tǒng)計特征的月份為58.3%左右。R2方面，模擬效果總體要優(yōu)于一階自相關系數，蘭州站NP模型各月份的R2基本位于盒箱范圍內，而NPL模型模擬效果不佳。根據表中數據計算，發(fā)現蘭州站NP模型控制在兩個均方差標準下的月份比例比NPL模型要高16.7%，龍門站兩模型基本相似，白馬寺站NPL模型控制在兩個均方差標準下的月份比例要高8.3%，與箱型圖結果基本一致?？傮w來看，對于各站NP和NPL模型二階自相關系數的模擬效果要優(yōu)于一階。

圖3 各站NP和NPL模型二階自相關系數R2

圖4說明了滯后1個月月徑流量與之前12個月的總徑流量之和之間的相關性。總的來看，在滯后一個月的情況下，NPL模型比NP模型較好地再現了這些相關性，如蘭州站的12月，龍門站的1月、2月，白馬寺站的10月等。但在白馬寺站部分月份和其余站點的個別月份，兩模型均有與實測徑流序列統(tǒng)計特征值差距較大的情況發(fā)生。滯后2月觀察到類似的結果，這表明在NPL模型模擬中，長期依賴性得到了恰當的表示。

圖4 各站NP和NPL模型滯時為1的月年互相關系數

圖5采用直觀地根據概率密度函數形態(tài)檢驗模擬徑流序列是否服從實測徑流序列分布，根據原序列概率密度分布圖選取具有代表性(如平頂性、多模態(tài)性和偏態(tài)性)的月份展示NP和NPL模型模擬效果。對比發(fā)現NPL模型在再現原序列概率分布的偏態(tài)性時更加貼近于實測序列，盡可能地再現了原徑流序列概率密度分布的不對稱性，如蘭州站的3月、白馬寺站的9月等；針對如蘭州站的8月、龍門站的3月以及白馬寺站的3月等表現的平頂性，NPL模型相比NP模型較好地重現了上述統(tǒng)計特性。但也發(fā)現在NP和NPL模型在模擬一些月份時效果不佳，這可能是因為高斯核函數的光滑性和對稱性而導致的。

圖5 各站NP和NPL模型代表性月份邊緣概率密度函數圖

4 小結

本文以黃河流域蘭州、龍門和白馬寺3個水文站點的天然月徑流序列資料為研究對象，建立了各站點月徑流的NP和NPL模型，對比分析兩模型的適用性，得到主要結論如下。

(1)對于均值、標準差而言，總體來看NP模型模擬效果要優(yōu)于NPL模型，各站點模擬序列的各月統(tǒng)計特征大多數要更接近于實測序列。

(2)對于偏態(tài)系數而言，NPL模型各站11月、12月的模擬效果較差，NP模型各站點的月徑流模擬序列的偏態(tài)系數基本能與實測序列相似。

(3)對于一階和二階自相關系數而言，總體來說對于各站NP和NPL模型二階自相關系數的模擬效果要優(yōu)于一階。

(4)在滯后1個月和2個月的月徑流量與之前12個月的總徑流量之和之間的相關關系方面，滯后1月和滯后2月的結果近似，NPL模型比NP模型較好地再現了這些相關性，表明在NPL模型模擬中，長期依賴性得到了恰當的表示。

(5)根據模擬序列密度函數形態(tài)檢驗徑流序列服從實測徑流序列分布結果，得到NPL模型在再現原序列概率分布的偏態(tài)性和平頂性時更加貼近于實測序列，較好地再現了原徑流序列概率密度分布的不對稱性。