基于不同目標函數(shù)GLUE方法參數(shù)不確定性研究

溫婭惠，李致家，2，霍文博，張漢辰，童冰星

(1.河海大學水文水資源學院，江蘇 南京 210098；2.河海大學水安全與水科學協(xié)同創(chuàng)新中心，江蘇 南京 210098)

0 引 言

水文系統(tǒng)的不確定性主要體現(xiàn)在自然現(xiàn)象、水文數(shù)據(jù)、模型結(jié)構(gòu)和模型參數(shù)4個方面。近年來，水文模型參數(shù)的不確定性問題得到了廣泛的研究，1992年Beven和Binley提出普適似然不確定性估計方法(Generalized Likelihood Uncertainty Estimation，GLUE)。GLUE方法認為：模型最后模擬結(jié)果的好壞由模型參數(shù)的綜合決定，而不是某一個模型參數(shù)所致[1]。1996年，F(xiàn)ree等將GLUE方法應用于TOPMODEL模型預報中[2]；1997年，F(xiàn)ranks等將GLUE方法應用于SWAT模型中[3]；2009年，衛(wèi)曉婧等人結(jié)合MCMC方法，對GLUE方法進行改進[4]；2010年，劉艷麗等人基于GLUE方法體系分析了大伙房模型參數(shù)的不確定性[5]。雖然GLUE方法在流域水文模型的不確定性分析中得到廣泛應用，但是在使用GLUE方法優(yōu)選參數(shù)時，對于選擇單目標似然函數(shù)(確定性系數(shù))的研究較多，而對于多目標似然函數(shù)的研究還需進一步發(fā)展。GLUE方法選取確定性系數(shù)為似然目標時，確定性系數(shù)只能描述水文過程的整體擬合特征，會使模擬誤差進行累積均化，造成模擬精度虛假增高，卻忽略了模型的其他判斷標準。比如說，洪峰誤差和峰現(xiàn)時間誤差，這些特征對于洪水預報都具有重要的參考價值及意義。新安江模型作為我國最重要的流域水文模型之一，已經(jīng)在洪水預報方面發(fā)揮了非常大的作用[6]。傳統(tǒng)預報過程中，不管通過人工調(diào)參或計算機自動率定，都是得出一組最優(yōu)參數(shù)，得到單一預報結(jié)果；又由于新安江模型參數(shù)之間具有較強的互補性，而新安江模型內(nèi)部的產(chǎn)匯流參數(shù)又有地區(qū)性的特征，因而給預報結(jié)果帶來很大的不確定性[7]。因此，本文以半濕潤的陳河流域研究區(qū)域，分別運用以確定性系數(shù)為目標函數(shù)的單目標GLUE方法和以綜合確定性系數(shù)、洪峰合格數(shù)、峰現(xiàn)時間合格數(shù)為似然目標函數(shù)的多目標GLUE方法，研究新安江模型參數(shù)的不確定性和敏感性以及對兩種方法篩選出的參數(shù)進行對比。

1 GLUE方法介紹

1.1 GLUE方法原理

GLUE方法基于Bayesian方程和Monte-Carlo隨機取樣方法。其在預先設定的參數(shù)分布取值空間內(nèi)，利用 Monte-Carlo隨機采樣方法獲取模型的參數(shù)值組合，運行模型。選定似然目標函數(shù)，計算模型模擬結(jié)果與觀測值之間的似然函數(shù)值，在所有的似然值中，設定一個臨界值，對高于臨界值的所有參數(shù)組似然值重新歸一化，按照似然值的大小，求出在某置信度下模型預報的不確定性范圍[8]。

1.2 GLUE方法步驟

1.2.1 定義似然判據(jù)

本文采用單目標GLUE方法和多目標GLUE方法分析研究。單目標GLUE方法：選擇用于反映模擬結(jié)果與實測結(jié)果吻合程度的確定性系數(shù)(Nash系數(shù))作為似然判據(jù)。即

(1)

多目標GLUE方法：本文采用綜合確定性系數(shù)、洪峰合格數(shù)、峰現(xiàn)時間合格數(shù)三個目標函數(shù)。多目標似然函數(shù)為

L(M(θI)/Y)=L(L1(θI)，L2(θI)，L3(θI))

(2)

其中：

(1)綜合確定性系數(shù)

(3)

式中，DYi為第i場洪水的確定性系數(shù)；α為平衡因子，表示模擬效果差的洪水對目標函數(shù)的影響，本文取4；n為洪水場次數(shù)。

(2)洪峰合格數(shù)

(4)

式中，l2i為第i場洪水是否合格，合格指洪峰誤差控制在20%之內(nèi)，記為1；不合格指洪峰誤差絕對值大于20%，記為0；n為洪水場次數(shù)。

(3)峰現(xiàn)時間誤差合格數(shù)

(5)

式中，l3i為第i場洪水的峰現(xiàn)時間誤差是否合格，合格，則峰現(xiàn)時間誤差小于3 h，記為1；不合格，則峰現(xiàn)時間誤差大于3 h，記為0；n為洪水場次數(shù)。

1.2.2 參數(shù)先驗分布和取值范圍確定

在新安江模型中，選擇蒸散發(fā)折算系數(shù)K、自由水蓄水容量SM、地下水出流系數(shù)KG、壤中流出水系數(shù)KI、地下水消退系數(shù)CG、壤中流消退系數(shù)CI以及河網(wǎng)水流消退系數(shù)CS這7個較為敏感的參數(shù)來研究其對模型模擬結(jié)果的不確定性。由于對參數(shù)缺乏足夠的信息和認知的局限性，因此假設參數(shù)的先驗分布服從均勻分布[9-11]，參數(shù)的取值范圍見表1。

表1 新安江次洪模型敏感參數(shù)取值范圍

1.2.3 計算過程

(1)根據(jù)Monte-Carlo隨機取樣方法，對以上的7個參數(shù)的取值范圍內(nèi)隨機取樣，生成10 000組參數(shù)組合。

(2)將10 000組參數(shù)依次帶入到流域的次洪模型中，計算預報流量與實際流量過程的確定性系數(shù)，統(tǒng)計洪峰合格數(shù)，峰現(xiàn)時間誤差合格數(shù)，根據(jù)單目標函數(shù)和多目標函數(shù)篩選原則，篩選出符合要求的參數(shù)組合。

(3)將篩選出的有效參數(shù)組合，帶入到需要驗證的洪水過程中。最后，將符合要求似然值的模擬流量按大小排序，估算出模型預報的置信區(qū)間為90%的不確定性的洪水序列。

(4)似然函數(shù)的更新，當有新的實測洪水數(shù)據(jù)加入計算時，利用之前模擬結(jié)果得到的參數(shù)組作為本次模擬的先驗分布，重新調(diào)整參數(shù)后驗概率分布和不確定性范圍。

2 流域介紹與模型率定

陳河流域位于北溫帶，屬于大陸性季風氣候[12]。該流域面積為1 350 km2，河長126 km，多年平均降水量700～900 mm，河流水量主要系雨水補給[13]，發(fā)生洪水的主要原因是局部暴雨[14-15]。流域平均徑流深度100～500 mm，徑流系數(shù)0.2～0.5，為相對較高產(chǎn)流區(qū)。陳河流域上所轄雨量站有兩河岔、厚畛子、沙梁子、老水磨、金井、南留、釣魚臺、麥場、板房子、小王澗(見圖1)。

本文選擇研究的參數(shù)有K、SM、KG、KI、CG、CI、CS這7個較為敏感的參數(shù)，其他不敏感參數(shù)數(shù)值(見表2)采用新安江模型的率定值。新安江模型參數(shù)率定采用人工優(yōu)選法和SCE-UA自動優(yōu)選法結(jié)合，計算時段長為1 h。

表2 陳河流域非敏感參數(shù)取值

3 參數(shù)敏感性分析

確定性系數(shù)反映洪水預報過程與實測過程之間的吻合程度，反映總體模擬的擬合好壞；因此，用確定性系數(shù)繪制參數(shù)與似然值散點分布圖，以陳河流域2003090319號洪水計算結(jié)果為例，繪制似然函數(shù)與散點分布圖(圖略)。散點圖反映了各個參數(shù)對模型不確定性的影響程度，體現(xiàn)了各個參數(shù)的敏感性。

參數(shù)似然值散點分布圖清晰地反映了參數(shù)值的選取與似然值之間的關(guān)系。從散點圖分析得，參數(shù)K變化平緩且分布比較均勻，說明該參數(shù)為非敏感參數(shù)，參數(shù)K在新安江日模型中為敏感參數(shù)，在次洪模型中為非敏感參數(shù)；參數(shù)SM在(50，200)的范圍內(nèi)有變化，參數(shù)SM數(shù)值越大，似然值呈偏小趨勢，但變化幅度不大，說明參數(shù)SM對不確定性有影響但在一定范圍內(nèi)影響不大；參數(shù)KG、KI、CG變化平緩，表明他們對不確定性影響不大；參數(shù)CI變化也較為平緩，但隨著參數(shù)值的增大，似然值呈現(xiàn)微弱幅度的增大；參數(shù)CS存在明顯的高峰區(qū)域，在(0，0.2)范圍內(nèi)變化較為平緩，在(0.2，3)范圍內(nèi)變化十分明顯，隨著參數(shù)CS數(shù)值的增大，似然值變小，說明CS對不確定性的影響比較顯著。

綜上所述，將研究的7個參數(shù)分為敏感參數(shù)(SM、CS)和非敏感參數(shù)(K、KG、KI、CG、CI)。

4 計算結(jié)果分析與比較

選擇陳河流域2003年～2012年的18場洪水對新安江模型參數(shù)進行不確定性研究。其中，前14場對參數(shù)進行先驗分布分析；后4場洪水進行后驗分布檢驗。將隨機采樣得到的10 000組參數(shù)分別代入陳河流域的每一場洪水中進行計算，得到10 000組模擬流量過程和目標函數(shù)值。結(jié)果發(fā)現(xiàn)，與實測流量相比，有些參數(shù)組合模擬的較好，有的較差；但流量過程線整體趨勢相同。而且，不同的參數(shù)組合可以得到相同的模擬過程。這恰恰印證了GLUE方法的中心思想，模型的結(jié)果由參數(shù)組合決定。

4.1 單目標GLUE方法篩選結(jié)果

根據(jù)單目標似然函數(shù)標準進行篩選，設定似然臨界值為0.7，不斷更新似然函數(shù)，參數(shù)組數(shù)不斷減少。在經(jīng)過14場洪水演算后篩選得到14個參數(shù)組，這些參數(shù)組滿足了使所有場次洪水的確定性系數(shù)均大于0.7，為最優(yōu)模型參數(shù)解的集合。單目標似然函數(shù)篩選結(jié)果(見表3)。其中，序號中的數(shù)值表示篩選的參數(shù)組合在隨機生成的10 000萬參數(shù)中的排列序號。

表3 單目標似然函數(shù)篩選參數(shù)組

圖2 單目標GLUE方法后四4場洪水后驗分布得到的不確定性范圍

將篩選出的14組參數(shù)對沒有用來更新似然分布的后4場洪水模型進行模擬，將模擬流量按大小排序，估算出90%置信度下新安江模型模擬的不確定性范圍。即，用累計似然分布的5%和95%兩個分位點作為不確定的上限和下限。陳河流域單目標GLUE方法后4場洪水檢驗的不確定性范圍(見圖2)，包括實測值，面平均雨量和90%置信區(qū)間范圍。

4.2 多目標GLUE方法篩選結(jié)果

根據(jù)多目標似然函數(shù)標準進行篩選，滿足綜合確定性系數(shù)L1(θi)≥(改為：L1(θi)≥；下同)0.8的有297組參數(shù)，滿足洪峰合格數(shù)L2(θi)≥11的有54組參數(shù)，滿足峰現(xiàn)時間合格數(shù)L3(θi)≥11的有82組參數(shù)。經(jīng)過篩選，同時滿足這三個似然函數(shù)的參數(shù)有32組次。表4列出了陳河流域多目標似然函數(shù)篩選結(jié)果。

表4 多目標似然函數(shù)篩選參數(shù)組

將篩選出的32組參數(shù)對沒有用來更新似然分布的后4場洪水模型進行模擬。即，陳河流域多目標GLUE方法后4場洪水后驗分布得到的不確定性范圍(見圖3)。

4.3 分析與比較

將單目標GLUE方法篩選的14組參數(shù)和多目標GLUE方法篩選的32組參數(shù)進行總和標準化計算，計算公式

(6)

式中，Xij′為總和標準化后的新值；Xij為原始值；i為參數(shù)類別；j為參數(shù)組號；n為參數(shù)組數(shù)。

總和標準化圖如圖4所示。

由圖4和表5中可知，參數(shù)K、KG、KI、CG、CI取值較為分散，參數(shù)SM和CS取值較為集中，同類參數(shù)在等效的似然值中的取值不同。這說明模型中每類參數(shù)不確定性程度和敏感度存在一定的差異。使用多目標GLUE方法篩選的參數(shù)范圍都大于單目標GLUE方法篩選的參數(shù)，非敏感參數(shù)最后固定的取值范圍與本文設定的先驗分布的取值范圍中基本不變，而敏感參數(shù)SM和CS取值相較于本文設定的參數(shù)先驗取值范圍，最后篩選出的數(shù)值基本上集中在一個小的范圍內(nèi)。值得注意的是，其中敏感參數(shù)CS基本在0.01～0.17之間，并沒有像其他敏感參數(shù)分布那樣差距很大；說明CS在模型計算有導向作用。這應該與模型以洪峰和徑流深合格為主要評價標準有關(guān)，也從側(cè)面反映了CS在對洪峰有主導控制作用。

圖3 多目標GLUE方法4場洪水后驗分布得到的不確定性范圍

圖4 參數(shù)標準化

單目標GLUE方法和多目標GLUE方法篩選出的參數(shù)對后四場洪水進行模擬分析，由圖2和圖3可知，在高流量區(qū)不確定性大，在低流量區(qū)不確定性小。統(tǒng)計檢驗洪水實測流量過程在90%置信區(qū)間所占的比例，即覆蓋率，計算結(jié)果見表6。由表6可知，用于檢驗的四場洪水，多目標GLUE方法得到的不確定性范圍相對于單目標GLUE方法得到的不確定性范圍大，且區(qū)間覆蓋率高，即實測流量被包入置信區(qū)間的范圍大，預報精度更高，得到的概率預報更具有實踐意義。

采用單目標GLUE方法共篩選出14組參數(shù)，采用多目標GLUE方法共篩選出32組參數(shù)。其中，兩種方法篩選出的重復參數(shù)組合有7組，參數(shù)序號分別為1 584，4 442，6 644，6 886，7 106，7 563，7 989。單目標GLUE方法選擇確定性系數(shù)為似然目標函數(shù)，選定0.7為似然值。即，滿足確定性系數(shù)大于等于0.7的參數(shù)組合才能被篩選出來，確定性系數(shù)表示預報過程與實測過程之間的吻合程度，會使模擬誤差進行累積均化。另外，對于一個流域，會存在模擬效果較差的洪水，即不同的參數(shù)組合都不能使該場洪水確定性系數(shù)大于等于0.7，以確定性系數(shù)為篩選標準時，則會使篩選出的參數(shù)大大減少或者沒有。相對于單目標GLUE方法，多目標GLUE方法在選擇似然函數(shù)時進行了優(yōu)化，采用綜合確定性系數(shù)、洪峰合格數(shù)、峰現(xiàn)時間合格數(shù)為似然目標函數(shù)，而且綜合確定性系數(shù)會通過平衡因子調(diào)解模擬效果較差的洪水，平衡因子取值越大，模擬效果差的洪水對綜合確定性系數(shù)影響更大；另外，洪峰相對誤差和峰現(xiàn)時間誤差也是半濕潤地區(qū)非常重要的參考指標。綜上，使用多目標GLUE方法篩選參數(shù)更合理。

表5 單目標方法和多目標方法篩選參數(shù)最大值和最小值

表6 兩種方法區(qū)間覆蓋率計算結(jié)果 %

5 結(jié) 論

以半濕潤的陳河流域為研究流域，運用基于單目標GLUE方法和多目標GLUE方法對新安江模型的參數(shù)進行的對比分析表明：

(1)在半濕潤的陳河流域，可將所研究的7個參數(shù)分為不敏感參數(shù)(KG、KI、CG、CI、CS)、敏感參數(shù)(K、SM)兩類。另外，以確定性系數(shù)的單目標GLUE方法和以綜合確定性系數(shù)、洪峰合格數(shù)、峰現(xiàn)時間合格數(shù)的多目標GLUE方法對參數(shù)不確定性進行對比分析討論的結(jié)果表明，多目標GLUE方法比單目標GLUE方法得到的區(qū)間覆蓋率更大，預報精度更高，有利于下一步的參數(shù)以及預報不確定性范圍的分析。

(2)實測流量并不能完全包含在90%置信度的流量界限內(nèi)，說明參數(shù)的取值范圍還不能覆蓋所有的值， GLUE不確定性分析方法還需在實踐中進一步完善和推廣；另外，還可將該方法與模糊數(shù)學結(jié)合，并對模型輸入不確定性和模型結(jié)構(gòu)不確定性進行研究，構(gòu)建更全面的不確定性評價體系。