周珍胜 王 林 冯夫健 谭 棉 何 兴 张再军
进化算法(Evolutionary Algorithm, EA)是一类模拟生物进化过程的自适应搜索算法,主要包括遗传算法(Genetic Algorithm, GA)、进化策略(Evolu-tion Strategy, ES)、遗传编程(Genetic Programming, GP)和进化规划(Evolutionary Programming, EP)等[1].EA为解决实际的复杂问题提供一种高效可行的思路,能在大规模、高维度的问题空间中搜索全局最优解,可广泛应用于组合优化、机器学习、数据挖掘等领域.EA计算时间分析可以通过首达时间刻画.首达时间是指算法首次找到最优解或满意解时的迭代(运行)次数[2].由于算法随机性的本质,算法运行一次并不能较好地刻画算法性能,因此,可借助平均首达时间评估和分析算法性能[3].
平均首达时间是指算法首次找到最优解或满意解时,对适应度函数的平均评估次数.平均首达时间的大小对应算法的全局搜索能力和收敛速度,值越小,收敛速度越快,越容易找到全局最优解.目前,进化算法在离散搜索空间和连续搜索空间优化问题上取得显著性成果,尤其连续型进化算法的计算时间分析引起学者的广泛关注[4].
在针对离散型进化算法的计算时间分析中,He等[5]将进化算法运行过程中种群状态建模成马尔可夫链,并引入漂移分析和距离函数,分析父代个体与子代个体之间的期望变化,得出平均首达时间上界.在此基础上,He等[6]提出进化算法时间分析的一般框架.Yu等[7]基于非齐次马尔可夫链,建立首达时间与收敛速度之间的桥梁,提出进化算法平均首达时间转换分析方法.Lengler等[8]结合漂移分析,假设父代个体优于子代个体,分析进化算法在线性函数上有更紧的平均首达时间上界.Yu等[9]提出转换分析方法,构造一条新的马尔可夫参考链,对比两条参考链之间每步计算时间的差异性,从而得出整体算法的计算时间复杂性.Wegener[10]通过种群个体适应值,将种群个体划分为从低级到高级的集合,分析种群从低级到高级的迁移概率界,得到计算时间复杂性适应度分层方法.Huang等[11]结合统计方法和平均增益模型,采集种群个体增益和曲面拟合技术,得到增益的概率分布函数,进一步推出估计进化算法的平均首达时间数学闭合解析式,并以此分析进化算法在球函数上的计算时间,有助于估计进化算法首达时间.冯夫健等[12]提出基于等同关系的演化算法的时间对比分析方法,不需要定义距离函数和满足漂移分析的特定条件,为算法改进提供理论上的依据.现有的乘法漂移[13-15]和可变漂移[16-18]证明平均首达时间下界的证明条件都强于时间上界的证明条件.结合马尔可夫链与漂移分析的模型是一种分析进化算法平均首达时间的有力工具.
综上所述,现有进化算法首达时间分析主要集中在对离散进化算法的分析,而连续型进化算法首达时间的分析相对较少,并且涉及复杂的矩阵计算.在连续型进化算法首达时间分析中,Jiang等[19]结合增益,分析(1+1)进化策略在超球体内使用均匀变异算子时,当种群个体到最优解距离减少到一半时,运行时间的下界为指数级.张宇山等[20]将首达时间视为一个停时,结合时齐马氏过程性质,提出(1+λ)进化策略平均首达时间停时理论模型.Agapie等[21]结合增长率,分析带高斯变异和均匀变异进化策略在球函数的首达时间是Ο(n),理论分析影响(1+1)进化策略的首达时间并非是变异算子,而是1/5成功选择的影响.Jägersküpper[22]研究进化策略在球函数上运行时间是如何取决于搜索空间维数,并进一步研究单模函数优化的计算时间,计算相对冗长.Agapie等[23]把父代个体与子代个体之间每步增长率建模成一个更新过程,分析带均匀变异算子的进化策略在倾斜平面上的计算时间复杂性,以此证明计算时间复杂性只与倾斜平面区域有关,而未关注计算时间复杂性与种群规模相关,导致得到的平均首达时间下界不紧致.张宇山等[24]基于平均增益,结合鞅论与勒贝格控制收敛定理,分析进化策略在球函数和二维倾斜平面的平均首达时间更紧的上界,而未给出相关理论下界分析,同时在上界理论分析中勒贝格控制收敛定理需要函数列(T0∧s)存在控制函数,这需要较强的假设.
总之,虽然连续型进化算法的首达时间分析在首达时间上界研究中已有一定成果,但是学者们较少关注其首达时间下界,以及相关的理论结果是针对具体案例展开,导致进化算法首达时间分析没有抽象为一般模型.
本文将进化算法首达时间视为更新过程的停时,引入鞅论和更新过程,结合瓦尔德不等式,推出文献[24]的下界闭合表达式,提出估计进化算法的平均首达时间分析的一般模型.本文的更新理论模型是建立在非负增长率随机过程上,不依赖算法的具体实现,有利于分离具体案例和算法,使模型更一般化.通过更新理论模型上界和下界,构建一个闭合区间,分析进化策略种群规模与平均首达时间的条件关系,得到平均首达时间与种群规模之间并非负相关.
本文重点讨论的是进化算法的随机过程模型,不失一般性,研究连续搜索空间的最优化问题,数学描述如下:
minf(x),
s.t.x∈S⊆Rn,
其中,f(x)表示目标函数,S表示搜索空间,n表示问题的规模.S中每个元素表示问题的一个解,对应算法中的每个个体,这些个体组成种群.记进化策略第t代种群
其中λ表示种群规模.一个典型ES求解最优化问题算法流程如算法1[25]所示.
算法1ES流程[25]
输入初始化种群
输出包含目前为止找到满意解
令t=0;
计算每个个体适应值,评估当代种群中个体优劣;while not stop do
fori=1∶λdo
通过交叉和变异产生新的种群;
end for
计算每个个体适应值,评估当代种群中个体优劣;
从λ个个体中选择ξ′t个体作为下一代的种群;
t=t+1;
end while
输出到目前为止找到的满意解.
算法1包括使用交叉、变异和选择,未对交叉和变异做任何限制.选择是采取非精英保存,即每代维持一个父代个体,通过变异产生λ个子代个体.再从λ个子代个体中选择个体作为下一代的种群.然而,本文给出的主要结果(即定理3)与这样的具体实现无关,几乎可以用于任何随机算法.
为了建立非负随机过程,本文假设在进化策略运行过程中,第t代种群中至少存在一个个体的适应值优于第t-1代种群个体的适应值.
由于算法很难找到最优解,当算法找到满意解
时,即
时,算法终止.
根据随机过程理论,设(Ω,H,P)为概率空间,
Ht=σ(d(ξ0),Y1,Y2,…,Yt)⊂H
为σ-代数流,H包含由信息d(ξ0),Y1,Y2,…,Yt生成的所有事件.鞅这个概念最初起源于公平赌博数学模型,后引入概率论中.鞅就是在掌握从初始到现在累积信息d(ξ0),Y1,Y2,…,Yt条件下预测未来与现在相同.
E(Xt+1|Ht)≥Xt,
E(Xt+1|Ht)≤Xt,
定义1适应值差函数 设函数d在S→R满足
定义2增长率[26]在ES运行过程中,增长率表示t-1代到t代种群之间到最优距离的减少量,即
质量增益表示父代与子代之间的适应值减少量[27].
则称Yt为进化策略在第t代增长率更新过程.
在随机算法中,研究者们通常感兴趣的是算法首次找到满意解时,算法对适应度函数的迭代次数[4], 即父代与子代之间到最优距离的减少量之和
Yt=d(ξ0)-ε.
首达时间可使用停时描述,下面给出停时的定义.
则称τ为Xt的一个更新过程停时.
与目前进化算法[28-29]研究不同,本文将ES从初始位置到满意解的渐进过程视为随机状态建立更新理论模型.在建立模型之前,需要给出定理1.
在定理1中,得到增长率更新过程期望与平均首达时间之间的等式.增长率更新过程的期望是平均首达时间与期望增长率的乘积.然而,期望增长率的概率密度函数与初始种群有关,推论1给出增长率更新过程的期望上界和下界.
E(Tτ|Ht)<∞,
μ1 则 根据推论1,可得到平均首达时间的更新定理,如定理2所示. μ1 若当dTτ→∞时,则对E(Tτ|Ht), 成立. 证明因为 dTτ+1 于是可得 而Tτ时刻时到达距离 在E(Tτ|Ht)之前到达距离dTτ的均值 由强大数定理可得,当E(Tτ|Ht)→∞时, 但是由于dTτ→∞时, E(Tτ|Ht)→∞, 所以dTτ→∞时, 又 类似可推得,当E(Tτ|Ht)→∞时, 由极限夹逼定理可得,当E(Tτ|Ht)→∞时, 证毕. μ1 证明对于下界证明,定义4中停时 等价于 根据瓦尔德不等式,有 则 而对于上界证明,记另一更新过程为 对于常数a>0, 从而 因此,由推论2可推出 即 故 证毕. 在进化策略运行中,当父代个体到最优解的适应值差函数值等于ε时,算法终止,即到达某种特定的停时距离: dTτ=Yt=d(ξ0)-ε. 则平均首达时间定义如定义5. 由于进化算法首达时间与初始种群的状态有关,为了刻画初始种群转移状态与首达时间的关系,E(Tτ|Ht)可简化为E(Tτ|d(ξ0)).推论3给出E(Tτ|d(ξ0))的估计式. 若存在2个关于初始种群的函数ψ1(ξ0)、ψ2(ξ0),使得 ψ1(ξ0) 则对E(Tτ|d(ξ0)), 成立. ∀t=1,2,…, Yt=d(ξ0)-ε=d(ξ0), 即ε=0,则在连续情形(对于所有t)下,μ1>0不存在.并且在连续搜索空间下,算法很难搜索到局部最优解.因此本文在收敛时间中使用瓦尔德不等式.基于推论3,可得到平均首达时间E(Tτ)的通用定理,如定理3所示. 为单调递增可积函数,若满足 则对Tτ, 成立. 证明令 当 时,存在一个 使得 即 假定 由推论3得 证毕. 类似地,可以证明E(Tτ|d(ξ0))上界成立. 本节运用更新理论模型,分析(1,λ)ES在倾斜平面和超平面问题上的平均首达时间上界和下界.(1,λ)ES带均匀变异算子用于解决最优化问题算法流程如算法2. 算法2(1,λ)ES流程 输入初始种群个体ξ0,步长Δ0(σ0) 令t=0; 计算每个个体适应值,评估当代种群中个体优劣; while not stop do fori=1∶λdo ξ′t=ξ′t-1+Δt(σ0)U; end for 计算每个个体适应值,评估当代种群中个体优劣; 从λ个个体中选择一个ξ′t个体作为下一代的种群; t=t+1; end while 输出到目前为止找到的满意解. 在(1,λ)ES算法中,每代以Δt(σ0)为步长,通过 ξ′t=ξ′t-1+Δt(σ0)U 生成λ个独立同分布子代个体,U表示均匀变异算子,满足均匀分布,概率密度函数为: 选择操作采取非精英保存,即每代维持一个父代个体.通过计算每个个体的适应值,评估种群个体适应值,从λ个子代个体中选择1个个体作为下一代的种群.如果(1,λ)ES算法收敛,即1/5变异子代个体适应值比父代个体适应值更小,则Δt(σ0)减小;否则Δt(σ0)增大.然而,在实际运行中,算法很难找到最优解,为了便于理论分析,假设步长Δt(σ)为一个定值1. 均匀分布变异是进化算法中常用的变异算子[21-22],这里使用均匀分布变异在倾斜平面和超平面上进行首达时间的分析. 本节主要分析(1,λ)ES带均匀变异算子在倾斜平面问题上的平均首达时间分析.倾斜平面问题是一个连续搜索空间的优化问题,为了降低计算量,本节主要讨论二维倾斜平面问题: f(y1,y2)=y2, (y1,y2)⊂S=[-c,c]×[0,c]. 由此可知,y2=b时,二维倾斜平面取得最大值.由于本文考虑最小值问题,因此,讨论的二维倾斜平面问题转化为 f=minf(y1,y2)=-y2, (y1,y2)⊂S=[-c,c]×[0,c]. 可见,y2=-c时二维倾斜平面取得最小值,即目标适应值f*=-c.由于二维倾斜平面问题只在y2上有分量,可以归结为一元问题.定义 为第t-1次迭代到最优解的适应值差,第t-1代到第t代的增长率为: (1,λ)ES的平均首达时间关键是推导增长率的概率密度函数,下面将推导带均匀变异(1,λ)ES的增长率分布函数. 定理4对于(1,λ)ES,如果变异算子服从U(-1,1),增长率δt-1的分布函数为: 证明由算法2可知,父代y2分量个体ξ′t-1通过均匀变异Uk(k=1,2,…,λ)产生λ个子代个体 ξ′t,k=(ξ′1,ξ′2,…,ξ′λ-1,ξ′λ),k=1,2,…,λ, 其个体之间相互独立,则 根据增长率定义,得 进而增长率δt的分布函数F(γ)表示如下. 1)当γ<0时, F(γ)=P(δt<γ)=0. 2)当γ=0时, 3)当0<γ≤1时, F(γ)=P(δt<γ)= P(maxUk<γ)= (P(δt<0)+P(δt=0)+P(0<γ<1))λ= 4)当γ>1时, F(γ)=1. 因此 证毕. 定理5在(1,λ)ES求解倾斜平面问题中,对于平均首达时间Tτ,有 成立. 证明根据定理4,可得 又因为 所以可得出算法平均增长率: 根据定理3,对于 可得 E(Tτ|d(ξ0))≤ 即 因此,估计平均首达时间的上下界: 证毕. 本节分析(1,λ)ES带均匀算子在超平面上的平均首达时间. 如果球函数的对称性在一个方向上被破坏,这就是RIDGE问题,相应的RIDGE函数表达式如下: 不失一般性,本节主要研究超平面问题,这是RIDGE函数中θ=0时的简化形式[31-32].相应的函数表达式为 其中x=(x1,x2,…,xn),考虑当x=(0,0,…,0),f(x)=a0为理想最小值. 定义 为第t-1次迭代到最优解的适应值差.那么第t-1代到第t代的增长率为: 虽然独立同分布的均匀分布具有可加性,但是计算量较复杂,积分不易计算.为了减少冗余的复杂计算,下面主要考虑n=5的情形,结合列维-林德伯格中心极限定理进行分析.增长率δt由变异算子及问题规模决定,因此增长率δt是一个随机变量,分布函数计算如定理6. 定理6设n=5时,对于(1,λ)ES,如果变异算子U服从U(-1,1),那么增长率δt的分布函数为: 则随机变量 其中~表示近似等于. 1)当γ<0时,即f(ξ′t)>f(ξ′t-1),则 F(γ)=0. 2)当γ=0时,即 则 3)当γ>0,即f(ξ′t) 时,则 因此 证毕. 不妨假定进化策略在初始化种群时,初始化种群到满意解的距离为d(ξ0),即有如下定理7成立. 成立. 证明根据定理4,可得 因此 假定算法初始化从原点出发,根据定理4,推导的首达时间Tτ下界 上界 证毕. 由定理5和定理7可得,(1,λ)ES在满足均匀分布变异算子时,理论分析可得到倾斜平面和超平面的首达时间的上下界. 本节实验目的是为了验证更新理论模型在分析(1,λ)ES在倾斜平面问题与超平面问题平均首达时间理论上下界的准确性. 定理5给出二维倾斜平面的平均首达时间上下界的闭合表达式,下面通过实验方式验证此理论结果.实验中设置如下实验参数: 固定误差 记理论分析下界 理论分析上界 表1给出(1,λ)ES在二维倾斜平面上实际平均首达时间与理论上的时间上界和下界.从表中可以看出,在选择种群小于λ=13时,随着种群规模的不断增大,平均首达时间不断减小,与理论分析时间的上下界一致.然而当种群大于λ=13时,随着种群的不断增大,平均首达时间和理论分析时间的上下界的结果并未减小,反而趋于一个稳定值.此结果表明种群规模与平均首达时间之间并非负相关. 表1 二维倾斜平面上理论与实际平均首达时间对比Table 1 Comparison of theoretical and actual first hitting time on a two-dimensional inclined plane 文献[24]引入鞅论和勒贝格控制收敛定理,分析平均首达时间的上界,得到紧致的闭合表达式,但是勒贝格控制收敛定理要求存在控制函数,这是一个很强的条件,而本文模型无需要求存在控制函数,削弱条件后得到更一般结果的时间上下界.文献[23]只给出一个渐近的表达式结果,而本文模型考虑种群规模对首达时间的影响,得出更紧的时间下界. 定理7给出超平面的平均首达时间的上下界的闭合表达式,下面验证此理论结果.实验中设置如下实验参数: 固定误差 记理论分析下界为: 理论分析上界为: 表2给出(1,λ)ES在超平面上实际平均首达时间与理论上的时间上下界.从表中可以看出,随着种群规模λ不断增长,实际运行平均首达时间不断减少,这与理论分析的结果一致.文献[30]通过加性漂移(Additive Drift)分析(1+1)EA在超平面上的计算时间上界,得到紧致的闭合表达式,但未关注首达时间下界.同时本文模型将(1+1)EA推广到(1,λ)ES,得到更通用的结果. 表2 五维超平面上理论与实际平均首达时间对比Table 2 Comparison of theoretical and actual first hitting time on a five-dimensional hyperplane 本文针对进化策略平均首达时间上下界估计困难的问题,提出更新理论模型予以解决.在倾斜平面和超平面问题上,分析(1,λ)ES带均匀分布变异的平均首达时间,得到估计进化策略首达时间的上下界闭合解析式,同时分析种群规模与平均首达时间之间的相关性.更新理论模型在连续搜索空间上呈现逐渐收敛的过程,并且在种群规模一定范围内,增加种群规模能减少进化策略的平均首达时间.数值实验表明,进化策略平均首达时间的理论估计与实际运行的首达时间一致.与此同时,本文方法可应用到其它进化算法,如微搜索算法、粒子群算法和差分进化算法等.今后将基于增长率,使用统计学的方法辅助理论方法估计进化策略的首达时间,尤其是带重组的计算时间分析.另外,本文假设进化策略变异步长是定值情况,因此,下一步将对可变步长展开分析,把理论结果推广到进化策略的实际应用中.3 进化策略的平均首达时间分析
3.1 (1,λ)ES在倾斜平面问题上的首达时间分析
3.2 (1,λ)ES在超平面问题上的首达时间分析
4 数值实验及结果分析
4.1 (1,λ)ES求解倾斜平面问题
4.2 (1,λ)ES求解超平面问题
5 结 束 语