改进萤火虫算法及其收敛性分析

张大力, 夏红伟, 张朝兴, 马广程, 王常虹,*

(1. 哈尔滨工业大学航天学院, 黑龙江 哈尔滨 150001; 2. 上海航天控制技术研究所, 上海 201109)

0 引 言

社会性生物的某些集体行为往往呈现出寻优的本质,通过对这类行为的公式化,得到了一类被称为群智能算法的优化算法,蚁群优化(ant colony optimization, ACO)算法、粒子群优化(particle swarm optimization, PSO)算法、人工蜂群(artificial bee colony, ABC)算法等均属此列。由于算法简单灵活、鲁棒性好,可用于求解许多传统方法无法解决的NP-hard问题,如旅行商问题、装箱问题等,已经广泛应用于人工智能、通信网络和工业生产等领域。

萤火虫算法(firefly algorithm, FA)是英国学者Yang于2008年提出的一种元启发式算法。该算法通过模拟自然界中萤火虫个体之间相互吸引的理想化行为达到寻优目的,一经提出便受到国内外研究人员的关注,广泛应用于计算机科学、复杂方程求解、结构和工程优化领域。但许多构造的和现实的优化问题具有大量的局部最优解,容易使算法因种群多样性的快速衰减而陷入早熟停滞,大幅降低优化效果。同时,大量局部最优对种群的吸引可能导致算法不断震荡而无法收敛。有学者通过引入自适应或随机机制针对上述问题进行改进。Gandomi等采用12种不同的混沌映射对参数进行自适应调节,其中高斯映射调节吸引度参数的改进算法取得了一定效果,但性能提升有限,寻优成功率不高。Yu等设计了一种自适应步长策略来提高算法的全局探索能力,每个个体的随机步长根据当前和历史信息分别计算,然而实际寻优效果与标准算法区别很小。Altay等将该改进算法用于调节数字水印技术中的比例因子,由于缺乏与标准算法的横向对比,不足以说明改进算法的优越性。Verma等在算法初始化阶段引入反向学习策略提高算法收敛速度,在迭代过程中采用基于维度的位置更新方法跳出局部最优,但与其他优化算法相比没有体现出明显的竞争力。Manju等参考量子PSO提出了量子FA,并利用基于二次插值的重组算子改进搜索过程,提高了算法的收敛速度和寻优能力,但无法处理复杂病态优化问题。柳长源等引入了一种驱散机制,若种群最优值在连续多次迭代后未发生变化,则对最优个体采用柯西分布进行变异,聚集在最优个体附近的个体采用均匀分布进行驱散,从函数优化结果来看改进效果并不显著。刘景森等提出了一种具有振荡、约束和自然选择机制的改进策略,提高了算法对高维问题的求解能力,但是增加了额外的调节参数,不利于推广应用,其收敛性证明仿照PSO算法证明的思路,简化到二维空间进行分析,无法向高维情况推广。部分学者选择从个体信息交互的角度对算法进行改进。Wang等提出了一种随机吸引模型,通过将个体与其邻域内有限个其他个体进行通信,降低了算法计算复杂度,但寻优性能没有明显提升。在此基础上,Dhal等进一步考虑粗糙集和局部搜索机制,提出了一种随机吸引粗糙FA,并成功应用于图像分割过程,但与标准算法及其他传统优化算法相比没有明显的竞争力。Shan等通过分布式并行技术将种群划分为多个子群,子群之间基于4种不同的通信策略进行信息交换,但每个子群的寻优能力偏弱,不能保证整体的寻优效果,且增加了额外的配置参数,适用性不强。Elakkiya等提出了一种基于社区启发的改进策略,个体只被由个较优个体组成的所谓社区中的个体吸引,并给出了社区规模的自适应更新机制,降低了计算复杂度,但无法保证寻优精度。还有一种改进策略是在FA中引入其他算法的思想。Zhao等借鉴ABC算法的思想将萤火虫种群分为领导者、开发者和跟随者3种角色,并为每个角色分配了一种学习策略,用来平衡算法的探索和挖掘能力,但没有取得显著的性能改进。Xing等为了使算法跳出局部最优,引入了遗传算法中的变异算子,同时仿照PSO算法将全局最优纳入到位置更新公式中,改进算法用于解决协同干扰资源分配问题,寻优结果与标准算法相比有少许提升。薛晗等结合文化算法和FA的优点,将迭代过程中获取的经验和知识组成信仰空间,并设计了一套规则来影响种群更新,但算法没有足够的跳出局部最优的能力。Gupta等将引力搜索和FA通过串行运行的方式结合起来,提出了一种混合引力搜索FA,计算复杂度较高,寻优能力没有明显提升。此外,Cheng等和Xue等采用混合编码的形式对个体进行编码,利用组合变异、协同进化等方式加快算法收敛速度,但该方法只适用于解决离散工程优化问题,缺乏适用性。

基于以上分析可以看出,目前FA仍然具有较大的改进潜力,并且理论分析比较薄弱,也在一定程度上制约了FA的进一步改进和发展。因此,本文提出一种改进FA(improved FA, IFA),以提升算法性能。首先,介绍FA的基本原理。其次,给出新的个体位置更新机制。围绕种群更新过程,引入位置置换变异策略和差分进化(differential evolution,DE)算法中的最优变异策略(DE/best/1),并给出算法流程。然后,给出改进算法收敛性的理论证明。最后,采用一组典型基准函数和装箱问题对改进算法有效性进行验证。

1 FA基本原理

FA受萤火虫群体的集聚性启发,其仿生原理如下:搜索空间中的单个元素代表萤火虫个体;拟优化的目标抽象为个体所在位置的亮度,衡量个体位置的优劣并决定其移动方向;元素搜索和更新的寻优过程模拟个体间的吸引和移动行为。

算法数学描述如下:

(1) 个体的位置记为=(1,2,…,),其亮度即为该位置对应的目标函数值()。

(2) 个体与间的吸引度更新公式定义为

(1)

式中:为初始吸引度,一般取值为1;为光强吸收系数,一般取值为1;为个体间的笛卡尔距离,有

(2)

(3) 个体被吸引并向其移动的位置更新公式为

(3)

式中:为算法迭代次数;为步长因子;为随机数,服从均匀分布或高斯分布。

2 IFA

威廉·奥卡姆曾言,“如无必要,勿增实体”,这就是著名的“奥卡姆剃刀定律”。从优化算法设计的角度来诠释,可理解为在效果相同或相近的前提下,尽量选择简单框架和规则。在这种思想的指导下,本文提出以下改进策略。

2.1 随机扰动策略

图1 吸引力项运动特征Fig.1 Motion characteristic of attraction term

由于式(3)第3项随机项的变化是独立的,没有利用个体和种群信息,因此其变化很难与个体的进化相匹配。一方面,如果取值过小,那么较小的搜索半径无法提供跳出局部最优的能力,如果取值较大,则会因为较大的搜索半径导致较慢的收敛速度和较低的收敛精度。另一方面,一旦在随机项中引入个体信息,那么可与改进后的第2项合并同类项,在功能上被其覆盖。

基于以上考虑,向原位置更新式(3)中引入一个随机扰动因子,使得吸引力项拥有一个随机权重,并去掉式(3)第3项,有

(4)

为保证式(4)在初始阶段的搜索能力,将吸引度更新为

(5)

式中:为吸引度调节尺度,根据调试经验建议取值范围为07～09。

随机扰动因子的引入提高了个体的搜索范围。一方面,在算法初始阶段,个体间距离较大,式(5)表示的随机权重简化为·,此时搜索步长较长,算法能够在全局范围内进行寻优,一定程度上提高了算法的全局探索能力;另一方面,在算法迭代后期,个体间距离大幅减小,随机权重简化为(+)·,而此时搜索步长较小,算法能够在当前位置的邻域深度挖掘。个体搜索范围的提升使得种群在运动过程中能够覆盖到更多搜索空间中的区域,如图2所示。

图2 种群运动趋势和个体分布情况Fig.2 Population movement trend and individual distribution

下面对随机扰动因子的选择进行分析。典型的随机分布有均匀分布、正态分布、柯西分布、列维分布等。为了更直观地区分不同分布的特征,给出上述分布100 000个样本的分布直方图,如图3所示。

图3 4种不同概率分布的100 000样本分布Fig.3 Four different probability distributions of 100 000 random samples

由图3可以看出,均匀分布产生的样本分布均匀且取值区间较小,选择均匀分布作为随机扰动因子可获得在各个方向上更加稳定的搜索能力。正态分布产生的样本呈现出左右对称、中间高两侧低的钟型分布特征,在原点附近取值的概率较高,离原点越远取值概率越低,选择正态分布作为随机扰动因子,存在初始阶段收敛慢的情况,不能保证算法收敛速度。柯西分布和列维分布产生的样本呈现出明显的厚尾特征,在原点附近取值的概率较高,偶尔会出现非常大的样本,选择柯西或列维分布作为随机扰动因子,算法搜索能力不稳定,且产生的大样本很有可能使个体更新超出搜索空间,精度和收敛速度均无法保证。基于以上分析,随机扰动因子首选均匀分布,在第4.1节中给出仿真，对以上分析做进一步说明。

2.2 双重变异策略

为进一步提升算法寻优性能,本文围绕个体和种群的更新过程增加两个变异机制。

221 位置置换变异策略

图4 位置置换变异原理和种群更新流程Fig.4 Principle of position substitution mutation and population renewal process

参数的取值是位置置换变异策略的关键。一方面,由于该操作在个体位置更新后进行,且嵌套于种群更新过程中,因此不能过大,否则会提高算法计算复杂度,降低收敛速度。另一方面,若次数太少,个体位置变化较小,无法达到跳出局部最优的目的,降低寻优精度。因此,置换变异至少2次,至多建议不超过5次,本文取=3。

222 DE/best/1变异策略

(6)

图5 DE/best/1/变异种群更新流程
Fig.5 DE/best/1 cross mutation population renewal process

DE/best/1变异操作的执行次数为。由于该部分没有循环嵌套,可以取大一些的数,从而充分利用历史最优信息指导搜索过程,提高算法收敛精度。参考差分进化算法的取值,可取种群规模。

两种变异策略得到的子代种群采取与父代种群竞争的模式选取进入下一轮循环的种群,在每轮迭代中保留较优个体,提升整体种群质量。

2.3 算法平衡性分析

全局探索与局部挖掘能力是群智能优化算法解决问题的核心,二者之间的平衡对算法性能至关重要。研究普遍认为,在算法开始阶段应以全局探索为主,在后期应以局部挖掘为主。改进算法能够很好地满足此规律。在算法迭代初始阶段,个体间相距较远,从位置更新式(4)和DE/best/1变异策略式(6)可直观地看出,此时算法有较大的搜索步长,能够在解空间内进行全局探索。算法迭代后期,种群不断聚集,搜索步长随之减小,此时算法的局部挖掘能力得到不断强化。第221节所述的位置置换变异策略是对个体自身信息的突变,该操作不受迭代过程的影响,无论在迭代前期还是迭代后期均能使个体有效地跳出当前位置,从而避免算法陷入局部最优。需要注意的是,在实际运行中,全局探索和局部挖掘能力是算法各个策略共同触发的,没有明确的界限。

2.4 算法流程

下面给出改进算法的流程。

初始化控制参数，，，位置变异次数为，DE/best/1变异执行次，种群规模为，搜索空间下界为，上界为，最大迭代次数为MaxGeneration。

构建服从均匀分布的初始种群，计算适应度和当前全局最优解g。

取种群中的一个个体，判断其与其他个体的吸引关系，根据式(2)、式(4)、式(5)更新该个体的位置并进行限位，同时更新全局最优值g，该步骤执行次。

根据第221节给出的策略对步骤3中选定的个体进行位置置换变异，该步骤执行次。

步骤3～步骤4连续执行，重复次数为次，得到子种群1。

根据第222节给出的策略对种群1进行变异操作，同时更新全局最优值g。该步骤执行次，得到子种群2。

父代种群和子代种群拼接并排序后截断保留 维数据，得到下一轮迭代的父代种群。

本文设置停止准则为迭代次数达到MaxGeneration，若满足条件，则终止判断，输出优化结果并可视化；若不满足条件，则重复步骤3～步骤7。

2.5 算法计算复杂度分析

结合算法流程进行分析,为了简便这里记为迭代次数,为适应度函数的计算量,那么算法语句总的执行次数可以表示为·[(5+)·+·(+4)·+(+7)·],此时改进算法的时间复杂度为(··)。显然,算法的计算复杂度与问题本身的复杂度直接相关,主要由适应度函数计算量和算法计算量决定。

3 收敛性证明

3.1 随机算法的收敛准则

Wets给出了一般随机优化算法收敛性判别准则,下面重述其主要结论。

假设有随机优化问题[,],其迭代方式为,第次迭代结果,则+1=(,)。其中,为适应度函数,为解空间,是次迭代产生的所有解。

在Lebesgue测度空间定义搜索的下确界:

=inf(∶[∈|()<]>0)

(7)

式中:[]表示集合上的Lebesgue测度。根据式(7)定义最优解区域为

(8)

式中:>0且足够小,>0且足够大。若算法找到了,中的一个点,则算法收敛。

下面给出算法收敛的充分条件:

全局搜索收敛

(9)

(10)

(11)

证毕

3.2 改进算法基本概念的数学描述

萤火虫位置状态和位置状态空间

萤火虫位置的状态=(,),由萤火虫的位置和全局最优值构成,所有可能的个体位置状态的集合称为萤火虫的位置状态空间,记为={=(,)|,∈,()<()}。其中,萤火虫个体在一次迭代过程中产生的两组子种群的个体位置记为1和2,其位置状态分别记为1=(1,)和2=(2,),为包含两子种群个体的位置状态空间。

群体状态和群体状态空间

种群位置状态集合称为群体状态,记为={,,…,},所有群体状态构成的群体状态空间记为={=(,,…,),∈,1≤≤},其中子群位置状态记为1=(11,12,…,1)和2=(21,22,…,2),且包含两子群状态空间。

3.3 IFA的马尔可夫模型建立

对∀=(,gb)∈,∀1=(1,gb)∈位置状态由转移到1记为1()=1。

萤火虫的位置状态由转移到1的转移概率为

(1 ()=1)=(→1′)·(gb→gb′)·(1′→1)·(gb′→gb″)

(12)

式中:(→1′)为式(4)对应的转移概率;(gb→gb′)为该步骤全局最优位置的转移概率;(1′→1)为步骤4对应的转移概率;(gb′→gb″)为该步骤更新全局最优的转移概率。

萤火虫位置状态由转移到1,意味着→1′,gb→gb′,1′→1,gb′→gb同时成立。作为子种群,有

(→1′)=1

(13)

全局最优位置gb的转移概率为

(14)

位置状态由1′→1的概率为

(15)

全局最优位置的再次转移概率为

(16)

证毕

对∀=(1,2,…,)∈,∀1=(11,12,…,1)∈,群体状态由转移到1,记作1()=1。

萤火虫位置群体状态由转移到1的转移概率为

(17)

若萤火虫位置群体状态由转移到1,那么表示中所有个体的位置状态同时转移到1中所有个体的位置状态,即1(1)=11,1(2)=12,…,1()=1同时成立。

证毕

在IFA算法迭代中,对∀1=(11,12,…,1)∈,∀2=(21,22,…,2)∈,记作2(1)=2。

萤火虫位置群体状态由1转移到2的转移概率为

(18)

若萤火虫位置群体状态由1转移到2,那么有1中所有个体位置状态同时转移到2中所有个体位置状态,即2(1)=21,2(1)=22,…,2(1)=2同时成立。在迭代过程中,第二次变异操作个体的状态转移与种群状态转移有关,根据式(7)有

(19)

证毕

对∀=(1,2,…,)∈,=(1,2,…,)∈,将萤火虫位置的群体状态转移到记作()=。

萤火虫位置群体状态由转移到的转移概率为

(20)

式中:为对父代种群和两子种群进行排序和截断操作,(gb‴→gb)为寻找种群中的最优个体的转移概率。

若萤火虫位置的群体状态由一步转移到,就意味着位置的群体状态发生了集体转移,由定理1～定理3可推出最终的转移概率公式。

证毕

萤火虫的位置群体状态序列{();≥0}是离散时空上的有限齐次马尔可夫链。

首先,对于任何一个优化算法来讲,其搜索空间都是有限的,任意个体状态=(,g)中的,g也是有限的,所以位置的状态空间是有限空间。又由于=(,,…,)由个个体组成,为大小有限的正整数,故萤火虫位置的群体状态也是有限的。

其次,∀(-1)∈,∀()∈,其状态转移概率(((-1))=())是由群体中所有个体的转移概率和两子群的转移概率决定。由定理1～定理3,可得

(1((-1))=1(-1))=((-1)→1′(-1))·(g(-1)→g′(-1))·(1′(-1)→1(-1))·(g′(-1)→g(-1))

(21)

(22)

(23)

则

(24)

通过以上分析,(((-1))=())仅与-1时刻的位置状态有关,故萤火虫位置的群体状态{();≥0}具有马尔可夫性,且算法迭代与时间无关,所以该过程是齐次的。

证毕

3.4 IFA的收敛性分析

定义最小化问题∶→,⊆

min() s.t.∈

式中:={,,…,}为维决策变量,所有满足约束的值构成的集合称为可行域。位置最优位置状态集为,其群体状态集为={=(,,…,)|∃∈,1≤≤}。

对闭集,∀∈,有∑∈=1,意味着从状态空间内的任何一点出发,始终无法达到以外的任意状态。

对于个体位置状态集{();>0},集合是空间上的一个闭集;对位置群体状态集{();≥0},集合是空间上的一个闭集。

同理,假如∀∈,∃∉,→成立,那么至少存在一个个体的位置由的内部转移至外部,又已证明是上的一个闭集,根据式(20),有(()=)=0,由闭集定义知,集合是空间上的一个闭集。

证毕

在萤火虫位置群体状态空间中,无法找到外的非空闭集,使其满足∩=。

假设在位置群体空间中能够找到一个闭集,满足∩=,那么∃∈,∀∈,有(gb)>(gb)。根据式(20),其满足各个组件不为0的条件,那么存在转移概率(()=)≠0,即不是闭集。定理得证。

证毕

假设〈〉是一组服从rand(0,1)的随机变量,且满足独立同分布,那么

随着萤火虫位置群体状态不断迭代更新,当次数趋于无穷时,萤火虫位置群体状态序列必然进入最优集。

由定理6和定理7以及引理1和引理2得证。

在迭代次数趋于无穷时,IFA能够以概率1收敛于全局最优。

一方面,在每次迭代过程中,IFA都保留了全局最优,从而保证了适应度的非增性,满足收敛条件一;另一方面,根据定理9,当迭代次数趋于无穷时,萤火虫位置的群体序列必然进入最优集,即找不到全局最优的概率为0,满足收敛条件二。故IFA以概率1收敛于全局最优。

证毕

4 仿真实验

由于实际应用中算法运行时间有限,故需通过固定时间的仿真实验对算法有效性进行验证。算法运行环境为Intel Core i5-3470处理器,CPU3.20 GHz,3.48 GB内存,Windows10,64位操作系统。

4.1 随机扰动因子性能对比

为验证第2.1节中的结论,对4种随机扰动因子下的改进算法RIFA, GIFA, CIFA, LIFA进行仿真试验和性能对比,R(平均)、G(正态)、C(Cauchy)、L(Levy)分别对应改进算法采用的随机分布。选用4个典型基准函数,如表1所示。

算法参数设置情况为=1,=08,=1,=3,=40,取种群规模,=30,～U(04,09,1,),最大迭代次数为2 000。初始种群由问题域内进行随机初始化产生。为提高结果的可信度,保证对比实验客观公正,对每一个测试,所有算法均独立执行50次,取平均最优值、时间和成功率(success rate, SR)的统计结果进行比较,最优结果由粗体标出,仿真结果如表2所示。

表1 基准函数Table1 Benchmarks functions

表2 4种随机扰动因子仿真结果Table 2 Four random disturbance factors simulation results

表2的统计结果显示,随机扰动因子选择均匀分布时,在寻优精度、收敛时间和SR上的表现最好,对所有基准函数均能找到最优解;选择正态分布时,算法收敛时间大幅提高,寻优精度和SR也出现了下降,特别是函数～的表现很差;柯西扰动因子的表现与正态分布类似,收敛时间比前者稍短,但寻优精度更加无法保证;列维扰动因子表现最差,运行效率最低,优化结果几乎不可用。因此,改进算法的随机扰动因子优先选择均匀分布。

4.2 与其他算法的性能对比

选取已在实践中检验过的优秀算法PSO、DE以及FA进行对比试验。各个算法的参数设置情况如下:

(1) PSO:惯性权重=1,阻尼比=0.99,加速因子=15,=2.0。

(2) DE(DE/rand/1):缩放因子F为范围在0.2～0.8的均匀随机数,变异概率=0.2。

(3) FA:=02,=0.98,=2,=1。

(4) IFA参数同第4.1节。

其他仿真条件同第4.1节,仿真结果见表3。

表3 固定迭代次数仿真结果Table 3 Simulation results of fixed iterations

表3的统计结果显示,对于函数～,改进算法均能够找到最优解、精度和SR均好于对比方法,在和尤为明显。但在平均仿真时间上,改进算法比其他算法长。这是由于PSO、DE和FA过早陷入局部最优,根据优化算法本身运行机制,如果没有足够的跳出局部最优的能力,那么会不断加速收敛,因此其运行时间短于改进算法。

函数是典型的复杂病态二次函数,其全局最优位于狭窄的、抛物形的山谷中,缺乏有价值的信息,优化难度较大,比较考验算法全局探索和局部挖掘能力。各方法对该函数进行优化的寻优曲线如图6所示。

图6 不同算法性能比较Fig.6 Performance comparison of different algorithms

由图6可以看出,PSO、DE和FA均出现了早熟收敛现象,而IFA在短暂收敛后,在迭代到151步时成功跳出局部最优并最终达到极高的优化精度。现给出算法迭代150步到151步时的个体分布情况。由于个体维数较高,因此通过评估每个个体与平均个体在搜索空间中的相对距离得到个体分布散点图,如图7所示。

图7 典型迭代时刻个体分布Fig.7 Individual distribution of typical iteration time

可以看出,种群在151步时部分个体突变至距离平均个体较远的位置,获得了更优解的同时,带动种群整体移动,整体的多样性比150步时有了明显的提升。通过以上仿真和分析,证明了所提改进策略在基准函数优化方面的有效性。

4.3 装箱问题实例

本文将IFA用于求解装箱问题。装箱问题是一个经典NP-hard组合优化问题,在运筹学、离散数学、计算机科学以及工业生产等方面有着广泛的应用。

设有足够数量的结构相同容量一致的箱子,,…,每个箱子的最大容量为。现有一组个物品,,…,，体积分别为(0<<,=1,2,…,)，需装入箱内,要求在满足约束的情况下打包所有物品,要求所需箱子数量尽可能的少。两个约束分别为:① 所有物体均不可分割,需保持完整并恰好放入一个箱子中;② 每个箱子所放物品体积总和均不能超过容量上限。假设物品形状对装箱没有影响,则该问题可描述为

(25)

其中,

假设30个物品体积分别为[6,3,4,6,8,7,4,7,7,5,5,6,7,7,6,4,8,7,8,8,2,3,4,5,6,5,5,7,7,12],箱子容量=30。固定迭代1 000步,种群规模为20,其余参数不变。优化结果与PSO、DE和FA比较,算法收敛情况如图8所示。

图8 装箱问题的性能比较Fig.8 Comparison on bin packing problem

由图8可以看出,IFA仅需12步即可达到其他3种算法的优化效果,并且经进一步迭代后能够跳出局部最优,从而得到更优的结果。表4给出了几种优化算法的总体运行时间(total time,TT)、收敛时间(convergence time,CT)和最优结果(Best)。比较可知,IFA在优化结果和收敛速度方面都有着明显的优势。

表4 装箱问题优化结果Table 4 Optimization results of bin packing problem

5 结 论

本文为克服FA易陷入局部最优的缺陷,提出了IFA。首先分析了个体位置更新机制的特点,提出了一种随机扰动策略,阐述了改进策略的有效性,并分析了不同扰动因子选取的可行性。通过引入元素置换变异策略和最优差分变异策略,提高了算法跳出局部最优的能力。基于随机算法收敛准则和马尔可夫过程的理论分析指出,改进算法在足够大的迭代次数下能够以概率1收敛到全局最优。对4种基准函数的仿真结果表明,随机扰动因子选择均匀分布时有最佳表现。在相同条件下与几种经典算法相比,改进算法的收敛精度和SR更高,验证了改进算法的优越性。将该方法对经典装箱问题进行求解,得到了较好的优化结果,验证了改进方法的有效性和实用价值。