马斯京根模型参数反演的改进粒子群算法

张新明，马艳

(哈尔滨工业大学 深圳研究生院，广东 深圳 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