馬斯京根模型參數反演的改進粒子群算法

張新明，馬艷

(哈爾濱工業大學 深圳研究生院，廣東 深圳 518055)

?

馬斯京根模型參數反演的改進粒子群算法

張新明，馬艷

(哈爾濱工業大學 深圳研究生院，廣東 深圳 518055)

摘要：針對傳統粒子群算法容易陷入局部最優(即“早熟”現象)的問題，將基于適應值共享原則的小生境策略與粒子群算法相結合，提出了一種改進的粒子群算法——小生境粒子群算法，并將之應用于4個典型測試函數的數值仿真以及基于馬斯京根模型的參數反演計算。數值模擬結果顯示，相比于傳統的粒子群算法，小生境粒子群算法具有精度高、收斂速度快的特點，但其抗噪性較差。為了進一步提高算法的抗噪性，將基于小波多分辨分析的多尺度反演策略和小生境粒子群算法相結合構造了多尺度小生境粒子群算法。帶有5%隨機噪聲的馬斯京根模型參數反演結果顯示，新提出的多尺度小生境粒子群算法能夠有效提升小生境粒子群算法的抗噪性，從而使反演結果的精度得到較大的改善。

關鍵詞：馬斯京根模型；參數反演；小生境粒子群算法；多尺度；隨機噪聲

眾所周知，由于其簡單實用性，馬斯京根模型是眾多洪水演算模型中應用最為廣泛的方法之一，它利用河段水量槽蓄方程代替復雜的水動力方程從而使計算過程極大簡化, 同時又能取得滿足實用的演算精度, 從而被國內外廣泛應用, 在洪水預報和防洪規劃中具有重要意義。然而，該方法在實際應用中的一個重要難點即是模型的參數優選問題。Gill[1]采用最小二乘方法求解非線性馬斯京根模型的三個參數值；Tung[2]將Hook-Jeeves模式搜索方法分別和線性回歸方法、共軛梯度方法、DFP方法相結合應用于馬斯京根模型的參數識別問題，得到了較好的結果；Yoon等[3]討論了馬斯京根模型參數估計的多種方法，如：基于最小二乘方法的線性及非線性回歸法、非線性迭代方法、線性規劃、二次規劃等。但由于馬斯京根模型本身的近似性和傳統方法的局限性, 上述方法較難得到全局最優解。而近些年來，隨著各種智能優化算法的發展，其在馬斯京根模型參數優化求解問題中的應用逐漸引起了眾多專家學者的關注。Mohan[4]采用遺傳算法研究了馬斯京根模型參數估計問題，結果顯示方法的全局尋優性明顯優于上述三種方法，并且不需要嚴苛的初始值猜測；Kim等[5]和程銀才等[6]分別將和聲搜索方法和混沌模擬退火方法用于相同問題的求解，確定參數x、n和K，有效地提高了收斂速度和模型精度；袁曉輝等[7]和魯帆等[8]針對傳統遺傳算法的不足，分別提出了相應的改進遺傳算法，用于非線性馬斯京根模型參數估計，得到了較好的結果；馬細霞等[9]和Chu等[10]分別提出采用粒子群算法對馬斯京根模型參數進行率定，研究結果顯示了粒子群算法較好的全局尋優能力。

粒子群優化算法[11]( particle swarm optimization, PSO )是由Kennedy 和Eberhart等提出的一種基于種群搜索的自適應進化計算技術。該算法具有并行處理和魯棒性好等特性。它不依賴于問題的具體領域, 以粒子群個體作為運算對象，直接以目標函數作為尋優搜索的基本信息，可以使用整個種群的信息，并且占用計算機內存少，尤其適用于求解一些非線性、多參數復雜系統的全局優化問題。但在實際應用中，其存在兩個主要的缺點，分別是容易陷入局部收斂和后期收斂速度慢。為了能夠有效的解決上述問題，本文提出將小生境思想[12]和粒子群算法結合構造小生境粒子群算法，并將其應用于馬斯京根模型的參數反演，同時為了提高方法的抗噪性，在反演過程中，加入了多尺度思想[13]，有效提高了反演結果的計算精度，改善了方法的抗噪性。

1馬斯京根模型

馬斯京根法是河道洪水演算中廣泛應用的方法, 其采用的基本方程如下：

(1)

式中： W為槽蓄量，t為時間，I、Q為入流量、出流量，K為槽蓄系數，x為流量比重因子。

該模型的解析解[8]為

Q2=C1I1+C2I2+C3Q1

(2)

式(1)中K值基本上反映的是河道穩定流的傳播時間，理論上應隨流量的增大而減小，不少實測資料也是這樣的；x值反映的是楔蓄在河槽調蓄中的影響，對于某河段，x在洪水漲落的過程中基本穩定；Δt的選取要求摘錄的洪水數值能比較真實地反映洪水的變化情況并且不要溜掉洪峰。

對于一個河段，只要確定參數K，x的值和演算時段Δt之后，就可以求出C1、C2、C3，代入式(2)，再根據上站流量過程和下站起始流量從而計算出下站的流量過程。

當Δt=2Kx時則C2=0，則式(2)就變成：

Q2=C1I1+C3Q1

(3)

式(3)計算簡便，又能獲得Δt的預見期。

2小生境粒子群算法

小生境是來自于生物學的一個概念，它是指在特定環境下的一種生存環境，即生物居住生活的小范圍或小棲息地。而應用于進化算法中的小生境技術就是把上述概念應用到進化算法中，產生小生境機制。主要的小生境機制有預選擇機制(preselection)、排擠機制(crowding)與共享機制(sharing)。小生境技術具有較強的局部搜索能力，能夠保持物種的多樣性，防止過早收斂。將其與粒子群算法相結合，構造小生境粒子群算法，能夠有效的促進粒子群算法跳出局部最優找到全局最優。 本文采用的是基于適應值共享的小生境技術[15]。

2.1基于適應值共享的小生境粒子群算法

小生境共享機制的基本思想為：利用共享函數去判斷個體之間相似程度，對適應值進行調整。也就是當適應度值減小時，代表一個個體與其他個體比較相似,反之，當適應度值增大時，代表一個個體與其他個體相似程度較差。按照上述方式進行就可以有效的控制相似個體復制過多,從而形成一種較好的小生境進化環境。

本文的適應度計算基于適應值共享原則，要想實現適應值共享，首先要定義一種距離度量方式。常用的距離度量方式有兩種：一種是在參數空間中被廣泛采用的歐幾里得距離，它也被稱作是表現型距離；另一種是在編碼空間中被廣泛采用的基因距離，它也被稱作是海明距離。采用第一種距離方式。

將粒子i與粒子j之間的距離記為符號dij，共享函數的表達式如下：

(4)

式中：α用來調整共享函數的形狀的常數，本文中選取α=1；σ為預先給定的小生境半徑。小生境半徑計算公式為[14]

(5)

因此，依據上述處理，假如在一個小生境中存在非常多的個體，那么在該小生境中基于適應值共享后的所有個體的適應值會大大的降低，這樣就會讓那些存在較少個體的小生境能夠存在并繁衍。

基于適應值共享原則的小生境粒子群算法優化的基本步驟為：

2)確定小生境種群個體。

首先，令i=1；

其次，計算兩個粒子個體的距離dij；

最后，讓上述距離dij小于小生境半徑σ，進而確定小生境群體；

3)按照粒子群算法對小生境群體進行速度和適應度更新，再對更新后的粒子更新其適應度；

4)計算每個粒子的適應值，保留最優的適應值及個體，檢查是否達到優化條件； 如果達到，則結束。否則，進入下一個粒子的小生境群體進行優化；

5)若沒有找到最優值， 則對每個粒子的小生境群體保留的最優個體組成新的群體空間，重復步驟2)～4)。

2.2多模態函數優化

為驗證小生境粒子群算法的有效性，選取表1所示的四個典型多模態函數作為測試函數，其中，F1為Rastrigrin 函數、F2為Griewangk's函數， F3為Branin函數、F4為Schaffer 函數，它們在自己所在的可行解范圍內的全局最優值均為0。

選取的參數設定如下：迭代次數為100，多模態函數優化時的參數范圍如表1所示。對上述函數進行優化時，在基于適應值共享原則的小生境粒子群算法(NPSO)中，選取粒子個數為30，粒子維數為2，小生境個數為4，小生境半徑可由式(5)計算求得。同時，采用基本粒子群算法(PSO)對上述4個函數也進行了計算，選取參數和小生境粒子群算法(NPSO)完全一致。為克服算法隨機性的影響，所有計算結果都是進行30次計算后的平均結果。迭代尋優結果的對比如表2和圖1所示。

從表2和圖1中，可以清晰地看到對于四個測試函數，NPSO在解的精確性上都明顯高于PSO，而在迭代次數和收斂時間上，前三個測試函數有明顯的提升，而對于函數F4則相差不大。總的來說，相比于基本粒子群算法，小生境粒子群算法在全局尋優方面更具有優勢。

表1　測試函數

表2　NPSO算法與PSO算法結果

圖1　NPSO與PSO的迭代尋優結果Fig.1　The iterative results of algorithm NPSO and algorithm PSO

3馬斯京根模型參數反演

首先給出馬斯京根模型參數反演的小生境粒子群算法的基本步驟:

算法1：基于適應值共享原則的小生境粒子群算法

2)利用模型(1)，計算該河段的流量數據并記為Q2(Δti,x,K)，其中Δti為i個不同的時段。

3)這樣便得到了不同時段上Q2的計算值與實測值，將他們作差后平方相加得出誤差的平方和，形式如下：

4)把J(x,K)作為小生境粒子群算法的適應值函數，適時調整參數x,K，使式(6)達到最小，即可得到參數的全局最優值。

為了檢驗上述方法的可行性，選取以下算例進行驗證。

3.1單參數反演

依據上述模型分別對槽蓄系數K、流量比重因子x進行了反演，其真值為K=18，x=0.1。

為了進一步提高算法的抗噪性，將多尺度反演策略和小生境粒子群算法相結合，提出了多尺度小生境粒子群算法。其主要求解思路如下：

算法2：多尺度小生境粒子群算法

2)設定待反演參數初始選擇區間(15≤K≤25和0≤x≤0.5)，以最粗糙尺度為起始尺度，采用小生境粒子群算法求解，目標函數為

3)將尺度加細，依據上一尺度的反演結果來確定下一尺度的待反演參數區間。具體來說即是：采用二分法將上一尺度的待反演參數區間劃分為兩等份，判斷上一尺度的反演結果落在哪一個子區間，以包含反演結果的子區間作為下一尺度的待反演參數區間；

4)在較細尺度下，采用小生境粒子群算法求解，目標函數為

5)重復上述步驟3)和4)直至最精確的尺度，得到最終反演結果。

在加噪5%的情形下，將多尺度小生境粒子群算法應用于馬斯京根模型的參數反演。尺度分解工具選用Matlab小波工具箱中的相關命令(尺度分解函數：wavedec，系數提取函數：appcoef，信號重構函數：wrcoef)，分解及重構基函數選用Daubechies小波，本文中分解為5個尺度。參數反演結果如表6和表7所示。表6是槽蓄系數K和流量比重因子x的最終反演結果。從表6中可以看到兩個參數的計算精度都有了較好的提升，和真值的相對誤差分別由21%提升到了4%和由36%提升到了接近10%。表7給出了在各個不同尺度下的反演結果。從表7中可以看到隨著尺度的減小，反演結果的精度在不斷的改善。但需要引起注意的是，在槽蓄系數K的多尺度反演過程中，最好的反演結果并沒有出現在原始尺度，而是出現在了第四層，因此在多尺度算法的實現過程中，要留意保留不同尺度下的數值計算結果，并通過最終比較給出最好的反演值。

表3　Q2的實測值

Δt(h)7.5.27.5.87.5.147.5.207.6.27.6.87.6.147.6.20I1/(m3·s-1)54153.4365457015.3265861235.7854264256.2326367365.3025470124.3265474136.6984177562.36254I2/(m3·s-1)54296.5894656324.1245855216.2165459126.0254162154.6540265847.3104270125.3698575584.26589Q1/(m3·s-1)46326.20314549325.6548754897.2364558125.6958762145.3542166859.0324870926.6529874326.56984Q(m)2i=Q2/(m3·s-1)51957.9849854620.7912357704.8946161038.8772164385.1320167969.3807772073.4487076072.67987

表4　參數反演結果(5≤K≤25,0≤x≤0.5)

表5　加噪5%參數反演結果(5≤K≤25,0≤x≤0.5)

表6　參數反演結果(多尺度小生境粒子群算法，加噪5%，15≤K≤25,0≤x≤0.5)

表7　不同尺度下參數K,x反演結果(加噪5%)

3.2雙參數反演

算法參數設定如下：粒子數量為30，迭代次數為100次，小生境個數為4個，粒子維數為2。同樣的，為克服算法隨機性的影響，所有計算結果也都是進行30次計算后的平均結果。同時對參數K和x進行反演，待反演參數初始范圍分別是5≤K≤25和0≤x≤0.2，反演結果如表8所示，從表中可以看到反演結果精度較高，相對誤差可以控制在4%之內。表9給出了在對實測值添加5%的噪聲情形下，采用小生境粒子群算法得到的反演結果，結果顯示反演的精度明顯降低，參數K的相對誤差達到了5%以上，參數x的相對誤差達到18%以上。表10給出了結合了多尺度策略后的反演結果，兩個參數反演值的精度都有所提高，和真值的相對誤差分別達到了3.08%和13.26%。反演結果顯示相比于單參數反演，雙參數反演的多尺度小生境粒子群算法的抗噪性有所欠缺，這還有待于進一步研究。

表8　馬斯京根模型雙參數反演結果

表9　雙參數反演結果(加噪5%)

表10　雙參數反演結果(多尺度小生境粒子群算法，加噪5%)

4結論

1)本文將基于適應值共享原則的小生境策略與傳統粒子群算法相結合，提出了一種改進粒子群算法—小生境粒子群算法，并應用于基于馬斯京根模型的參數反演中。通過對槽蓄系數和流量比重因子進行單參數和雙參數反演計算實驗，結果顯示：如果不考慮噪聲的話，小生境粒子群算法具有理想的計算精度和收斂速度。單參數的反演值和真值的相對誤差可以達到10-6，和真值幾乎完全吻合；雙參數的反演結果的相對誤差也可以控制在4%之內。但是，該方法的抗噪性卻是不甚理想。

2)為進一步提高算法的抗噪性，將基于小波多分辨分析的多尺度反演策略和小生境粒子群算法相結合構造了多尺度小生境粒子群算法，并將之應用于帶有5%隨機噪聲的馬斯京根模型參數反演，使反演結果的精度得到了明顯的改善，有效的提升了小生境粒子群算法的抗噪性。

3)數值模擬實驗顯示提出的小生境粒子群算法能夠有效的改善傳統粒子群算法的全局尋優性，而多尺度反演策略的加入能夠有效的改進算法的抗噪性。因此在進行實際河道洪水演進計算時，兩者結合所構造的多尺度小生境粒子群算法能夠為馬斯京根模型的參數估計問題提供一種高效的計算方法。

參考文獻：

[1]GILL M A. Flood routing by the Muskingum method[J]. Journal of Hydrology, 1978, 36(3/4): 353-363.

[2]TUNG Y K. River flood routing by nonlinear Muskingum method[J]. Journal of hydraulic engineering, 1985, 111(12): 1447-1460.

[3]YOON J, PADMANABHAN G. Parameter estimation of linear and nonlinear Muskingum models[J]. Journal of water resources planning and management, 1993, 119(5): 600-610.

[4]MOHAN S. Parameter estimation of nonlinear Muskingum models using genetic algorithm[J]. Journal of hydraulic engineering, 1997, 123(2): 137-142.

[5]KIM J H, GEEM Z W, KIM E S. Parameter estimation of the nonlinear Muskingum model using harmony search[J]. Journal of the American Water Resources Association, 2001, 37(5): 1131-1138.

[6]程銀才, 李明華, 范世香. 非線性馬斯京根模型參數優化的混沌模擬退火法[J]. 水電能源科學, 2007, 25(1): 30-33.

CHENG Yincai, LI Minghua, Fan Shixiang. Application of Chaotic simulated annealing algorithm to parameter optimization of nonlinear Muskingum model[J]. Water resources and power, 2007, 25(1): 30-33.

[7]袁曉輝, 張雙全, 張勇傳, 等. 非線性馬斯京根模型參數率定的新方法[J] . 水利學報, 2001, 32(5): 77-81.

YUAN Xiaohui, ZHANG Shuangquan, ZHANG Yongchuan, et al. Parameter estimation of nonlinear Muskingum Model using mixed genetic algorithm[J]. Journal of hydraulic engineering, 2001, 32(5): 77-81.

[8]魯帆, 蔣云鐘, 王浩, 等. 多智能體遺傳算法用于馬斯京根模型參數估計[J]. 水利學報, 2007, 38(3): 289-294.

LU Fan, JIANG Yunzhong, WANG Hao, et al. Application of multi-agent genetic algorithm to parameter estimation of Muskingum model[J]. Journal of hydraulic engineering, 2007, 38(3): 289-294.

[9]馬細霞, 舒丹丹, 黃渝桂. 基于PSO的非線性馬斯京根模型參數率定新方法[J]. 鄭州大學學報: 工學版, 2007, 28(4): 122-125.

MA Xixia, SHU Dandan, HUANG Yugui. Parameter estimation method of nonlinear Muskingum model based on PSO[J]. Journal of Zhengzhou University: engineering science, 2007, 28(4): 122-125.

[10]CHU H J, CHANG L C. Applying particle swarm optimization to parameter estimation of the nonlinear Muskingum model[J]. Journal of hydrologic engineering, 2009, 14(9): 1024-1027.

[11]KENNEDY J, EBERHART R. Particle swarm optimization[C]//Proceedings of the IEEE International Conference on Neural Networks. Perth, WA: IEEE, 1995: 1942-1948.

[12]HUTCHINSON G E. Concluding remarks[J]. Cold spring harbor symposium on quantitative biology, 1957, 22: 415-427.

[13]ZHANG Xinming, ZHOU Chaoying, LIU Jiaqi, et al. Multiparameter identification of fluid-saturated porous medium with the wavelet multiscale method[J]. Journal of porous media, 2009, 12(3): 255-264.

[14]GOLDBERG D E, RICHARDSON J. Genetic Algorithms with sharing for multimodal function optimization[C]//Proceedings of the 2nd International Conference on Genetic Algorithms on Genetic Algorithms and Their Application. Hillsdale: Lawrence Erlbaum Associates, 1987: 41-47.

Improved particle swarm optimization for parameter inversion of Muskingum model

ZHANG Xinming, MA Yan

(Shenzhen Graduate School, Harbin Institute of Technology, Shenzhen 518055, China)

Abstract：In this study, the parameter inversion problem of the Muskingum model is considered. To overcome the "premature" phenomenon of particle swarm optimization, a new niche particle swarm optimizaiton (NPSO) is presented. NPSO combines traditonal particle swarm optimization with the fitness-sharing principle. By applying four test functions and parameter inversion based on the Muskingum model and comparing these with traditional particle swarm optimization, the efficiency in convergent speed and precision of this method are verified. However, inversion results are not good because of stochastic noise. To improve the antinoise capability of NPSO, a multiscale NPSO is constructed by combining the multiscale strategy with NPSO and by applying the parameter inversion of the Muskingum model with 5% stochastic noise. Inversion resultsverify the effectiveness of the improved algorithm; the antinoise performance of the NPSO has been increased, and the precision of the parameter inversion result is significantly improved.

Keywords：Muskingum model; parameter inversion; niche particle swarm optimization; multiscale; stochastic noise

中圖分類號：X522

文獻標志碼：A

文章編號：1006-7043(2016)02-0271-06

doi:10.11990/jheu.201407078

作者簡介：張新明(1979-)，男，副教授，博士.通信作者：張新明，E-mail:xinmingxueshu@gmail.com.

基金項目：國家自然科學基金資助項目(41004052).

收稿日期：2014-07-31.網絡出版日期：2015-12-29.

網絡出版地址：http://www.cnki.net/kcms/detail/23.1390.u.20151229.1711.010.html