基于Kriging模型的自适应多阶段并行代理优化算法

乐春宇，马义中

(南京理工大学 经济管理学院，江苏 南京 210094)

0 引言

现代工程优化设计中，常采用高精度仿真模型获取数据，如有限元分析和流体动力学等[1]，如何在优化过程中尽可能少地调用高精度仿真模型，以提高优化效率，显得尤为重要。由于优化过程的黑箱特性无法获取过程的梯度信息，限制了基于梯度信息的高效优化算法的使用，使基于非梯度的启发式优化算法得到快速发展[2-3]。在已有启发式算法中，遗传算法(Genetic Algorithm, GA)和粒子群优化(Particle Swarm Optimization, PSO)算法是解决黑箱问题最常用的算法，然而启发式算法的最优解需要通过迭代在大量个体中竞争产生，不适用于需要使用昂贵的仿真模型来获取数据的工程优化设计。

代理优化(Surrogate-Based Optimization, SBO)可以用少量样本点构建代理模型，并对真实过程进行探索，逐步发现过程的最优解，是解决高耗时黑箱工程优化的有效方法[4-5]。相对于多项式响应曲面、径向基函数、支持向量回归和神经网络模型等，Kriging模型在给出预测值的同时还能给出预测偏差[6]，结合预测值和预测偏差能够生成多种加点准则(Infill Sampling Criteria, ISC)，通过优化ISC选择新的样本点，能够自适应地更新模型并选择下一个样本点，因此Kriging模型成为最具代表性的代理模型。在常用的ISC中，JONES等[7]提出高效全局优化(Efficient Global Optimization, EGO)算法，其中使用的EI(expected improvement)准则兼具全局开发(exploration)和局部探索(exploitation)的特性，是最经典的ISC；WATSON等[8]提出的WB2准则在EI准则的基础上增加了预测均值的负数项，更加注重局部开发，其寻优效率和精度更好[9-10]；SCHONLAU等[11]提出的可行性概率(Probability of Feasibility, PoF)能够计算未知点落在约束内的概率，是代理优化中最常用的约束处理准则。在每次迭代过程中，EGO算法只使用EI值最大的一个点更新模型，在计算资源充足和并行仿真技术越来越完善的情况下，同时仿真多个点的数据和仿真一个点的数据所花费的时间资源相同。因此，如何同时选取多个点作为新的样本点更新代理模型，以减少迭代次数，提高优化效率，受到广泛关注。

代理优化中，在迭代中加入多个样本点的方法称为多点填充或并行优化，其主要目的是使用并行计算资源，同时获取多个新样本点的真实响应值，从而压缩优化周期。常用的多点填充策略有3类：

(1)同时使用多个方法产生多个样本点。HAMZA等[12]为了避免部分优化算法求解EI值时陷入局部最优解，采用多种优化算法优化EI值来产生多个样本点；HUTTER等[13]通过一个代理模型设置不同参数来产生多个样本点；VIANA等[14]、YE等[15]用多个代理模型产生多个点。这类策略操作简单，但加点数量受限，不够灵活。

(2)多次更新同一个准则产生多个点。GINSBOURGER等[16]提出的KB(Kriging believer)和CL(constant liar)算法采用“假”响应值替代真实响应值，并不断更新Kriging模型，能够一次产生多个样本点；ZHAN等[17]提出伪装期望改进算法，采用新的样本点构造影响函数，再用影响函数更新EI并产生下一个样本点。这类策略可以加入任意数量的样本点，然而如果同时加入的样本点较多，则需多次更新ISC，会降低其精确度。

(3)结合多目标优化算法，同时优化多个准则，产生Pareto解集，再从解集中选取特定数量的点。对于无约束优化，FENG等[18]将EI准则的两部分同时作为优化目标，产生Pareto解集；LI[19]采用Kriging模型的预测均值和预测偏差作为优化的两个目标；WANG等[20]综合考虑新加入点的失败率和提升率，将改进概率(Probability of Improvement，PI)和EI两个准则作为多目标优化的两个目标；对于约束优化，PARR等[21]直接将PoF准则和EI准则作为两个目标同时优化；HABIB等[22]在约束范围内进行局部开发和全局探索，将EI准则的两部分同时乘以PoF准则作为多目标优化的两个目标；HABIB等[23]考虑加入样本点的空间填充性，将距离以及EI准则和PoF准则的乘积值作为两个优化目标。这类多点填充策略研究最为广泛，其做法是选取多个合理的ISC作为多目标优化的目标。

已有的多点填充算法重复使用某种策略选取样本点，直至达到终止条件，其模式比较单一。每种准则都有各自优势和不同适用场景，注重全局探索的准则更适合探索整个区域，注重局部开发的准则更适合在已知最优解附近加入样本点。为了充分利用每种ISC在不同阶段的优势，进一步提高代理优化的效率，本文提出一种基于Kriging模型的自适应多阶段并行代理优化算法(Parallel Surrogate-based optimization algorithm based on Kriging Model using an adaptive multi-phases strategy, PSKM)，该算法是一种针对约束问题的多点填充代理优化算法，在全局探索阶段采用EI和PoF准则探索可行域和约束边界，在局部开发阶段采用WB2和PoF准则开发约束问题的最优解，同时制定了自适应的两阶段切换策略，使算法能够根据加点的情况选择更好的加点策略。两个阶段均采用文献[21,24-25]的多目标优化框架应对约束并产生备选点集，再选取多个样本点完成多点填充。最后通过3个数值案例和一个工程实例，将所提算法与已有的多点填充算法进行对比，验证了所提算法的性能。

1 有关概念界定

1.1 单目标约束优化

本文需要解决的单目标约束优化问题可以用以下数学式表达：

s.t.

gi(x)≤0,i=1,…,m;

x∈∈d。

(1)

式中：x为d维自变量向量，y(x)为其响应值；gi(x)为约束函数，m为约束的个数。文中假设y(x)和gi(x)均为复杂的黑箱过程，需要调用昂贵的仿真模型来获取数据。

1.2 Kriging模型

假设X=[x1,x2,…,xn]为n个已知样本点，其响应值Y=[y1,y2,…,yn]。Kriging代理模型需要通过已知的n个点预测样本空间中任意位置点的响应值，其预测值为

(2)

式中C=[c1,c2,…,cn]为已知样本点在预测点中所占的权重。Kriging模型将响应函数看作某个正态随机过程的实现，该随机过程定义为

Y(x)=u+Z(x)。

(3)

(4)

式中：r为已知点和预测点之间的相关系数向量；R为已知点之间的相关系数矩阵；F是长度为n的列向量，F=[1,1,…,1]T。模型具体的推导以及参数说明，可参考文献[6,26]。

1.3 加点准则

ISC是由模型产生、用于指导寻找新样本点的一种标准，本节简单介绍文中需要使用的3个ISC。

(1)EI准则

该准则计算未知点给最优解带来改进的值为

I(x)=max(ymin-y(x),0)。

(5)

(6)

式中Φ和φ分别为标准正态分布的累积概率函数和概率密度函数。EI准则相加的两部分分别体现局部开发和全局探索的能力。在已知样本点处的s2(x)=0会导致该点的EI准则值为0，因此用EI准则选点能够保证已知点不会被重复选中。另外，EI准则能够更好地寻找可能存在的最优解，具有更好的探索性。

(2)WB2准则

WB2准则更加注重局部开发，是在EI准则基础上改良的一种ISC，其表达式为

(7)

由式(7)可见，WB2准则是在EI准则的基础上加上预测均值的负数项，其兼具局部开发和全局探索的特性。因为预测均值负数项的加入增加了局部开发的权重，所以WB2准则能够集中在已知最优解处加点，从而提高寻优精度。

(3)PoF准则

PoF准则能计算样本点落在可行区域的概率值，在假设约束相互独立的条件下，其表达式如下：

(8)

1.4 自适应代理优化框架

EGO算法是最经典的自适应代理优化算法，现有的代理优化和EGO类似，均具有相同的优化框架，其主要分为两个阶段：①采用适当的抽样方法，抽取少量初始样本点，建立Kriging模型；②采用Kriging模型生成ISC，优化ISC并找到使ISC达到最大的点，获取该点响应值将其加入样本点中并更新代理模型，重复执行直到达到终止条件。ISC通过代理模型产生，优化ISC不需要调用昂贵的仿真模型，可减少仿真次数，提高优化效率。

2 PSKM分析及实现

PSKM内容主要包括3方面：①在约束下的全局探索和局部开发两个阶段，每个阶段如何选择样本点并完成多点填充；②如何确定两个阶段的切换条件，以保证两个阶段自适应切换；③根据前面两点的内容，将所有子步骤整合成一个完整的并行代理优化算法。下面从这3方面对PSKM进行分析和说明。

2.1 约束下的全局探索和局部开发

EI准则表达式中，标准正态分布的累积概率函数和概率密度函数可以看作为局部开发和全局探索的权重。WB2准则在EI准则的基础上增加了预测值的负数项，相当于增加了局部开发的权重。EI和WB2均兼具全局探索和局部开发的特性，其中EI准则的全局探索性能更好，能更快地在最优解附近加点，但是最优解的精度不高；WB2的局部开发性能更好，寻找最优解的速度慢，然而精度较高。因此，PSKM在全局探索和局部开发两个阶段，分别采用EI和WB2准则。

关于约束处理问题，两个阶段都采用PoF准则。EI准则的数值会慢慢降低至0，WB2准则的数值量级由真实问题决定，而PoF取值在0～1之间。为了避免某一个准则占据主导地位，导致探索约束边界能力不足，两个阶段均使用多目标优化框架。约束下的全局探索阶段将EI和PoF作为两个目标同时优化，局部开发将WB2和PoF作为两个目标同时优化，两种方法的数学表达式如下：

x∈∈d；

(9)

x∈∈d。

(10)

采用NSGA-Ⅱ可以求解出相应准则的折中解，对样本点进行第一层筛选，将备选点从整个取值范围缩小到Pareto解集上，然后在Pareto上进行第二层筛选，挑选出指定数量的样本点。假设每次迭代需要加入q个新的样本点，可以对Pareto解集进行聚类，将其分割成互不相交的q个子集，再用特殊的选点准则从每个子集中挑选出一个样本点即可，本文采用K均值聚类算法对解集进行分割。根据指定的多点填充数量q，将生成的Pareto解集P分割成q个子集P1,P2,…,Pq，q个子集互不相交且和为P，即:

Pi∩Pj=∅,i≠j；

P1∪P2∪…∪Pq=P。

(11)

然后用式(12)的选点标准进行第二层筛选，在每一个子集中选出一个点组成xnew，同时用仿真模型获取q个点的响应值，将其作为新的样本点加入已有样本点中更新模型，从而实现多点填充。

xnew=[xnew1,xnew2,…,xnewq]。

(12)

选点标准使用的是EI准则值和PoF准则值的乘积。可以从每个子集中选出尽可能满足约束且EI准则值尽可能大的点。因为EI准则不会重复选择已知样本点，所以能够保证算法的收敛性。最后，将以上步骤整合为约束下两阶段并行优化的两个子算法，分别命名为Cglobal和Clocal，算法步骤如表1所示。

表1 Cglobal和Clocal算法步骤表

2.2 两阶段自适应切换算法

表2 G2L和L2G算法步骤表

2.3 PSKM流程

3 算例分析

本章用3种不同类型的数值算例和1个实际工程实例验证PSKM的性能，并与已有的几种经典代理优化算法进行对比。单点填充算法采用文献[11]中经典的CEI(constrained expected improvement)算法。对于多点填充算法，第二类多点填充算法选取文献[16]的KB算法，为使该算法能够处理约束，本文将其中的EI准则换成EI准则值乘以PoF准则值，将这种算法简称为KBMP(Kriging believer multi-points)；第三类多点填充的算法选取文献[21]中的双目标多点填充算法，简称为DOMP(double-objective multi-points)。

3.1 算例介绍

3个经典的数值算例是来自文献[25]中的CAM(camel)，SAS(sasena)，GOM(gomez)。CAM是一个单约束多极值优化问题，该算例在约束边界的拐角处有多个局部最优解和一个全局最优解；SAS是一个多约束多极值问题，其有3个约束条件，2个有效约束和1个无效约束；GOM是单约束且可行域不连续的约束优化问题，其可行域的范围很小，会出现初始样本点不在可行域中的情况。工程实例来自文献[27]的拉伸弹簧设计问题TSD(tension spring design)，设计变量为线圈的平均直径、弹簧的线径和弹簧的有效线圈数，设计目标为：在线圈偏转度和线圈强度等约束下，使弹簧的重量最低。

3.2 试验设置

(1)设置初始样本点 所有试验的初始样本点均采用最优拉丁方设计抽取，数量为8d，其中d为自变量的维度。采用最优拉丁方设计时，设置不同的随机种子会得到不同的样本点。如果初始样本点中存在离最优点比较近的点，则将导致算法更快地找到全局最优解。为了屏蔽初始样本点的影响，采用20个不同的随机种子产生初始样本点，每个算法对每个算例均进行20次试验。

(2)设置NSGA-Ⅱ 种群数量设置为80q，q为多点填充的数量，最大迭代次数为100，交叉分布指数为5，交叉概率为0.7，变异分布指数为10，变异概率为0.2。

(3)设置终止条件 对于3种多点填充算法，设置每次迭代加入3个新的样本点。所有试验都以加入点的总数T作为终止条件，设CAM，SAS，GOM的终止条件T=30，CEI算法迭代30次，多点填充算法迭代10次；设TSD的T=45，CEI算法迭代45次，多点填充算法迭代15次。

3.3 比较指标

为了更好地对比说明算法性能，本文从以下3个指标对试验结果进行分析：

(1)算法的绝对误差 对于每个算例，计算每个算法20次试验中每次迭代产生的最优解的均值f(xbest)，并计算其与真实最优解f(x*)之间的绝对误差，以此衡量算法收敛到全局最优解的效率。算法的绝对误差

(13)

(2)算法的最优解性能 对于每个算例，统计每种算法20次试验最优解的性能，包括最值、均值和标准差，并通过对20次试验结果做t检验来衡量每种算法寻找最优解的能力，即算法寻优的有效性。

(3)算法的稳健性 统计20次试验结果的数据，采用绝对误差以10为底对数形式的负值制箱线图，通过观察每种方法的波动情况来衡量算法寻优的稳健性。

3.4 结果分析

根据4.3节提到的比较标准对试验数据进行处理，结果如下：

(1)算法的收敛图

由于每个算例结果的跨度不同，为了更好地观察结果，对AE值取以10为底的对数，并取负数形式，绘制4个算例的算法收敛图，如图2所示。

图中横坐标为新加入样本点的个数，纵坐标为AE以10为底对数形式的负值，其值越大表示算法求解的最优解和全局最优解的差值越小。CEI算法每次迭代只加入一个新的样本点，因此其收敛的曲线更加密集。由图可见，PSKM每次迭代可以加入多个样本点，相比单点填充CEI，相同精度下PSKM的迭代次数更少，说明了PSKM的多点填充策略在加快代理优化算法收敛速度方面的有效性。考虑总共消耗的计算资源，即加入相同样本点数观察算法的收敛情况，对于算例CAM和SAS，PSKM的收敛曲线均在其他算法的左上方，说明PSKM的收敛速度更快，最终求解的最优解更好；对于算例GOM和TSD，PSKM虽然在前期寻优速度较慢，但是因为在局部开发阶段采用WB2准则能找到更精确的最优解，所以随着加点的进行，其收敛曲线会超过其他算法。

(2)算法最优解性能表

对于每种算例中每种算法的20次试验，计算20次试验结果的最小值、最大值、均值和方差，并对每个算例下算法两两之间的试验结果进行t检验，检验置信度为95%时其均值是否较小，结果如表3所示。

表3 算法最优解性能表

4个算例中，PSKM的最小值、均值和最大值均优于其他3个算法。因为PSKM不但能在局部开发阶段用Clocal充分开发最优解区域，寻找精度更好的结果，而且能自动切换状态，在充分开发某一区域后，用Cglobal算法继续探索其他区域，所以对于拥有局部最优解的CAM和SAS，PSKM也能及时跳出局部最优，寻找更好的最优解。从结果的标准差来看，对于GOM算例，DOMP算法的所有结果较大，导致标准差偏小，因此其总体结果不如PSKM；对于其他3个算例，PSKM算法的标准差更小。综上所述，从最值、均值、标准差和t检验的结果来看，PSKM的结果都不比其余3种算法差，说明PSKM求解问题更有效，所求解的精确度更高。

(3)算法最优解箱线图

绘制每个算例中每种算法20次结果的箱线图，如图3所示，图中纵坐标为AE以10为底对数形式的负值。可见4种算例的20次试验结果，PSKM的箱线图更靠近上方，说明PSKM算法屏蔽初始样本点影响的能力更好，算法的稳健性更好。

4 结束语

在计算资源充足且并行计算技术越来越成熟的情况下，在更新代理模型过程中一次加入多个样本点，可以将充足的计算资源转换成横向的时间资源，从而缩短优化周期，提高优化效率。现有的多点填充代理优化算法都是重复使用某种策略加点，模式较单一，为了充分利用每种ISC在不同阶段的优势，进一步提高代理优化的效率，本文提出PSKM，算法内部可以自适应地选择全局探索或局部开发，并在不同阶段采用不同子算法寻找新的样本点。算例的结果表明，PSKM不但收敛速度快，而且求解的精确性和稳健性更好。