环境模型参数优化方法的比较

摘要：模型参数优化是通过极小化目标函数使得模型输出和实际观测数据之间达到最佳的拟合程度. 由于环境模型本身的复杂性,常规优化算法难以达到参数空间上的全局最优. 近年来,随着计算机运算效率的快速提高,直接优化方法得到了进一步开发与广泛应用. 本文比较了CRS、SCE UA、SA 和Annealing2Simplex 等4 种算法应用于环境模型参数优化的结果和计算效率。 关键词：参数优化 环境模型 CRS 算法 SCE UA 算法 Simulated2Annealing 算法 Annealing2Simplex 算法 人们对环境系统的深入研究是建立在环境模型的广泛应用基础上的. 为了更加精确地刻画环境系统的行为,环境模型在近10 年里表现出了强烈的复杂化趋势;不同空间尺度、不同时间过程模型的耦合,进一步加剧了这一过程. 环境模型的复杂性导致了模型结构和参数可识别性问题的提出,并成为当今环境建模理论研究的热点[1 ,2 ] . 其中在不确定性的框架下,模型参数的优化是研究的一个重要方面. 解决优化问题的难度主要取决于模型参数的空间维数和模型本身的非线性特征. 一般来说,参数越多、非线性越强,优化时间和精度就越差,同时也越不能够保证优化算法是否收敛到整体最优. 传统经验表明,求解优化问题的困难主要体现为[4 ,5 ] : ①全局搜索可能收敛到多个不同的吸引域; ②每一个吸引域可能包含一个或多个局部最小值; ③目标函数在n 维参数空间上不连续; ④参数及相互间存在高度灵敏性和显著非线性干扰; ⑤在最优解的附近,目标函数往往不具有凸性。 优化算法可以分为直接算法和间接算法2大类. 间接算法(如牛顿法以及各种以牛顿法为基础的改进算法) 的局限性主要在于要求目标函数在相关值域上必须是可微的;而直接算法仅涉及目标函数值的计算,不需要计算目标函数的导数. 因此尽管后者的计算效率相对较低,但在环境问题的实际应用中,它可以有效和简洁地解决由于模型复杂性所衍生的不可微函数的优化问题. 同时,与求解问题所耗费的时间相比,通常更为关注解的结果能够在多大程度上描述系统的行为[1 , 2 ] . 所有这些使得直接算法在近年来得到了迅速发展和广泛的应用. 本文的目的在于分析和比较几种近年来渐为接受的参数直接优化算法的计算效率、计算精度及其算法稳定性. 由于直接算法本质上的随机性,以及对于不同优化问题所表现出的算法特性上的差异,因此本文仅以经典测试函数和环境水文模型为例,对几种算法进行详细的比较研究. 1参数优化算法 模型参数的优化即是寻求一组参数,使模型的输出与实际观测数据之间按给定目标函数的度量方式达到最佳拟合,即: Min f ( ys , yobv ,θ) = f ( ys , yobv ,θ3 ) θ ∈ S (1) 式中, f ( y , θ) 为目标函数; ys 表示模型输出变量, yobv为系统观测值;θ3表示参数可行域S 上的最优参数向量. 最早提出的直接算法均是简单随机方法,这些算法除了计算效率低外,主要缺陷在于要求目标函数在邻域上必须是不相关的. 控制随机搜索算法(CRS) [10 , 11 ] 有效地克服了简单随机算法的主要缺点,在计算过程中保存指定数目的参数样本及其对应的目标函数值,引入几何学中“重心”的概念,即考虑了新点产生的随机性,又在一定程度上保证了搜索的整体性. 复合形混合演化算法( SCE UA) [4 , 5 ]所采用的竞争演化和复合形混合的概念继承了CRS 算法中全局搜索和复合形演化的思想,该方法将生物自然演化过程引入到数值计算中,模拟了生物进化的过程,提高了计算效率和全局搜索整体最优的能力. 模拟退火算法(SA) [7 , 12 ]则假设优化问题的解及其目标函数分别与固体物质的微观状态及其能量所对应,采用蒙特卡洛(Monte Carlo) 随机方法模拟固体稳定“退火”的过程,并假设优化过程中递减目标函数值的控制参数t 与“退火”过程中的温度T 所对应,对于控制参数t 的每一个值,算法持 续进行“产生新解2判断2接受/ 舍弃”的迭代过程,应用该算法的关键在于确定合理的冷却进度表. 退火单纯形算法(AS) [8 , 9 ]综合了下山单纯形方法和模拟退火法2 种优化算法,充分利用单纯形的形变信息,以一定温度T 所对应的概率接受准则来指导单纯形的映射、压缩和扩展等过程,提高了计算效率和算法稳定性. 运用直接算法进行参数优化首先要正确选取算法的内部控制参数. 一般根据具体问题的特点,采取数值实验的方法不断测试参数“性能”,从而确定既保证算法收敛到满意的极值,又不至于使计算过于耗费时间而导致算法失效. 选择控制参数本身实际就是一个“优化”过程. 2测试函数的优化 首先采用一个广泛使用的、具有多个局部极小值的经典函数,对上述几种优化算法的全局搜索能力进行较为全面的测试,其函数形式为: 尽管测试函数仅有x 、y 2 个变量,但其等值线却具有复杂的空间结构,即多个局部极小、值域空间上的不连续和全局最小值周围的香蕉状狭长低谷(图1) . 现已知参数全局最小值为(5 ,5) ,所对应的目标函数值为10.0。 如前所述,为了保证直接算法的优化精度和避免算法失效,在算法的运行中需要在如下３个优化终止准则中进行选择和组合: ①前m组最佳估计参数收敛到可接受的参数空间,即最优参数的不确定性小于ε1 ; ②最佳目标函数值的优化改进速率小于给定值ε2 ; ③搜索过程中参数采样总数不超过给定的计算次数n. 表1 给出了4 种直接算法分别基于终止准则1 和2 ,以及不同初值假设条件下的测试函数参数最优值和计算效率. 为了尽可能消除算法随机性对于算法比较的影响,表格中数据除了计算次数的偏差外,均为相同条件下5 次运算结果的算术平均值. 由于AS 算法中的下山单纯形方法具有很强的方向性,因而在两个终止准则约束下,算法均可迅速达到全局最优化,且结果具有较高精度. 尽管SCE UA 算法引入竞争演化思想后提高了样本空间的搜索效率,但是不同种群从不同方向上向全局最优点逼近,以及同一种群内部仍然采用的单纯形控制随机搜索方法,使得该算法在计算效率上远不及AS 算法. 然而SCEUA 算法在各种初值条件下均可得到十分精确的全局优化结果,可见该算法的稳定性和可靠性比AS 算法具有更大的优势. 与前2 种算法相比,CRS 算法和SA 算法在计算效率和精度方面均表现不佳. 特别是当参数初值发生变化时,必须相应地调整控制参数才能得到较好的计算结果,这说明2 种算法的适应能力较差. 表1 4 种直接算法的参数优化结果1) Table 1 Parameter optimization results of four direct optimal algorithms 优化算法 编号 假设初值 参数优化结果 目标函数值 计算次数 准则1 准则2 x y x y x y 准则1 准则2 准则1 准则2 偏差/ % CRS 1 0.40 0.60 0.5000 0.5001 0.5000 0.5001 10.0001 10.0001 565 584 -3.25 2 0.10 0.20 0.5000 0.4999 0.5000 0.5000 10.0001 10.0000 806 841 -4.16 3 0.20 0.80 0.5001 0.5001 0.5000 0.5000 10.0002 10.0000 585 618 -5.34 4 0.70 0.70 0.4998 0.5001 0.5000 0.4999 10.0002 10.0000 541 571 -5.25 5 0.80 0.85 0.5000 0.5000 0.5000 0.5000 10.0000 10.0000 595 583 2.06 SCE UA 1 0.40 0.60 0.5000 0.5000 0.5000 0.5000 10.0000 10.0000 479 356 34.55 2 0.10 0.20 0.5000 0.5000 0.5000 0.5000 10.0000 10.0000 494 350 41.14 3 0.20 0.80 0.5000 0.5000 0.5000 0.5000 10.0000 10.0000 501 359 39.55 4 0.70 0.70 0.5000 0.5000 0.5000 0.5000 10.0000 10.0000 431 305 41.31 5 0.80 0.85 0.5000 0.5000 0.5000 0.5000 10.0000 10.0000 480 335 43.28 SA 1 0.40 0.60 0.5001 0.5002 0.5000 0.5000 10.0013 10.0000 230 349 -34.10 2 0.10 0.20 0.5041 0.4973 0.5000 0.5000 10.1300 10.0000 408 504 -19.05 3 0.20 0.80 0.4999 0.5002 0.5000 0.4999 10.0002 10.0000 468 579 -19.17 4 0.70 0.70 0.5002 0.4999 0.5000 0.4999 10.0005 10.0001 468 582 -19.59 5 0.80 0.85 0.5001 0.5000 0.5001 0.5000 10.0000 10.0000 490 544 -9.93 AS 1 0.40 0.60 0.4999 0.5001 0.4999 0.5001 10.0000 10.0000 73 63 15.87 2 0.10 0.20 0.5000 0.5001 0.5000 0.5001 10.0001 10.0001 150 242 -38.02 3 0.20 0.80 0.5000 0.5000 0.5000 0.5000 10.0000 10.0000 108 135 -20.00 4 0.70 0.70 0.5001 0.4999 0.5001 0.4999 10.0000 10.0000 99 191 -48.17 5 0.80 0.85 0.4999 0.5001 0.5001 0.4999 10.0000 10.0001 157 157 0.00 1) 计算次数偏差的计算公式为( n准则1 – n准则2) / n准则2 表2 目标函数计算次数的统计特性 Table 2 Statistical characteristics of the calculatingtimes of objective function 计算次数的统计量 CRS SCE- UA SA AS 准则1的标准方差 106.89 27.36 106.65 35.46 的标准方差 114.05 22. 15 96.24 66.52 准则1和准则2偏差均值/% -3.19 39.97 –20.37 –18.06 可以采用各种初始条件下对应同一终止准则计算次数的标准方差来表征算法的随机性.尽管基于控制性随机搜索思想构造的直接优化算法不可能完全消除计算结果的随机性,但是算法结构的精心设计将在一定程度上减小随机性影响,同时也就意味着算法稳定性的增大. 表2 中结果显示SCE UA 算法的计算次数在各种初始条件下的标准方差远小于其他3 种直接算法,再次表明了算法所具有的良好稳定性. 这一结果的产生是由于SCE UA 算法结构中引入的复合形概念及其竞争混合算法大大削弱了计算结果的随机性. 对应准则1 和2 的计算次数偏差,表征了算法对于选用不同终止准则的敏感性. 根据表2 结果可以得到4 种算法的敏感性排序为(从大到小) : SCE UA > SA > AS > CRS ,特别是CRS 算法的计算次数偏差的最大值和统计平均值分别为- 5.34 %和- 3.19 %,远小于其他3种算法. 更进一步发现,算法结构复杂性与算法对于终止准则选用的敏感性之间高度相关,表现为本例中算法敏感性排序与复杂性顺序完全一致. 其原因在于,直接算法在结构复杂性所形成的内部约束与目标函数的外部约束的共同作用下趋向全局最优,其结果必然是一个不断协调和折衷的过程. 算法越复杂,内部约束就越强,其计算结果就越发表现出与终止准则之间的高度相关性. 因此简单算法对于外部终止准则控制性约束的适应能力较之复杂算法更为灵活.