基于LH-OAT的馬斯京根模型參數(shù)全局敏感性分析

劉 春 燁

(太原理工大學(xué)水利科學(xué)與工程學(xué)院，太原 030024)

馬斯京根法是河道洪水流量演算中廣泛應(yīng)用的一種水文學(xué)方法，該法在建立馬斯京根槽蓄方程的基礎(chǔ)上，與水量平衡方程聯(lián)立，推出上下游流量關(guān)系方程，進行河段洪水的演算[1]。馬斯京根法應(yīng)用的一個前提是模型參數(shù)的合理率定，通常對參數(shù)K、x進行率定，再由K、x率定值計算參數(shù)C0、C1、C2。參數(shù)K、x的率定方法有分析法和試算法，但這2種方法較為復(fù)雜且率定結(jié)果精確性較差。因此，部分學(xué)者采用了直接優(yōu)化參數(shù)C0、C1、C2的方法，主要有POS算法[1]、遺傳算法、最小二乘法、SCE-UA算法[2]等。

對模型參數(shù)進行敏感性分析是參數(shù)率定的基礎(chǔ)。敏感性分析即分析參數(shù)對模型的影響程度，根據(jù)參數(shù)的影響程度，篩選出重要的影響參數(shù)，減少參數(shù)的率定次數(shù)和調(diào)參過程中的不確定性，提高率定效率。水文模型敏感性分析方法有局部敏感性分析法、全局抽樣敏感性分析法。局部敏感性分析法原理簡單但無法完全考慮參數(shù)間的相互影響，不能滿足提高模型精確性的要求。全局敏感性分析法可綜合分析多個參數(shù)的變化對模型運行結(jié)果的影響，且該方法已在各種水文模型中得到了應(yīng)用，如RSA方法、GLUE方法[3]和LH-OAT方法[4]等。目前，國內(nèi)外直接對馬斯京根模型參數(shù)C0、C1、C2敏感性分析較少。本研究引入Latin-Hypercube抽樣算法，并結(jié)合OAT 方法,采用LH-OAT全局敏感性分析方法直接對馬斯京根模型參數(shù)C0、C1、C2進行全局敏感性分析。在敏感性分析的基礎(chǔ)上，確定參數(shù)的敏感性大小順序，采用假定最優(yōu)參數(shù)法進行馬斯京根模型參數(shù)率定，發(fā)現(xiàn)按敏感性順序進行參數(shù)的率定，可提高馬斯京根模型參數(shù)的率定效率。

1 LH-OAT全局敏感性分析方法

現(xiàn)階段處理參數(shù)敏感性的方法有隨機法、模糊法、區(qū)間分析法。前2種方法需要較多的數(shù)據(jù)，而區(qū)間估計法只需已知參數(shù)的上下界，方法較為簡單。LH-OAT全局敏感性分析方法屬于區(qū)間估計法，是采用無量綱參數(shù)來表示全局敏感性的一種方法，該法的適用條件為已知模型參數(shù)的取值范圍。LH-OAT全局敏感性分析方法的理論依據(jù)分別為OAT敏感性分析法和LH抽樣法。OAT敏感性分析法在Monte Carlo抽樣方法的基礎(chǔ)上進行重復(fù)試驗，每次對其中一個參數(shù)進行微小的擾動，保持其他參數(shù)不變，帶入敏感性公式，即可得到每個參數(shù)的敏感度，具有較高的精確性。但Monte Carlo抽樣法在參數(shù)取值范圍內(nèi)隨機抽樣，會產(chǎn)生大量的樣本，使OAT敏感性分析方法的運算量變大。LH抽樣方法，假設(shè)抽樣符合均勻分布，將每個參數(shù)區(qū)間均分為n層，分別從每個分層中隨機抽樣一次，生成一個LH抽樣參數(shù)組，減少了采樣頻次，能提高敏感性分析的效率。LH-OAT全局敏感性分析法結(jié)合了LH抽樣法的高效性和OAT算法的精確性。

LH-OAT全局敏感性分析法具體步驟如下：①采用LH抽樣法，將整個參數(shù)空間均分為N層，分別從每個分層中隨機抽樣一次，組成一個抽樣參數(shù)組。②假設(shè)有P個模型參數(shù)，按照OAT方法，對每個LH 抽樣參數(shù)組進行P次參數(shù)的微小擾動并且每次只擾動其中一個參數(shù)，計算每次擾動前后目標(biāo)函數(shù)的變化值。③根據(jù)LH-OAT方法定義，模型總共需運行N(P+1)次。每個參數(shù)的相對敏感度由以下公式計算得到:

(1)

式中：M為敏感性分析的目標(biāo)函數(shù)；ei,k為第i個參數(shù)在第k層中抽樣值；Δei,k為ei,k的擾動變化值；Si,k為參數(shù)ei,k在k層的相對敏感度。

由參數(shù)相對敏感度Si,k可計算參數(shù)的全局敏感度GSi,k：

(2)

式中：N為參數(shù)空間的總層數(shù)；GSi,k為全局敏感度。

2 馬斯京根模型

現(xiàn)階段河道洪水流量演算有水文學(xué)和水力學(xué)2種方法，在水文學(xué)中應(yīng)用最多的為馬斯京根法,其公式為：

Q2=C0I2+C1I1+C2Q1

(3)

式中：Q1為計算時段初始時刻的出流量；Q2為時段末的河段出流量；I1、I2分別為計算時段始、末的河段入流量；C0、C1、C2均為流量演算系數(shù)，計算公式如下：

(4)

(5)

(6)

式中：K為蓄量常數(shù)，即蓄量流量關(guān)系曲線的坡度；x為流量比重系數(shù)；Δt為計算時段長。

馬斯京根法模擬河道流量過程中不可避免地存在著不確定性，對模型參數(shù)進行全局敏感性分析是研究參數(shù)不確定性的分析方法之一。在傳統(tǒng)的馬斯京根模型參數(shù)分析過程中，是通過公式(4)、(5)、(6)來計算C0、C1、C2的值，本文提出直接對C0、C1、C2進行敏感性分析，避免了由K、x取值不準(zhǔn)確對參數(shù)C0、C1、C2造成的影響。

3 實例分析

3.1 模型參數(shù)敏感性計算

為定性了解參數(shù)C0、C1、C2對馬斯京根模型模擬洪水流量的影響程度，提出采用LH-OAT全局敏感性分析方法對模型參數(shù)進行全局敏感性分析。以文獻[5]中數(shù)據(jù)為例，分析C0、C1、C2的敏感性，河段上斷面入流、 下斷面出流情況見表1，所選河道穩(wěn)定，計算時段長為 Δt=6 h。其中時段1-4為連續(xù)時段，時段5-8為連續(xù)時段。由公式(4)、(5)、(6)利用多元函數(shù)求極值的方法得到C0、C1、C2參數(shù)的取值范圍[6]分別為(-9/11，3/13)，(1/4，19/21)，(1/21,1)。

表1 不同時段流量 m3/s

首先，將馬斯京根模型的3個參數(shù)C0、C1、C2在其取值范圍內(nèi)均分成10份，共11層；然后利用Latin-Hypercube抽樣法進行取樣組合，可確定11種組合，且C0、C1、C2應(yīng)滿足公式：C0+C1+C2=1，以避免組合參數(shù)的不合理；按照OAT分析計算方法，對每種組合參數(shù)微小擾動10%(保持2個參數(shù)不變，對另一個參數(shù)進行微小擾動)后利用公式(1)進行計算，每組數(shù)據(jù)需計算11×(3+1)=44次，4個連續(xù)時段可以確定3組計算數(shù)據(jù)，有2組不同的連續(xù)時段，共需計算44×3×2=264次。利用公式(2)對相對敏感性計算結(jié)果進行均值處理，最終得到6組反映3個參數(shù)C0、C1、C2敏感性的無量綱數(shù)值，具體計算結(jié)果見表2。從表1中隨機抽取一組數(shù)據(jù)，分別把參數(shù)C0、C1、C2的擾動幅度變?yōu)?0%、30%，每改變擾動幅度一次需計算44次，變2次擾動幅度，故需計算88次，分析變幅后的敏感性，計算結(jié)果見表3。

表2 參數(shù)的全局敏感度Tab.2 Global sensitivity of parameters

表3 不同擾動幅度參數(shù)敏感度Tab.3 Sensitivity of different disturbance amplitude parameters

3.2 參數(shù)敏感性結(jié)果分析

以徑流量作為敏感性評價目標(biāo)，從模型參數(shù)上來看，1～4時段(即河道漲水階段)馬斯京根模型參數(shù)針對不同流量會產(chǎn)生相應(yīng)波動，但參數(shù)敏感性總趨勢為C0>C1>C2；5-8時段(即河道落水階段)參數(shù)敏感性趨勢為C0VC1>VC2。

圖1 不同時段參數(shù)敏感性比較Fig.1 Comparison of parameter sensitivity

圖2 不同擾動幅度參數(shù)敏感度趨勢Fig.2 Sensitivity trend of different disturbance amplitude parameters

以洪峰流量作為敏感性評價目標(biāo)，當(dāng)出現(xiàn)洪峰時(3～4時段)，C0、C1、C2雖遵循漲水階段的規(guī)律，但C0、C1的敏感性相對前一時段的值減小，C2的敏感性相對增大，3個參數(shù)的敏感性相差不大，無主要敏感性參數(shù)。分析敏感性變化的原因，是由于洪峰流量較河道普通流量大，當(dāng)流量過大時，河道的穩(wěn)定性減弱，影響因素相應(yīng)發(fā)生變化且影響因素增加。

采用假定最優(yōu)參數(shù)法分析參數(shù)敏感性對率定效率和率定結(jié)果的影響，首先選定一組參數(shù)作為最優(yōu)參數(shù)，分別采用2種方法進行參數(shù)率定，第1種方法是根據(jù)參數(shù)的敏感性大小順序進行參數(shù)率定，第2種方法是根據(jù)擾亂參數(shù)敏感性順序進行率定。以率定結(jié)果與最優(yōu)參數(shù)的偏離程度、率定效率來評價2種方法的優(yōu)劣，結(jié)果見表4、圖3。方法1率定11次，精度為99.8%，方法2率定16次，精度為99.6%。因此，按敏感性由大到小進行參數(shù)率定可提高69%的率定效率。以洪峰流量為率定目標(biāo)時，由于沒有主要敏感性參數(shù)，參數(shù)的敏感性差別不大，率定效率會有所下降。

表4 率定結(jié)果 %

圖3 率定結(jié)果比較Fig.3 Comparison of the calibration results

4 結(jié) 語

采用LH-OAT全局敏感性分析方法分析馬斯京根模型參數(shù)的敏感性，并通過實例進行演算，在分析演算過程中得出以下結(jié)論。

(1)直接采用LH-OAT敏感性分析方法對C0、C1、C2進行敏感性分析，確定了各參數(shù)對馬斯京根模型的影響程度，對河道資料較少、無法獲得K、x值的地區(qū)的河道演算具有較好的適用性。

(2)以徑流量為評價目標(biāo)，當(dāng)河道漲水時，參數(shù)C0較敏感，C1、C2一般敏感；當(dāng)河道落水時，C0、C1一般敏感，C2較為敏感，因此，建議在漲水和落水階段分別分析參數(shù)的敏感性。以洪峰流量作為評價目標(biāo)時，C0、C1的敏感性相對減弱，C2的敏感性相對增加，但3個參數(shù)敏感性差別較小，無主要敏感性參數(shù)。

(3)按敏感性大小順序進行參數(shù)率定，可以大幅度提高馬斯京根模型參數(shù)的率定效率。在實際工程計算時需針對研究河道實際情況研究參數(shù)C0、C1、C2的敏感性，以期最接近河道流量的真實情況，敏感性分析成果為馬斯京根模型參數(shù)的合理率定奠定了基礎(chǔ)。

□

[1] 李 匡,胡宇豐,梁犁麗，等.PSO 算法在馬斯京根法參數(shù)率定中的應(yīng)用[J].水電站機電技術(shù),2016,39(8):88-91.

[2] 李鴻雁,趙 娟,王玉新，等.擴域搜索遺傳算法優(yōu)化馬斯京根參數(shù)及其應(yīng)用[J].吉林大學(xué)學(xué)報,2011,41(3):861-865.

[3] 李 丹,張 翔,張 揚.水文模型參數(shù)敏感性的區(qū)間分析[J].水利水電科技進展,2011,31(1):29-32,41.

[4] 劉慧如,吳建華,孫 毅，等.基于LH-OAT的雙曲正切模型參數(shù)敏感性分析[J].水電能源科學(xué),2016,34(9):84-86.

[5] 林三益.水文預(yù)報[M]. 北京：中國水利水電出版社,2001.

[6] 李 匡,付 力,胡宇豐，等.馬斯京根法參數(shù)C0、C1、C2取值范圍的確定[J].南水北調(diào)與水利科技,2012,10(5):43-45,55.