JADE: Adaptive Differential Evolution with Optional External Archive

摘要(Abstract)

本文提出了一种新的差分进化(DE)算法——JADE,旨在通过实现一种新的变异策略 “DE/current-to-pbest” 并结合可选的外部归档,以及以自适应方式更新控制参数,来提升优化性能。该“DE/current-to-pbest”策略是经典“DE/current-to-best”的泛化形式,而可选归档操作则利用历史数据来提供进化方向的信息。这两种操作均能增加种群多样性并改善收敛性能。参数自适应机制能自动将控制参数更新至合适数值,避免了用户需要预先了解参数设置与优化问题特性之间关系的先验知识,从而有助于提高算法的鲁棒性。仿真结果表明,对于一组20个基准问题,JADE在收敛性能上优于或至少相当于其他经典或自适应DE算法、经典粒子群优化(PSO)以及文献中的其他进化算法。带有外部归档的JADE在处理相对高维问题时显示出良好的结果。此外,本文结果明确表明,不存在适用于各种问题、甚至同一问题不同优化阶段的固定控制参数设置。

关键词 :自适应参数控制,差分进化,进化优化,外部归档。


I. 引言(INTRODUCTION)

差分进化(DE)已被证明是一种简单而高效的进化算法,适用于许多现实世界中的优化问题[1]-[5]。然而,根据实验研究[6]和理论分析[7],其性能在很大程度上仍依赖于控制参数(如变异因子和交叉概率)的设置。尽管已有关于参数设置的建议[2], [6], [8],但参数设置与优化性能之间的相互作用仍然复杂且未被完全理解。这主要是因为不存在一种固定的参数设置能够适用于各种不同的问题,甚至无法适用于同一问题在不同进化阶段的需求。

即使对于在整个进化搜索过程中固定其参数的算法(例如经典的 DE/rand/1/bin),通过试错法来调整控制参数通常也需要繁琐的优化试验。受此考虑启发,近年来引入了不同的自适应或自适应机制[9]-[18],以动态更新控制参数,而无需用户事先了解参数设置与优化问题特性之间的关系。此外,如果设计得当,参数自适应能够提高算法的收敛性能。

上述参考文献中提出的参数自适应机制可以根据控制参数的改变方式进行分类。根据 Angeline [19] 和 Eiben 等人 [20], [21] 引入的分类方案,我们可以将参数控制机制解释为以下三类:

  1. 确定性参数控制:控制参数通过某些确定性规则进行改变,而不考虑进化搜索的任何反馈。一个例子是 Holland [22] 提出的随时间变化的变异率。

  2. 自适应参数控制:利用来自进化搜索的反馈来动态改变控制参数。这类方法的例子包括 Rechenberg 的“1/5 规则”[23] 以及 [11], [12] 中的模糊逻辑自适应 DE。几个新近提出的 DE 算法,如 SaDE [13]、jDE [15]、SaNSDE [18] 以及本文提出的算法,也可归入此类。:其中一些算法在原始论文中被归类为“自适应性”。然而,尽管种群中的每个个体都与其自身的控制参数相关联,但这些参数本身在优化过程中并不经历变异和交叉。因此,根据 Eiben 等人的定义,这些算法可以被视为自适应方案。

  3. 自适应性参数控制采用“进化的进化”方法来进行控制参数的自适应。控制参数直接与个体关联,并经历变异和重组/交叉。由于更好的参数值倾向于产生更有可能存活的个体,这些值可以传播给更多的后代。用于多目标优化(MOO)的 SPDE [9] 和 DESAP [10] 属于此类。

如果设计得当,自适应或自适应性参数控制可以通过将参数动态适应于不同适应度地形(landscapes)的特性来增强算法的鲁棒性。因此,它适用于各种优化问题而无需试错。此外,如果对于特定问题的不同进化阶段,控制参数能够被调整到合适的值,则可以提高收敛速度。

对于许多基准问题,自适应和自适应性 DE 算法已显示出比没有参数控制的经典 DE 算法更快、更可靠的收敛性能[13]–[17]。因此,寻找最合适的 DE 变体并将其作为引入参数自适应操作的基础方案是很有意义的。一些自适应算法[10], [14], [18]是基于经典的 DE/rand/1/bin 开发的,该策略以鲁棒性强但收敛效率较低而闻名。其他算法[13], [14], [16]则同时实现了 DE/rand/1/bin 和一个或多个贪婪的 DE 变体(如 DE/current-to-best/1/bin 和 DE/best/1/bin),并动态更新使用每种变体生成后代的概率。

迄今为止,还没有方法是完全基于利用当前种群中最优解信息的贪婪 DE 变体(如 DE/current-to-best/1/bin 和 DE/best/1/bin)开发的。原因似乎很简单:贪婪变体通常可靠性较差,并可能导致诸如早熟收敛等问题,尤其是在解决多峰问题时[8], [24]。然而,我们注意到,一个设计良好的参数自适应方案通常有利于增强算法的鲁棒性。因此,在自适应 DE 算法中,贪婪 DE 变体的可靠性问题变得不那么关键,而其快速收敛的特性则更具吸引力

除了当前种群中的最优解,历史数据是另一个可用于提高收敛性能的信息来源。一个例子是粒子群优化(PSO)[25],其中利用历史上探索过的最优解来指导当前种群(粒子群)的移动。在本文中,我们感兴趣的不是先前探索的最优解,而是一组近期探索过的较差解,并认为它们与当前种群的差异是指向最优解的有希望的方向。这些解也被用来提高种群多样性。

鉴于以上考虑,我们引入了一种新的贪婪变异策略——“DE/current-to-pbest”及其可选的外部归档,并将其作为我们的参数自适应算法 JADE 的基础。JADE 算法的早期版本已在[26]中发表。作为 DE/current-to-best 的泛化,DE/current-to-pbest 不仅利用最优解的信息,还利用其他优质解的信息。具体来说,在 DE/current-to-pbest 中,可以随机选择前 100p% (p ∈ (0, 1]) 中的任何一个解来扮演在 DE/current-to-best 中专门为单一最优解设计的角色。此外,归档的较差解与当前种群之间的差异可以被纳入变异操作中。尽管具有贪婪特性,所提出的策略能够使种群多样化,从而可以缓解诸如早熟收敛等问题。算法的可靠性通过自适应参数控制得到进一步提高。

论文的其余部分组织如下:第二部分描述了差分进化的基本操作流程,并介绍了有助于在第三部分回顾各种自适应或自适应性 DE 算法的符号和术语。新的算法 JADE 在第四部分详细阐述,并详细解释了带可选归档的 DE/current-to-pbest 和自适应参数控制机制。第五部分给出了仿真结果,用于比较 JADE 与其他进化算法的性能。最后,第六部分总结了结论性意见。


II. DE 的基本操作(BASIC OPERATIONS OF DE)

在本节中,我们将描述差分进化的基本操作,并介绍必要的符号和术语,以便于后续解释不同的自适应 DE 算法。

差分进化遵循进化算法的通用流程。初始种群 {xi,0=(x1,i,0,x2,i,0,…,xD,i,0)∣i=1,2,…,NP}\{x_{i,0} = (x_{1,i,0}, x_{2,i,0}, \ldots, x_{D,i,0})| i = 1, 2, \ldots, NP\}{xi,0=(x1,i,0,x2,i,0,,xD,i,0)i=1,2,,NP} 根据均匀分布随机生成,其范围满足 xjlow≤xj,i,0≤xjupx_{j}^{low} \leq x_{j,i,0} \leq x_{j}^{up}xjlowxj,i,0xjup,其中 j=1,2,…,Dj = 1, 2, \ldots, Dj=1,2,,D,这里 DDD 是问题的维度,NPNPNP 是种群大小。初始化之后,DE 进入一个包含进化操作的循环:变异交叉选择

变异:在每一代 ggg,此操作基于当前的父代种群 {xi,g∣i=1,2,…,NP}\{x_{i,g}| i = 1, 2, \ldots, NP\}{xi,gi=1,2,,NP} 创建变异向量 vi,gv_{i,g}vi,g。以下是文献中常用的几种不同变异策略:

  1. “DE/rand/1”
    vi,g=xr0i,g+Fi⋅(xr1i,g−xr2i,g)(1) v_{i,g} = x_{r_0^i,g} + F_i \cdot (x_{r_1^i,g} - x_{r_2^i,g}) \tag{1} vi,g=xr0i,g+Fi(xr1i,gxr2i,g)(1)

  2. “DE/current-to-best/1”
    vi,g=xi,g+Fi⋅(xbest,g−xi,g)+Fi⋅(xr1i,g−xr2i,g)(2) v_{i,g} = x_{i,g} + F_i \cdot (x_{\text{best},g} - x_{i,g}) + F_i \cdot (x_{r_1^i,g} - x_{r_2^i,g}) \tag{2} vi,g=xi,g+Fi(xbest,gxi,g)+Fi(xr1i,gxr2i,g)(2)

  3. “DE/best/1”
    vi,g=xbest,g+Fi⋅(xr1i,g−xr2i,g)(3) v_{i,g} = x_{\text{best},g} + F_i \cdot (x_{r_1^i,g} - x_{r_2^i,g}) \tag{3} vi,g=xbest,g+Fi(xr1i,gxr2i,g)(3)

其中,索引 r0ir_0^ir0i, r1ir_1^ir1ir2ir_2^ir2i 是从集合 {1,2,…,NP}\{i}\{1, 2, \ldots, NP\}\backslash\{i\}{1,2,,NP}\{i} 中均匀随机选择的不同整数,xr1i,g−xr2i,gx_{r_1^i,g} - x_{r_2^i,g}xr1i,gxr2i,g 是一个差分向量,用于对相应的父向量 xi,gx_{i,g}xi,g 进行变异,xbest,gx_{\text{best},g}xbest,g 是当前第 ggg 代中的最优向量,FiF_iFi 是变异因子,通常取值范围在区间 (0, 1+] 内。在经典 DE 中,Fi=FF_i = FFi=F 是一个固定参数,用于在所有世代生成所有变异向量;而在许多自适应 DE 算法中,每个个体 iii 都关联其自己的变异因子 FiF_iFi

上述变异策略可以通过引入多个差分向量(而不仅仅是 xr1i,g−xr2i,gx_{r_1^i,g} - x_{r_2^i,g}xr1i,gxr2i,g)来进行推广。根据所采用的差分向量数量 kkk,生成的策略被命名为 “DE/–/k”。

需要注意的是,试验向量的某些分量可能会违反预定义的边界约束。解决此问题的方法包括重置方案、罚函数方案等。然而,约束问题并非本文重点。因此,我们遵循一种简单的方法,将违反边界的分量设置为被违反的边界与父代个体相应分量的中点值 [2],即:

vj,i,g={(xjlow+xj,i,g)/2,if uj,i,g<xjlow(xjup+xj,i,g)/2,if uj,i,g>xjup v_{j,i,g} =\begin{cases} (x_{j}^{\text{low}}+x_{j,i,g})/2, & \text{if } u_{j,i,g}<x_{j}^{\text {low}} \\ (x_{j}^{\text{up}}+x_{j,i,g})/2, & \text{if } u_{j,i,g}>x_{j}^{\text {up}} \end{cases} vj,i,g={(xjlow+xj,i,g)/2,(xjup+xj,i,g)/2,if uj,i,g<xjlowif uj,i,g>xjup

其中 uj,i,gu_{j,i,g}uj,i,gxj,i,gx_{j,i,g}xj,i,g 分别表示第 ggg 代时试验向量 ui,g\mathbf{u}_{i,g}ui,g 和父向量 xi,g\mathbf{x}_{i,g}xi,g 的第 jjj 个分量。当最优解位于边界附近或边界上时,此方法尤其有效。

交叉:变异之后,一个二项式交叉操作形成最终的试验/后代向量 ui,g=(u1,i,g,u2,i,g,…,uD,i,g)\mathbf{u}_{i,g}=(u_{1,i,g},u_{2,i,g},\ldots,u_{D,i,g})ui,g=(u1,i,g,u2,i,g,,uD,i,g)

uj,i,g={vj,i,g,if rand(0,1)≤CRi or j=jrand,xj,i,g,otherwise(4) u_{j,i,g}=\begin{cases} v_{j,i,g}, & \text{if rand}(0,1)\leq CR_{i}\text{ or } j=j_{\text{rand}},\\ x_{j,i,g}, & \text{otherwise} \end{cases} \tag{4} uj,i,g={vj,i,g,xj,i,g,if rand(0,1)CRi or j=jrand,otherwise(4)

其中 rand(0,1) 是在区间 [0,1] 上均匀分布的随机数,且为每个 jjj 和每个 iii 独立生成;jrand=randint(1,D)j_{rand}=\text{randint}(1,D)jrand=randint(1,D) 是一个从 1 到 DDD 之间随机选择的整数,且为每个 iii 重新生成;交叉概率 CRi∈[0,1]CR_{i}\in[0,1]CRi[0,1] 大致对应于从变异向量继承的向量分量的平均比例。在经典 DE 中,CRi=CRCR_{i}=CRCRi=CR 是一个固定参数,用于在所有世代生成所有试验向量;而在许多自适应 DE 算法中,每个个体 iii 都关联其自己的交叉概率 CRiCR_{i}CRi

选择:选择操作根据其适应度值 f(⋅)f(\cdot)f(),从父向量 xi,g\mathbf{x}_{i,g}xi,g 和试验向量 ui,g\mathbf{u}_{i,g}ui,g 中选择较好的一个。例如,对于一个最小化问题,被选中的向量由下式给出:

xi,g+1={ui,g,if f(ui,g)≤f(xi,g)xi,g,otherwise(5) \mathbf{x}_{i,g+1}=\begin{cases} \mathbf{u}_{i,g}, & \text{if } f(\mathbf{u}_{i,g}) \leq f(\mathbf{x}_{i,g})\\ \mathbf{x}_{i,g}, & \text{otherwise} \end{cases} \tag{5} xi,g+1={ui,g,xi,g,if f(ui,g)f(xi,g)otherwise(5)

并且它将作为下一代的父向量。

如果试验向量 ui,g\mathbf{u}_{i,g}ui,g 优于父向量 xi,g\mathbf{x}_{i,g}xi,g,即改进或进化进展 Δi,g=f(xi,g)−f(ui,g)\Delta_{i,g}=f(\mathbf{x}_{i,g})-f(\mathbf{u}_{i,g})Δi,g=f(xi,g)f(ui,g) 为正,则 (5) 式中的操作称为成功更新。相应地,用于生成 ui,g\mathbf{u}_{i,g}ui,g 的控制参数 FiF_{i}FiCRiCR_{i}CRi 分别称为成功变异因子成功交叉概率

上述的一对一选择程序在不同的 DE 算法中通常是固定的,而交叉操作除了 (4) 式中的二项式操作外,也可能有其他变体。因此,一个经典的差分进化算法在历史上被命名为,例如,DE/rand/1/bin,意味着它采用 DE/rand/1 变异策略和二项式交叉操作。


III. 自适应 DE 算法(ADAPTIVE DE ALGORITHMS)

本节简要回顾了一些近期的自适应和自适应性 DE 算法,这些算法随着进化搜索的进行而动态更新控制参数。这为后续比较这些参数自适应机制与 JADE 中提出的机制提供了便利。

A. DESAP

文献[10]提出了一种名为 DESAP 的算法,它不仅动态自适应变异和交叉参数 η\etaηδ\deltaδ(这通常是其他自适应或自适应性算法中的情况),还自适应种群大小 π\piπ。DESAP 的基础策略与经典的 DE/rand/1/bin 略有不同,而与文献[9]中提出的方案更为相似——虽然 δ\deltaδπ\piπ 的含义分别与 CRCRCRNPNPNP 相同,但 η\etaη 指的是在交叉后应用附加的正态分布变异的概率。实际上,普通的变异因子 FFF 在 DESAP 中是保持固定的。

在 DESAP 中,每个个体 iii 都关联其自身的控制参数 δi\delta_{i}δiηi\eta_{i}ηiπi\pi_{i}πi。这些参数会经历交叉和变异,方式与相应的向量 xi\mathbf{x}_{i}xi 类似。如果 f(ui)<f(xi)f(\mathbf{u}_{i}) < f(\mathbf{x}_{i})f(ui)<f(xi),这些控制参数的新值将与 ui\mathbf{u}_{i}ui 一起存活下来。尽管其原理简单,但 DESAP 的性能并不令人满意。在五个 De Jong 测试问题中,它仅在一个问题上优于传统 DE,而其他结果则非常相似。事实上,正如作者所述,DESAP 主要用于演示通过自适应更新种群大小以及其他控制参数来进一步减少控制参数的可能性。

B. FADE

由 Liu 和 Lampinen [11] 提出的模糊自适应差分进化(FADE),是一种新的 DE 变体,它使用模糊逻辑控制器来自适应变异和交叉操作的控制参数 FiF_{i}FiCRiCR_{i}CRi。类似的方法也被独立地提出用于多目标 DE 优化[12]。与许多其他自适应 DE 算法(DESAP [10] 除外)类似,FADE 假设种群大小是预先调优的,并在整个进化过程中保持固定。这种模糊逻辑控制方法在一组 10 个标准基准函数上进行了测试,并在问题维度较高时显示出比经典 DE 更好的结果。

C. SaDE

Qin 和 Suganthan 提出了 SaDE [13],它同时实施两种变异策略“DE/rand/1”和“DE/current-to-best/1”。它根据这两种策略在过去 50 代中的成功比率来调整使用任一种策略生成后代的概率。据信,这种自适应过程可以针对所考虑的问题,在不同学习阶段逐渐演化出最合适的变异策略。这与文献[27, 28]中提出的方案类似,该方案同时采用多种竞争启发式方法(包括几种 DE 变体、单纯形法和进化策略),并动态调整使用任何这些启发式方法生成后代的概率。

在 SaDE 中,变异因子在每个世代独立生成,遵循均值为 0.5、标准差为 0.3 的正态分布,并被截断到区间 (0, 2]。文中指出,这种方案可以在整个进化过程中保持局部(具有小 FiF_{i}Fi 值)和全局(具有大 FiF_{i}Fi 值)搜索能力,以生成潜在良好的变异向量。

交叉概率根据一个独立的正态分布随机生成,该分布均值为 CRmCR_mCRm,标准差为 0.1。与 FiF_iFi 相反,CRiCR_iCRi 值在下次重新生成之前会保持五代不变。均值 CRmCR_mCRm 初始化为 0.5。为了使 CRCRCR 适应到合适的值,作者每 25 代基于自上次更新 CRmCR_mCRm 以来记录的成功 CRCRCR 值来更新一次 CRmCR_mCRm

为了加速 SaDE 的收敛,作者进一步在 200 代后对一些优良个体应用了局部搜索程序(拟牛顿法)。

虽然 SaDE 最初是为解决无约束优化问题而提出的,但后来通过在求解一组约束优化问题时实现五种候选变异策略而得到了扩展[14]。

D. jDE

Brest 等人[15]提出了一种新的自适应 DE,jDE,它基于经典的 DE/rand/1/bin。与其他方案类似,jDE 在优化过程中固定种群大小,同时自适应与每个个体关联的控制参数 FiF_iFiCRiCR_iCRi。初始化过程将每个个体的 FiF_iFi 设为 0.5,CRiCR_iCRi 设为 0.9。与 SaDE 随机化 FiF_iFiCRiCR_iCRi 并通过记录的成功值更新 CRmCR_mCRm 的方法不同,jDE 在每个世代以概率 r1=r2=0.1r_1 = r_2 = 0.1r1=r2=0.1FiF_iFiCRiCR_iCRi 重新生成新值,分别根据 [0.1, 1] 和 [0, 1] 上的均匀分布。据信,更好的参数值倾向于产生更有可能存活的个体,因此这些值应该传播给下一代。

实验结果表明,jDE 的性能显著优于经典的 DE/rand/1/bin、FEP 和 CEP [29]、自适应的 LEP 和 Best Levy [30] 以及 FADE 算法 [11]。jDE 通过自适应两种变异策略得到了进一步扩展,新算法命名为 jDE-2 [17],在一组 25 个基准函数上显示出与原始 SaDE(在 200 代后实施了局部拟牛顿搜索程序)具有竞争力的结果。

IV. JADE

在本节中,我们提出一种新的差分进化算法JADE。该算法采用"DE/current-to-pbest"变异策略(带有可选存档机制),并通过自适应方式调整控制参数F和CR。JADE沿用了第II节公式(4)和(5)所描述的交叉与选择操作。

A. DE/Current-to-pbest

DE/rand/1是差分进化算法中最早提出的变异策略[1][2],也被认为是文献中最为成功且应用最广泛的方案[31]。然而,文献[6]指出DE/best/2可能具有优于DE/rand/1的特性,文献[32]则表明DE/best/1在多数技术问题的研究中表现更优。此外,文献[33]的作者认为融入最优解信息是有益的,并在其算法中采用了DE/current-to-best/1策略。相较于DE/rand/k类策略,诸如DE/current-to-best/k和DE/best/k这类贪婪策略通过融入最优解信息,在进化搜索中获得了更快的收敛速度。然而,这种最优解信息也可能导致种群多样性下降,进而引发早熟收敛等问题。

针对现有贪婪策略(式2和式3)收敛速度快但可靠性不足的问题,本文提出一种新型变异策略——“DE/current-to-pbest with optional archive”,作为本文提出的自适应DE算法JADE的核心基础。
在这里插入图片描述
图1. JADE 采用的 DE/current-to-pbest/1 变异策略示意图。 虚线表示优化问题的等高线。vᵢ 是使用相关联的变异因子 Fᵢ 为个体 xᵢ 生成的变异向量。

如图1所示,在DE/current-to-pbest/1策略(无存档机制)中,变异向量按如下方式生成:
在这里插入图片描述
其中,xbest,gpx_\mathrm{best,g}^pxbest,gp 是从当前种群中前100p%100p\%100p%p∈(0,1]p\in(0,1]p(0,1])的精英个体中随机选取的,FiF_iFi是与个体 xix_ixi 相关联的变异因子,其值通过后文式(10)提出的自适应机制在每一代重新生成。DE/current-to-pbest本质上是DE/current-to-best策略的推广形式——通过随机选取前100p%100p\%100p%精英个体中的任一解,取代DE/current-to-best中单一全局最优解的角色。

新近探索的劣势解相较于当前种群,能够为进化方向提供额外的有效信息。设A为存档的劣势解集合,P为当前种群。在带存档机制的DE/current-to-pbest/1策略中,变异向量生成方式如下:
在这里插入图片描述
其中,xi,gx_{i,g}xi,g𝑥𝑟1,𝑔𝑥_{𝑟1,𝑔}xr1,g𝑥pbest,𝑔𝑥_{pbest,𝑔}xpbest,g的选取方式与式(6)中相同(均从当前种群P中选择),而 x~r2,g\tilde{\mathbf{x}}_{r2,g}x~r2,g 则是从当前种群P与存档A的并集𝑃∪𝐴中随机选取的。

存档操作的设计力求简洁,以避免产生显著的计算开销。存档初始化为空集,在每代进化过程中,根据式(5)的选择操作被淘汰的父代解会被添加至存档。若存档规模超过预设阈值(例如NP),则通过随机移除部分解以维持存档容量为NP。显然,当设置存档规模为零时,式(7)即退化为式(6)的特例。

该档案库不仅提供了关于进化方向的信息,还能够有效提升种群的多样性。此外,如后续参数自适应调整公式(10)-(12)所示,我们通过鼓励使用较大的F值来增强种群多样性。因此,尽管算法会偏向潜在最优解(可能为局部极小值)的进化方向,但所提出的DE/current-to-pbest/1算法仍能有效避免陷入局部极小。简而言之,DE/rand/1算法在无特定方向偏好的情况下,仅在相对较小的区域内进行搜索;而结合档案库机制的DE/current-to-pbest/1算法则在更大的搜索区域内开展探索,且其搜索方向会偏向于有潜力的进化方向。

JADE的伪代码如表I所示。接下来介绍F和C R的参数自适应。
在这里插入图片描述

B. Parameter Adaptation

在第g代时,每个个体xi的交叉概率CRi均独立地根据均值为 μCRμ_{CR}μCR 、标准差为0.1的正态分布生成。
在这里插入图片描述
然后截断为[0,1]。将SCRS_{CR}SCR表示为第g代所有成功交叉概率CRiCR_iCRi的集合。均值μCRμ_{CR}μCR被初始化为0.5,并在每代结束时按如下方式更新:
在这里插入图片描述
其中c是0和1之间的正常数,meanA(⋅)meanA(·)meanA()是通常的算术平均值。

类似地,在第g代时,每个个体xix_ixi的变异因子FiF_iFi均独立地根据位置参数μFμ_FμF、尺度参数0.1的柯西分布生成。
在这里插入图片描述
Fi≥1F_i≥1Fi1时将其截断为1,若Fi≤0F_i≤0Fi0则重新生成。令SFS_FSF表示第g代中所有成功变异因子的集合。柯西分布的位置参数μFμ_FμF被初始化为0.5,并在每代结束时按如下方式更新:
在这里插入图片描述
其中,meanL(·) 表示莱默平均(Lehmer mean)。
在这里插入图片描述

C. Explanations of the Parameter Adaptation

μCRμ_{CR}μCR的自适应调整遵循以下核心原则:较优的控制参数值能够产生具有更高生存概率的个体,因此这些参数值应当被传递到后续进化代际中。其核心机制在于记录近期成功的交叉概率值,并利用这些历史信息指导新一代CRiCR_iCRi的生成。公式(8)中设定较小标准差的原因在于:若标准差过大,参数调整效率将显著降低——例如在标准差趋近无穷大的极端情况下,截断正态分布将与μCRμ_{CR}μCR的实际取值完全无关。

CRCRCR的调整相比,μFμ_FμF的适配包含两项独特操作:首先,FiF_iFi的生成采用截断柯西分布。与正态分布相比,柯西分布能更有效地实现变异因子的多样化,从而避免在贪婪变异策略(如DE/best、DE/current-to-best和DE/current-to-pbest)中经常出现的早熟收敛现象——这种现象通常发生在变异因子高度集中于某个特定值时。

其次,μFμ_FμF 的适配过程通过采用Lehmer均值(如公式12所示)而非μCRμ_{CR}μCR 调整中使用的算术均值,赋予较大成功变异因子更高权重这种Lehmer均值有助于传播更大的变异因子,从而提升进化进度率。相反,SFS_FSF 的算术均值往往会小于变异因子的最优值,导致μFμ_FμF偏小并引发后期早熟收敛。μFμ_FμF 偏小的趋势主要源于进化搜索中成功概率与进度率之间的差异:采用较小FiF_iFi 的DE/current-to-pbest策略本质上类似于(1+1)进化策略(ES)[23],两者都在基向量的小邻域内生成后代。对于(1+1)ES而言,已有研究表明(在corridor and sphere 函数上严格证明[23][34]),变异方差越小通常成功概率越高。但趋近于0的变异方差显然会导致进化停滞。因此,通过赋予较大成功变异因子更高权重来实现更快进化进度,是一种简洁而高效的解决方案。

关于公式(9)和(11)中的常数c:当c=0时,系统不进行参数自适应调整;反之,成功的CRiCR_iCRiFiF_iFi 的"生命周期"约为1/c代——具体而言,经过1/c代后,当c趋近于零时,原有的μCRμ_{CR}μCRμFμ_FμF值会按(1−c)1/c(1−c)^{1/c}(1c)1/c→1/e≈37%的比例衰减。

D. Discussion of Parameter Settings

经典差分进化算法(DE)包含两个需要用户调整的控制参数F和CR。这些参数具有问题依赖性,因此通常需要繁琐的试错过程才能为特定问题找到合适的取值。与之相反,JADE算法新引入的两个参数c和p预计对不同问题具有不敏感性——这源于它们在JADE中的特定作用:c控制参数自适应速率,p决定变异策略的贪婪程度。如第五节所示,JADE通常在1/c∈[5,20](即μCRμ_{CR}μCRμFμ_FμF值的生命周期为5-20代)且p∈[5%,20%]时表现最优,这意味着变异过程中我们仅考虑前5%-20%的高质量解。
好的,这是论文《JADE: Adaptive Differential Evolution with Optional External Archive》中“V. SIMULATION”部分的翻译。


V. 仿真实验(SIMULATION)

在本节中,JADE 被应用于最小化一组 13 个可扩展的基准函数(维度为 D=30D=30D=30 或 100)[29, 35] 和一组低维度(D=2∼6D=2\sim 6D=26)的 Dixon-Szego 函数 [36],如表 II 和表 III 所示。

JADE 与两种最近的自适应 DE 算法(jDE 和 SaDE)、经典的 DE/rand/1/bin 以及经典粒子群优化(PSO)算法 [37] 进行了比较。同时,还与 rand-JADE 和 nona-JADE 这两个 JADE 的变体进行了比较,以识别其不同组件带来的益处。为公平比较,我们在所有仿真中将 JADE 的参数固定设置为 p=0.05p=0.05p=0.05c=0.1c=0.1c=0.1。我们遵循 jDE [15] 和 SaDE [13] 原论文中的参数设置,但在 SaDE 中禁用了拟牛顿局部搜索。对于 PSO,参数值选自 [37],因为根据我们研究的其他文献,这些值通常效果更好。DE/rand/1/bin 的参数设置为 F=0.5F=0.5F=0.5CR=0.9CR=0.9CR=0.9,这与 [1, 8, 15, 38] 中使用或推荐的值一致。

表 II 总结了这 13 个可扩展的基准函数。所有这些函数的最优值均为 f∗=0f^{*}=0f=0,它们的一些不同特性简要总结如下:f1f_{1}f1-f4f_{4}f4 是连续单峰函数。f5f_{5}f5 是 Rosenbrock 函数,在 D=2D=2D=2 和 3 时是单峰的,但在高维情况下可能具有多个极小值 [39]。f6f_{6}f6 是一个不连续的阶梯函数,f7f_{7}f7 是一个带噪声的四次函数。f8f_{8}f8-f13f_{13}f13 是多峰函数,其局部极小值的数量随问题维度呈指数级增长 [29]。此外,f8f_{8}f8 是本文研究的唯一一个带边界约束的函数。多峰 Dixon-Szego 函数的特性在表 III 中简要描述,感兴趣的读者可参阅 [36] 了解详情。

在所有仿真中,当 D≤10D\leq 10D10=30=30=30=100=100=100 时,我们分别将种群大小 NPNPNP 设置为 30、100 和 400。本节报告的所有结果均基于 50 次独立运行得到。

为清晰起见,最佳次佳算法的结果分别用粗体斜体标出;如果并非所有或大多数算法产生相同结果。

表 II 维度为 D 的测试函数。每个函数的全局最小值均为 0。它们分别是 Sphere、Schwefel 2.22、Schwefel 1.2、Schwefel 2.21、Rosenbrock、Step、Noisy Quartic、Schwefel 2.26、Rastrigin、Ackley、Griewank 以及两个 Penalized 函数 [35]。
在这里插入图片描述
表 III 维度为 D 的测试函数。每个函数的全局最小值均为 0。它们分别是 Sphere、Schwefel 2.22、Schwefel 1.2、Schwefel 2.21、Rosenbrock、Step、Noisy Quartic、Schwefel 2.26、Rastrigin、Ackley、Griewank 以及两个 Penalized 函数 [36]。
在这里插入图片描述

A. JADE 与其他进化算法的比较(Comparison of JADE With Other Evolutionary Algorithms)

各算法对 f1f_{1}f1-f13f_{13}f13 所得结果的均值和标准差总结在表 IV 和表 V 中。这些统计数据仅针对 f1f_{1}f1-f4f_{4}f4f7f_{7}f7 在优化结束时计算。对于其他函数,也报告了中期结果,因为不同算法获得的最终结果可能由于 MATLAB 的精度问题而相同为零或达到某个残余误差(参见图 3 中 f10f_{10}f10f13f_{13}f13 的误差平台)。在这些情况下,仅比较中期结果,并可能以粗体或斜体标出。

表 IV 30 维问题 f1–f13f_1–f_{13}f1f13 的实验结果(基于 50 次独立运行的平均值)
在这里插入图片描述
表 V 100 维问题 f1–f13f_1–f_{13}f1f13 的实验结果(基于 50 次独立运行的平均值)
在这里插入图片描述
在表 VI 和表 VII 中,我们总结了各算法的成功率(SR)成功运行下的平均函数评估次数(FESS)。如果找到的最佳解满足足够的精度(对于噪声函数 f7f_{7}f710−210^{-2}102,对于所有其他函数为 10−810^{-8}108),则认为实验成功。FESS 和 SR 分别用于比较不同算法在成功运行中的收敛速度和可靠性。

表 VI 30 维问题 f1–f13f_1–f_{13}f1f13 的实验结果(基于 50 次独立运行的平均值)。
注意:在 JADE 及其他竞争算法中(不包括 JADE 的两个变体),最佳和次佳结果分别用粗体和斜体标出。
在这里插入图片描述
表 VII 100 维问题 f₁–f₁₃ 的实验结果(基于 50 次独立运行的平均值)
注:在 JADE 及其他竞争算法中(JADE 的两个变体除外),最佳和次佳结果分别以粗体和斜体标出。
在这里插入图片描述
为便于说明,我们在图 2 和图 3 中分别绘制了一些 30 维和 100 维问题的收敛曲线图。注意,在这些图中我们绘制的是中值曲线(而非表格中报告的均值),因为当一个算法可能只是偶尔导致假收敛时,这些曲线与箱线图一起能提供更多信息。箱线图在特定代数处为最佳和次佳算法绘制,有助于说明 50 次独立运行结果的分布情况。
在这里插入图片描述
图2. 测试函数 f1, f3, f5 和 f7 (维度 D=30) 的收敛曲线图(中值曲线与箱线图)。横轴为进化代数,纵轴为50次独立运行的函数值中位数。箱线图同时为JADE算法和其他算法中表现最佳者绘制。
在这里插入图片描述
图3. 测试函数 f9, f10, f11 和 f13 (维度 D=100) 的收敛曲线图(中值曲线与箱线图)。横轴为进化代数,纵轴为50次独立运行的函数值中位数。箱线图同时为JADE算法和其他算法中表现最佳者绘制。

从图 2、图 3 以及表 IV 和表 V 呈现的结果中,可以得出关于不同算法收敛速度和可靠性的重要观察结果。首先,这些结果表明,对于问题集 f1f_{1}f1-f13f_{13}f13,带或不带归档的 JADE 的总体收敛速度分别是最佳和次佳的:不带归档的 JADE 在相对低维问题(D=30D=30D=30)上效果最好,而带归档的 JADE 在高维问题(D=100D=100D=100)情况下实现了最佳性能。这可能是因为对于各种 100 维问题,种群大小 400 并不足够,而归档可以通过在变异操作中引入更多多样性来在一定程度上缓解这个问题。在两种情况下(带和不带归档),JADE 在优化 f1f_{1}f1-f7f_{7}f7f10f_{10}f10-f13f_{13}f13 时收敛最快,这些函数包括单峰或多峰、无噪声或有噪声、连续或不连续。在 f8f_{8}f8f9f_{9}f9 的情况下,jDE 表现最佳,而 JADE 取得了非常有竞争力的结果。SaDE 和 PSO 的收敛速度通常比 JADE 和 jDE 差,但比 DE/rand/1/bin 好。

其次,比较不同算法的可靠性也很重要,如表 VI 和表 VII 所示。PSO 和 DE/rand/1/bin 表现最差,因为前者饱受频繁的早熟收敛之苦,而后者在有限的函数评估次数内通常收敛非常缓慢。SaDE 优于 PSO 和 DE/rand/1/bin,但对于某些问题,尤其是在问题维度较高时,其可靠性明显不尽如人意。正如预期的那样,jDE 表现非常好,因为其基础的变异策略“DE/rand/1”比贪婪策略更鲁棒。有趣的是,对于所有 13 个可扩展函数,JADE 的可靠性与 jDE 相似。JADE 的高可靠性和快速收敛源于其自适应参数控制以及归档辅助的“DE/current-to-pbest/1”变异策略所提供的方向信息和多样性改善。

在低维 Dixon-Szego 函数的情况下,表 VIII 中的仿真结果表明没有明显更优的算法。这与文献 [15] 中的观察类似,即对于这组函数,自适应算法并未显示出相对于经典 DE/rand/1/bin 的明显性能提升。经过进一步研究,这一观察结果并不令人惊讶。首先,对 DE/rand/1/bin 使用从 {0.1, 0.3, 0.5, 0.7, 0.9}\{0.1,\,0.3,\,0.5,\,0.7,\,0.9\}{0.1,0.3,0.5,0.7,0.9} 中选择的 FFFCRCRCR 进行的各种仿真表明,参数设置 F=0.5F=0.5F=0.5CR=0.9CR=0.9CR=0.9 对于所有 Dixon-Szego 函数始终能产生最佳或接近最佳的性能。其次,仿真结果也表明,在优化这些低维问题所需的少量代数内,参数自适应无法有效发挥作用。
表 VIII
Dixon-Szegö 函数 f₁₄–f₂₀ 的实验结果(基于 50 次独立运行的平均值)
在这里插入图片描述
在表 IX 中,JADE 进一步与文献中报告结果的其他进化算法进行了比较:文献[30, Table III]中的自适应 LEP 和 Best Levy 算法,文献[40, Tables I-3]中的 NSDE 以及文献[15, Table III]中的 jDE。很明显,JADE 在大多数情况下表现最佳或次佳,并且实现了比其他竞争算法更好的整体性能。

表 IX 函数 f₁, f₃, f₅, f₈–f₁₃ (D=30) 及 Dixon-Szegö 测试函数 f₁₈–f₁₉ 的实验结果:其中带归档与不带归档的 JADE 结果基于 50 次独立运行计算;自适应 LEP 与 BestLevy 的结果取自 [30, Table III];NSDE 的结果取自 [40, Tables I–3];jDE 的结果取自 [15, Table III]
在这里插入图片描述

B. JADE 各组件的益处(The Benefit of JADE Components)

我们有兴趣识别不带归档的 JADE 的两个组件(“DE/current-to-pbest”变异策略和参数自适应)带来的益处。为此,我们考虑了 JADE 的两个变体,即 rand-JADE 和 nona-JADE。它们与不带归档的 JADE 在以下方面不同:rand-JADE 实施 DE/rand/1(而不是 DE/current-to-pbest)变异策略,而 nona-JADE 不采用自适应参数控制 [即,我们不根据 (9) 和 (11) 式更新 μCR\mu_{CR}μCRμF\mu_{F}μF,而是在整个优化过程中在 (8) 和 (10) 式中将 μCR=μF=0.5\mu_{CR}=\mu_{F}=0.5μCR=μF=0.5 设为固定值]。

表 VI 和表 VII 总结了 JADE、rand-JADE、nona-JADE 和 DE/rand/1/bin 的仿真结果。FESS 性能表明,nona-JADE(在大多数情况下)和 rand-JADE(在某些情况下)比 DE/rand/1/bin 收敛更快,这分别归因于它们相对贪婪的变异策略和自适应参数控制。然而,rand-JADE 和 nona-JADE 都对某些函数存在频繁的早熟收敛问题。

与它的两个变体相比,JADE 在收敛速度和可靠性方面都取得了显著更好的性能。这表明贪婪策略“DE/current-to-pbest”和参数自适应之间存在互利合作。实际上,贪婪变异策略可能从两个相反的方向影响种群多样性:它倾向于通过将个体移向少数几个最优解来降低多样性,但也可能通过加快在前进方向上的优化搜索来增加多样性。在经典的不带参数自适应的 DE 算法中,贪婪变异策略通常会导致早熟收敛,因为多样性降低的前者效应起着关键作用。然而,参数自适应方案能够将参数调整到合适的值,从而提高“DE/current-to-pbest”在有希望方向上的进展速度。因此,多样性增加的后者效应变得能够平衡前者效应,从而提高了算法的收敛性能。

C. JADE 中 μF\mu_FμFμCR\mu_{CR}μCR 的演化(Evolution of μF\mu_FμF and μCR\mu_{CR}μCR in JADE)

在经典 DE 中,两个控制参数 FFFCRCRCR 需要针对不同的优化函数通过试错法进行调整。在 JADE 中,它们由自适应参数 μF\mu_FμFμCR\mu_{CR}μCR 控制,这些参数随着算法进程而演化。μF\mu_FμFμCR\mu_{CR}μCR 的演化过程以均值曲线和误差棒的形式绘制在图 4 中。误差棒暗示了 μF\mu_FμFμCR\mu_{CR}μCR 在 50 次独立运行中清晰的演化趋势。例如,在球形函数 f1f_1f1 的优化中,μF\mu_FμFμCR\mu_{CR}μCR 变化很小。这是合理的,因为在不同进化阶段,适应度地形的形状是相同的。对于椭球函数 f3f_3f3μF\mu_FμFμCR\mu_{CR}μCR 在经历明显的初始变化后达到稳定状态(即初始分布在……的种群)。这与 f5f_5f5 的情况类似,当种群落入 Rosenbrock 函数的狭窄山谷后,μF\mu_FμFμCR\mu_{CR}μCR 达到稳定状态。对于 f9f_9f9,适应度地形随着算法进程呈现出不同的形状,因此 μF\mu_FμFμCR\mu_{CR}μCR 也相应地演化到不同的值。

这些观察结果与直觉一致;即,不存在固定的 FFFCRCRCR 参数设置能够适用于各种问题(即使是像球形函数 f1f_{1}f1 和椭球函数 f3f_{3}f3 这样可线性变换的问题),或者适用于同一问题在不同进化阶段的需求。

我们进一步研究了 μF\mu_{F}μFμCR\mu_{CR}μCR 初始值对参数自适应过程的影响。根据广泛的实验研究(由于篇幅限制未报告结果),μF\mu_{F}μF 的初始值对算法性能影响很小,而建议使用中等至较大的 μCR\mu_{CR}μCR 初始值,特别是对于不可分函数。总的来说,μCR=μF=0.5\mu_{CR}=\mu_{F}=0.5μCR=μF=0.5 的初始设置对所有测试函数都效果良好,因此被认为是 JADE 的标准设置

D. JADE 的参数值(Parameter Values of JADE)

JADE 引入了自身的参数 cccppp,它们分别决定了 μCR\mu_{CR}μCRμF\mu_{F}μF 的自适应速率以及变异策略的贪婪程度。据信,根据它们在 JADE 中的作用,这两个参数对问题不敏感;因此,与经典 DE 中依赖于问题的参数(FFFCRCRCR)选择相比,这是一个优势。然而,寻找适用于不同问题的这些参数的合适范围仍然是有意义的。如图 5 所示,绘制了 JADE 在不同参数组合下的最佳值中位数:1/c∈{1,2,5,10,20,50,100}1/c\in\{1,2,5,10,20,50,100\}1/c{1,2,5,10,20,50,100}p⋅NP∈{1,3,5,10,20,30,50}p\cdot NP\in\{1,3,5,10,20,30,50\}pNP{1,3,5,10,20,30,50},其中 NP=100NP=100NP=100。正如预期,较小的 1/c1/c1/c 值(例如,1/c=11/c=11/c=1)或 ppp 值(例如,p⋅NP=1p\cdot NP=1pNP=1)在某些情况下可能导致不太令人满意的结果。前者由于缺乏足够的信息来平滑更新 CRCRCRFFF 而导致假收敛,而后者则过于贪婪而无法维持种群的多样性。然而,很明显,JADE 在很大的 cccppp 参数范围内(与图 2 和图 3 中其他算法的结果相比)优于或与其他算法具有竞争力。具体来说,它表明 JADE(带或不带归档)在 1/c∈[5,20]1/c\in[5,20]1/c[5,20]p∈[5%, 20%]p\in[5\%,\,20\%]p[5%,20%] 的取值范围内工作最佳。


VI. 结论(CONCLUSION)

参数自适应通过在进化搜索过程中自动将控制参数更新至合适的值,有利于提高进化算法的优化性能。很自然地,可以将参数自适应与贪婪变异策略(如“DE/current-to-pbest”)相结合,以提高收敛速度,同时将算法的可靠性维持在较高水平。这促使我们提出了 JADE,一种参数自适应的“current-to-pbest”差分进化算法。在“current-to-pbest”中,我们利用了多个最优解的信息,以平衡变异的贪婪性与种群的多样性。JADE 的参数自适应是通过基于变异因子和交叉概率的历史成功记录来演化这些参数而实现的。我们还引入了外部归档来存储近期探索过的较差解,并利用它们与当前种群的差异作为指向最优解的有希望的方向。

JADE 在一组取自文献的经典基准函数上进行了测试。与其他经典和自适应 DE 算法以及经典 PSO 算法相比,JADE 在收敛速度和可靠性方面显示出更好或至少具有竞争力的优化性能。根据文献中报道的结果,JADE 也显示出优于其他进化算法的性能。

带归档的 JADE 在优化高维问题时显示出了良好的结果。与其他技术(如协同协同进化 [41, 42])相比,JADE 在提高高维问题收敛性能方面有着相当不同的侧重点。期望 JADE 可以作为一个合适的基础方案,将协同协同进化融入其中以获得更好的性能 [43, 44]。此外,将当前工作扩展到多目标优化(MOO)也很有意义,因为参数自适应多目标优化已被评论为“一个非常有趣但很少有研究者在专业文献中涉及的课题 [45]”。尽管我们的初步研究已显示出有希望的结果 [46],但在将参数自适应方案纳入多目标进化优化方面,仍存在许多悬而未决的问题。

Matlab Code

JADE.m

JADE算法。

%**************************************************************************************************
%Reference:  J. Zhang and A. C. Sanderson, "JADE: adaptive differential evolution
%                     with optional external archive," IEEE Trans. Evolut. Comput., vol. 13,
%                     no. 5, pp. 945-958, 2009.
%
% Note: We obtained the MATLAB source code from the authors, and did some
%           minor revisions in order to solve the 25 benchmark test functions,
%           however, the main body was not changed.
%**************************************************************************************************

clc;
clear all;
tic;

format long;
format compact;

'JADE'

% Choose the problems to be tested. Please note that for test functions F7
% and F25, the global optima are out of the initialization range. For these
% two test functions, we do not need to judge whether the variable violates
% the boundaries during the evolution after the initialization.
problemSet = [1 : 6 8 : 24];
for problemIndex = 1 : 23

    problem = problemSet(problemIndex)

    % Define the dimension of the problem
    n = 30;

    popsize = 100;

    switch problem

        case 1

            % lu: define the upper and lower bounds of the variables
            lu = [-100 * ones(1, n); 100 * ones(1, n)];
            % Load the data for this test function
            load sphere_func_data
            A = []; M = []; a = []; alpha = []; b = [];

        case 2

            lu = [-100 * ones(1, n); 100 * ones(1, n)];
            load schwefel_102_data
            A = []; M = []; a = []; alpha = []; b = [];

        case 3

            lu = [-100 * ones(1, n); 100 * ones(1, n)];
            load high_cond_elliptic_rot_data
            A = []; a = []; alpha = []; b = [];

            if n == 2, load elliptic_M_D2,
            elseif n == 10, load elliptic_M_D10,
            elseif n == 30, load elliptic_M_D30,
            elseif n == 50, load elliptic_M_D50,
            end

        case 4

            lu = [-100 * ones(1, n); 100 * ones(1, n)];
            load schwefel_102_data
            A = []; M = []; a = []; alpha = []; b = [];

        case 5

            lu = [-100 * ones(1, n); 100 * ones(1, n)];
            load schwefel_206_data
            M = []; a = []; alpha = []; b = [];

        case 6

            lu = [-100 * ones(1, n); 100 * ones(1, n)];
            load rosenbrock_func_data
            A = []; M = []; a = []; alpha = []; b = [];

        case 7

            lu = [0 * ones(1, n); 600 * ones(1, n)];
            load griewank_func_data
            A = []; a = []; alpha = []; b = [];

            if n == 2, load griewank_M_D2,
            elseif n == 10, load griewank_M_D10,
            elseif n == 30, load griewank_M_D30,
            elseif n == 50, load griewank_M_D50,
            end

        case 8

            lu = [-32 * ones(1, n); 32 * ones(1, n)];
            load ackley_func_data
            A = []; a = []; alpha = []; b = [];

            if n == 2, load ackley_M_D2,
            elseif n == 10, load ackley_M_D10,
            elseif n == 30, load ackley_M_D30,
            elseif n == 50, load ackley_M_D50,
            end

        case 9

            lu = [-5 * ones(1, n); 5 * ones(1, n)];
            load rastrigin_func_data
            A = []; M = []; a = []; alpha = []; b = [];

        case 10

            lu = [-5 * ones(1, n); 5 * ones(1, n)];
            load rastrigin_func_data
            A = []; a = []; alpha = []; b = [];
            if n == 2, load rastrigin_M_D2,
            elseif n == 10, load rastrigin_M_D10,
            elseif n == 30, load rastrigin_M_D30,
            elseif n == 50, load rastrigin_M_D50,
            end

        case 11

            lu = [-0.5 * ones(1, n); 0.5 * ones(1, n)];
            load weierstrass_data
            A = []; a = []; alpha = []; b = [];
            if n == 2, load weierstrass_M_D2, ,
            elseif n == 10, load weierstrass_M_D10,
            elseif n == 30, load weierstrass_M_D30,
            elseif n == 50, load weierstrass_M_D50,
            end

        case 12

            lu = [-pi * ones(1, n); pi * ones(1, n)];
            load schwefel_213_data
            A = []; M = []; o = [];

        case 13

            lu = [-3 * ones(1, n); 1 * ones(1, n)];
            load EF8F2_func_data
            A = []; M = []; a = []; alpha = []; b = [];

        case 14

            lu = [-100 * ones(1, n); 100 * ones(1, n)];
            load E_ScafferF6_func_data
            if n == 2, load E_ScafferF6_M_D2, ,
            elseif n == 10, load E_ScafferF6_M_D10,
            elseif n == 30, load E_ScafferF6_M_D30,
            elseif n == 50, load E_ScafferF6_M_D50,
            end
            A = []; a = []; alpha = []; b = [];

        case 15

            lu = [-5 * ones(1, n); 5 * ones(1, n)];
            load hybrid_func1_data
            A = []; M = []; a = []; alpha = []; b = [];

        case 16

            lu = [-5 * ones(1, n); 5 * ones(1, n)];
            load hybrid_func1_data
            if n == 2, load hybrid_func1_M_D2,
            elseif n == 10, load hybrid_func1_M_D10,
            elseif n == 30, load hybrid_func1_M_D30,
            elseif n == 50, load hybrid_func1_M_D50,
            end
            A = []; a = []; alpha = []; b = [];

        case 17

            lu = [-5 * ones(1, n); 5 * ones(1, n)];
            load hybrid_func1_data
            if n == 2, load hybrid_func1_M_D2,
            elseif n == 10, load hybrid_func1_M_D10,
            elseif n == 30, load hybrid_func1_M_D30,
            elseif n == 50, load hybrid_func1_M_D50,
            end
            A = []; a = []; alpha = []; b = [];

        case 18

            lu = [-5 * ones(1, n); 5 * ones(1, n)];
            load hybrid_func2_data
            if n == 2, load hybrid_func2_M_D2,
            elseif n == 10, load hybrid_func2_M_D10,
            elseif n == 30, load hybrid_func2_M_D30,
            elseif n == 50, load hybrid_func2_M_D50,
            end
            A = []; a = []; alpha = []; b = [];

        case 19

            lu = [-5 * ones(1, n); 5 * ones(1, n)];
            load hybrid_func2_data
            if n == 2, load hybrid_func2_M_D2,
            elseif n == 10, load hybrid_func2_M_D10,
            elseif n == 30, load hybrid_func2_M_D30,
            elseif n == 50, load hybrid_func2_M_D50,
            end
            A = []; a = []; alpha = []; b = [];

        case 20

            lu = [-5 * ones(1, n); 5 * ones(1, n)];
            load hybrid_func2_data
            if n == 2, load hybrid_func2_M_D2,
            elseif n == 10, load hybrid_func2_M_D10,
            elseif n == 30, load hybrid_func2_M_D30,
            elseif n == 50, load hybrid_func2_M_D50,
            end
            A = []; a = []; alpha = []; b = [];

        case 21

            lu = [-5 * ones(1, n); 5 * ones(1, n)];
            load hybrid_func3_data
            if n == 2, load hybrid_func3_M_D2,
            elseif n == 10, load hybrid_func3_M_D10,
            elseif n == 30, load hybrid_func3_M_D30,
            elseif n == 50, load hybrid_func3_M_D50,
            end
            A = []; a = []; alpha = []; b = [];

        case 22

            lu = [-5 * ones(1, n); 5 * ones(1, n)];
            load hybrid_func3_data
            if n == 2, load hybrid_func3_HM_D2,
            elseif n == 10, load hybrid_func3_HM_D10,
            elseif n == 30, load hybrid_func3_HM_D30,
            elseif n == 50, load hybrid_func3_HM_D50,
            end
            A = []; a = []; alpha = []; b = [];

        case 23

            lu = [-5 * ones(1, n); 5 * ones(1, n)];
            load hybrid_func3_data
            if n == 2, load hybrid_func3_M_D2,
            elseif n == 10, load hybrid_func3_M_D10,
            elseif n == 30, load hybrid_func3_M_D30,
            elseif n == 50, load hybrid_func3_M_D50,
            end
            A = []; a = []; alpha = []; b = [];

        case 24

            lu = [-5 * ones(1, n); 5 * ones(1, n)];
            load hybrid_func4_data
            if n == 2, load hybrid_func4_M_D2,
            elseif n == 10, load hybrid_func4_M_D10,
            elseif n == 30, load hybrid_func4_M_D30,
            elseif n == 50, load hybrid_func4_M_D50,
            end
            A = []; a = []; alpha = []; b = [];

        case 25

            lu = [2 * ones(1, n); 5 * ones(1, n)];
            load hybrid_func4_data
            if n == 2, load hybrid_func4_M_D2,
            elseif n == 10, load hybrid_func4_M_D10,
            elseif n == 30, load hybrid_func4_M_D30,
            elseif n == 50, load hybrid_func4_M_D50,
            end
            A = []; a = []; alpha = []; b = [];

    end

    % Record the best results
    outcome = [];

    %Main body which was developed by the authors

    time = 1; % 当前运行次数

    % The total number of runs
    totalTime = 1;

    while time <= totalTime

        rand('seed', sum(100 * clock)); % 设置随机种子

        % Initialize the main population
        popold = repmat(lu(1, :), popsize, 1) + rand(popsize, n) .* (repmat(lu(2, :) - lu(1, :), popsize, 1));

        valParents = benchmark_func(popold, problem, o, A, M, a, alpha, b); % 计算适应度值 NP*1

        c = 1/10;
        p = 0.05;

        CRm = 0.5;  % 初始化CR均值
        Fm = 0.5; % 初始化F均值

        Afactor = 1;  % 用于控制外部存档和种群大小的比例

        archive.NP = Afactor * popsize; % the maximum size of the archive 外部存档的大小
        archive.pop = zeros(0, n); % the solutions stored in te archive 劣解
        archive.funvalues = zeros(0, 1); % the function value of the archived solutions 劣解的适应度值

        %% the values and indices of the best solutions 升序排序,并得到对应的索引,最好的值在最前面
        [valBest, indBest] = sort(valParents, 'ascend');

        FES = 0;  % 函数评估次数
        while FES < n * 10000 %& min(fit)>error_value(problem)

            pop = popold; % the old population becomes the current population

            if FES > 1 && ~isempty(goodCR) && sum(goodF) > 0 % If goodF and goodCR are empty, pause the update
                CRm = (1 - c) * CRm + c * mean(goodCR);
                Fm = (1 - c) * Fm + c * sum(goodF .^ 2) / sum(goodF); % Lehmer mean
            end

            % Generate CR according to a normal distribution with mean CRm, and std 0.1
            % Generate F according to a cauchy distribution with location parameter Fm, and scale parameter 0.1
            [F, CR] = randFCR(popsize, CRm, 0.1, Fm, 0.1);  % 每个个体对应一个F和CR

            r0 = [1 : popsize];
            popAll = [pop; archive.pop];  % PUA
            [r1, r2] = gnR1R2(popsize, size(popAll, 1), r0);  % 生成参数r1和r2

            % Find the p-best solutions
            pNP = max(round(p * popsize), 2); % choose at least two best solutions
            randindex = ceil(rand(1, popsize) * pNP); % select from [1, 2, 3, ..., pNP]
            randindex = max(1, randindex); % to avoid the problem that rand = 0 and thus ceil(rand) = 0
            pbest = pop(indBest(randindex), :); % randomly choose one of the top 100p% solutions

            % == == == == == == == == == == == == == == == Mutation == == == == == == == == == == == == ==
            vi = pop + F(:, ones(1, n)) .* (pbest - pop + pop(r1, :) - popAll(r2, :));

            vi = boundConstraint(vi, pop, lu);  % 越界处理

            % == == == == = Crossover == == == == =
            mask = rand(popsize, n) > CR(:, ones(1, n)); % mask is used to indicate which elements of ui comes from the parent
            % choose one position where the element of ui doesn't come from the parent
            rows = (1 : popsize)';cols = floor(rand(popsize, 1) * n)+1;
            jrand = sub2ind([popsize n], rows, cols); % 使用 sub2ind 将子脚标转换为线性索引
            mask(jrand) = false;  % 确保每个个体至少有一个维度来自突变向量
            ui = vi; ui(mask) = pop(mask);

            valOffspring = benchmark_func(ui, problem, o, A, M, a, alpha, b); % 计算适应度值

            FES = FES + popsize;  % 函数评估次数+popsize

            % == == == == == == == == == == == == == == == Selection == == == == == == == == == == == == ==
            % I == 1: the parent is better; I == 2: the offspring is better
            [valParents, I] = min([valParents, valOffspring], [], 2);
            popold = pop;

            archive = updateArchive(archive, popold(I == 2, :), valParents(I == 2));  % 更新外部存档

            popold(I == 2, :) = ui(I == 2, :);

            goodCR = CR(I == 2);
            goodF = F(I == 2);

            [valBest indBest] = sort(valParents, 'ascend');

        end

        outcome = [outcome min(valParents)];  % 记录本次运行的最优值

        time = time + 1;

    end

    sort(outcome)
    mean(outcome)
    std(outcome)

end
toc;

randFCR.m

为种群中的每个个体生成F和CR。

function [F,CR] = randFCR(NP, CRm, CRsigma, Fm,  Fsigma)

% this function generate CR according to a normal distribution with mean "CRm" and sigma "CRsigma"
%           If CR > 1, set CR = 1. If CR < 0, set CR = 0.
% this function generate F  according to a cauchy distribution with location parameter "Fm" and scale parameter "Fsigma"
%           If F > 1, set F = 1. If F <= 0, regenrate F.
%
% Version: 1.1   Date: 11/20/2007
% Written by Jingqiao Zhang (jingqiao@gmail.com)

%% generate CR
CR = CRm + CRsigma * randn(NP, 1);
CR = min(1, max(0, CR));                % truncated to [0 1]

%% generate F
F = randCauchy(NP, 1, Fm, Fsigma);
F = min(1, F);                          % truncation

% we don't want F = 0. So, if F<=0, we regenerate F (instead of trucating it to 0)
pos = find(F <= 0);
while ~ isempty(pos)
    F(pos) = randCauchy(length(pos), 1, Fm, Fsigma);
    F = min(1, F);                      % truncation
    pos = find(F <= 0);
end

% Cauchy distribution: cauchypdf = @(x, mu, delta) 1/pi*delta./((x-mu).^2+delta^2)
function result = randCauchy(m, n, mu, delta)

% http://en.wikipedia.org/wiki/Cauchy_distribution
result = mu + delta * tan(pi * (rand(m, n) - 0.5));

updateArchive.m

更新存档并输入解决方案
步骤1:将新解决方案添加至存档
步骤2:移除重复项
步骤3:如有需要,随机删除部分解决方案以维持存档容量

function archive = updateArchive(archive, pop, funvalue)
% Update the archive with input solutions
%   Step 1: Add new solution to the archive
%   Step 2: Remove duplicate elements
%   Step 3: If necessary, randomly remove some solutions to maintain the archive size
%
% Version: 1.1   Date: 2008/04/02
% Written by Jingqiao Zhang (jingqiao@gmail.com)

if archive.NP == 0, return; end

if size(pop, 1) ~= size(funvalue,1), error('check it'); end

% Method 2: Remove duplicate elements
popAll = [archive.pop; pop ];
funvalues = [archive.funvalues; funvalue ];
[dummy IX]= unique(popAll, 'rows');
if length(IX) < size(popAll, 1) % There exist some duplicate solutions
  popAll = popAll(IX, :);
  funvalues = funvalues(IX, :);
end

if size(popAll, 1) <= archive.NP   % add all new individuals
  archive.pop = popAll;
  archive.funvalues = funvalues;
else                % randomly remove some solutions
  rndpos = randperm(size(popAll, 1)); % equivelent to "randperm";
  rndpos = rndpos(1 : archive.NP);
  
  archive.pop = popAll  (rndpos, :);
  archive.funvalues = funvalues(rndpos, :);
end

gnR1R2.m

选择两个随机个体r1和r2。

function [r1, r2] = gnR1R2(NP1, NP2, r0)

% gnA1A2 generate two column vectors r1 and r2 of size NP1 & NP2, respectively
%    r1's elements are choosen from {1, 2, ..., NP1} & r1(i) ~= r0(i)
%    r2's elements are choosen from {1, 2, ..., NP2} & r2(i) ~= r1(i) & r2(i) ~= r0(i)
%
% Call:
%    [r1 r2 ...] = gnA1A2(NP1)   % r0 is set to be (1:NP1)'
%    [r1 r2 ...] = gnA1A2(NP1, r0) % r0 should be of length NP1
%
% Version: 2.1  Date: 2008/07/01
% Written by Jingqiao Zhang (jingqiao@gmail.com)

NP0 = length(r0);

r1 = floor(rand(1, NP0) * NP1) + 1;
for i = 1 : inf
    pos = (r1 == r0);
    if sum(pos) == 0
        break;
    else % regenerate r1 if it is equal to r0
        r1(pos) = floor(rand(1, sum(pos)) * NP1) + 1;
    end
    if i > 1000, % this has never happened so far
        error('Can not genrate r1 in 1000 iterations');
    end
end

r2 = floor(rand(1, NP0) * NP2) + 1;
for i = 1 : inf
    pos = ((r2 == r1) | (r2 == r0));
    if sum(pos)==0
        break;
    else % regenerate r2 if it is equal to r0 or r1
        r2(pos) = floor(rand(1, sum(pos)) * NP2) + 1;
    end
    if i > 1000, % this has never happened so far
        error('Can not genrate r2 in 1000 iterations');
    end
end

boundConstraint.m

越界处理。

function vi = boundConstraint (vi, pop, lu)

% if the boundary constraint is violated, set the value to be the middle
% of the previous value and the bound
%
% Version: 1.1   Date: 11/20/2007
% Written by Jingqiao Zhang, jingqiao@gmail.com

[NP, D] = size(pop);  % the population size and the problem's dimension

%% check the lower bound
xl = repmat(lu(1, :), NP, 1);
pos = vi < xl;  % NP*D找到所有违反约束下界的值的位置
vi(pos) = (pop(pos) + xl(pos)) / 2;

%% check the upper bound
xu = repmat(lu(2, :), NP, 1);
pos = vi > xu;  % NP*D找到所有违反约束上界的值的位置
vi(pos) = (pop(pos) + xu(pos)) / 2;

参考(Reference)

[1] R. Storn and K. Price, “Differential evolution a simple and efficient heuristic for global optimization over continuous spaces,” J. Global Optimization, vol. 11, no. 4, pp. 341–359, 1997.
[2] K. V. Price, R. M. Storn, and J. A. Lampinen, Differential Evolution: A Practical Approach to Global Optimization, 1st ed. New York: SpringerVerlag, Dec. 2005.
[3] R. Joshi and A. C. Sanderson, “Minimal representation multisensor fusion using differential evolution,” IEEE Trans. Syst., Man Cybern. Part A, vol. 29, no. 1, pp. 63–76, Jan. 1999.
[4] J. Zhang, V. Avasarala, and R. Subbu, “Evolutionary optimization of transition probability matrices for credit decision-making,” Eur. J. Oper. Res., to be published.
[5] J. Zhang, V. Avasarala, A. C. Sanderson, and T. Mullen, “Differential evolution for discrete optimization: An experimental study on combinatorial auction problems,” in Proc. IEEE World Congr. Comput. Intell., Hong Kong, China, Jun. 2008, pp. 2794–2800.
[6] R. Gamperle, S. D. Muller, and P. Koumoutsakos, “A parameter study for differential evolution,” in Proc. Advances Intell. Syst., Fuzzy Syst., Evol. Comput., Crete, Greece, 2002, pp. 293–298.
[7] J. Zhang and A. C. Sanderson, “An approximate Gaussian model of differential evolution with spherical fitness functions,” in Proc. IEEE Congr. Evol. Comput., Singapore, Sep. 2007, pp. 2220–2228.
[8] E. Mezura-Montes, J. Velázquez-Reyes, and C. A. Coello Coello, “A comparative study of differential evolution variants for global optimization,” in Proc. Genetic Evol. Comput. Conf., Seattle, WA, Jul. 2006, pp. 485–492.
[9] H. A. Abbass, “The self-adaptive pareto differential evolution algorithm,” in Proc. IEEE Congr. Evol. Comput., vol. 1. Honolulu, HI, May 2002, pp. 831–836.
[10] J. Teo, “Exploring dynamic self-adaptive populations in differential evolution,” Soft Comput.: Fusion Found., Methodologies Applicat., vol. 10, no. 8, pp. 673–686, 2006.
[11] J. Liu and J. Lampinen, “A fuzzy adaptive differential evolution algorithm,” Soft Comput.: Fusion Found., Methodologies Applicat., vol. 9, no. 6, pp. 448–462, 2005.
[12] F. Xue, A. C. Sanderson, P. P. Bonissone, and R. J. Graves, “Fuzzy logic controlled multiobjective differential evolution,” in Proc. IEEE Int. Conf. Fuzzy Syst., Reno, NV, Jun. 2005, pp. 720–725.
[13] A. K. Qin and P. N. Suganthan, “Self-adaptive differential evolution algorithm for numerical optimization,” in Proc. IEEE Congr. Evol. Comput., vol. 2. Sep. 2005, pp. 1785–1791.
[14] V. L. Huang, A. K. Qin, and P. N. Suganthan, “Self-adaptive differential evolution algorithm for constrained real-parameter optimization,” in Proc. IEEE Congr. Evol. Comput., Jul. 2006, pp. 17–24. [15] J. Brest, S. Greiner, B. Boskovic, M. Mernik, and V. Zumer, “Selfadapting control parameters in differential evolution: A comparative study on numerical benchmark problems,” IEEE Trans. Evol. Comput., vol. 10, no. 6, pp. 646–657, Dec. 2006.
[16] J. Brest, V. Zumer, and M. S. Maucec, “Self-adaptive differential evolution algorithm in constrained real-parameter optimization,” in Proc. IEEE Congr. Evol. Comput., Vancouver, BC, Jul. 2006, pp. 215–222.
[17] J. Brest, B. Boskovic, S. Greiner, V. Zumer, and M. S. Maucec, “Performance comparison of self-adaptive and adaptive differential evolution algorithms,” Soft Comput.: Fusion Found., Methodologies Applicat., vol. 11, no. 7, pp. 617–629, 2007.
[18] Z. Yang, K. Tang, and X. Yao, “Self-adaptive differential evolution with neighborhood search,” in Proc. IEEE Congr. Evol. Comput., Hong Kong, China, Jun. 2008, pp. 1110–1116.
[19] P. J. Angeline, “Adaptive and self-adaptive evolutionary computations,” in Computational Intelligence: A Dynamic Systems Perspective. 1995, pp. 152–163.
[20] A. E. Eiben, R. Hinterding, and Z. Michalewicz, “Parameter control in evolutionary algorithms,” IEEE Trans. Evol. Comput., vol. 3, no. 2, pp. 124–141, Jul. 1999.
[21] A. E. Eiben and J. E. Smith, Introduction to Evolutionary Computing. Natural Computing. New York: Springer, 2003.
[22] J. H. Holland, Adaptation in Natural and Artificial Systems. Ann Arbor, MI: The University of Michigan Press, 1975.
[23] T. Back, Evolutionary Algorithms in Theory and Practice: Evolution Strategies, Evolutionary Programming, Genetic Algorithms. New York: Oxford University Press, 1996.
[24] R. Mendes, I. Rocha, E. C. Ferreira, and M. Rocha, “A comparison of algorithms for the optimization of fermentation processes,” in Proc. IEEE Congr. Evol. Comput., Vancouver, BC, Jul. 2006, pp. 2018–2025.
[25] R. C. Eberhart, Y. Shi, and J. Kennedy, Swarm Intelligence. 1st ed., San Mateo, CA: Morgan Kaufmann, Mar. 2001.
[26] J. Zhang and A. C. Sanderson, “JADE: Self-adaptive differential evolution with fast and reliable convergence performance,” in Proc. IEEE Congr. Evol. Comput., Singapore, Sep. 2007, pp. 2251–2258.
[27] J. Tvrdik, I. Krivy, and L. Misik, “Evolutionary algorithm with competing heuristics,” in Proc. MENDEL 2001, Int. Conf. Soft Computing, Brno, Czech, Jun. 2001, pp. 58–64.
[28] J. Tvrdík, L. Misik, and I. Krivy. “Competing heuristics in evolutionary algorithms,” Intell. Technol. Theory Applicat., pp. 159–165, 2002.
[29] X. Yao, Y. Liu, and G. Lin, “Evolutionary programming made faster,” IEEE Trans. on Evol. Comput., vol. 3, no. 2, pp. 82–102, Jul. 1999.
[30] C. Y. Lee and X. Yao, “Evolutionary programming using mutations based on the Lévy probability distribution,” IEEE Trans. Evol. Comput., vol. 8, no. 1, pp. 1–13, Feb. 2004.
[31] B. V. Babu and M. M. L. Jehan, “Differential evolution for multiobjective optimization,” in Proc. IEEE Congr. Evol. Comput., Dec. 2003, pp. 2696–2703.
[32] U. Pahner and K. Hameyer, “Adaptive coupling of differential evolution and multiquadrics approximation for the tuning of the optimization process,” IEEE Trans. Magnetics, vol. 36, no. 4, pp. 1047–1051, Jul. 2000.
[33] E. Mezura-Montes, J. Velazquez-Reyes, and C. A. Coello Coello, “Modified differential evolution for constrained optimization,” in Proc. IEEE Congr. Evol. Comput., Vancouver, BC, Jul. 2006, pp. 25–32.
[34] H. G. Beyer, Theory of Evolution Strategies. New York: Springer-Verlag, Apr. 2001.
[35] X. Yao, Y. Liu, K.-H. Liang, and G. Lin. “Fast evolutionary algorithms,” in Proc. Advances Evol. Computing: Theory Applicat., New York, 2003, pp. 45–94.
[36] L. C. W. Dixon and G. Szegö, “The global optimization problem: An introduction,” in Proc. Toward Global Optimization 2, Amsterdam, Netherlands: North-Holland, 1978, pp. 1–15.
[37] I. C. Trelea, “The particle swarm optimization algorithm: Convergence analysis and parameter selection,” Inform. Process. Lett., vol. 85, no. 6, pp. 317–325, 2003.
[38] J. Vesterstroem and R. Thomsen, “A comparative study of differential evolution, particle swarm optimization, and evolutionary algorithms on numerical benchmark problems,” in Proc. IEEE Congr. Evol. Comput., Jun. 2004, pp. 1980–1987.
[39] Y.-W. Shang and Y.-H. Qiu, “A note on the extended rosenbrock function,” Evol. Comput., vol. 14, no. 1, pp. 119–126, 2006.
[40] Z. Yang, J. He, and X. Yao, “Making a difference to differential evolution,” in Proc. Advances Metaheuristics Hard Optimization, Dec. 2007, pp. 397–414.
[41] Z. Yang, K. Tang, and X. Yao, “Large scale evolutionary optimization using cooperative coevolution,” Inform. Sci., vol. 178, no. 15, pp. 29852999, Feb. 2008.
[42] O. Olorunda and A. P. Engelbrecht, “Differential evolution in highdimensional search spaces,” in Proc. 2007 IEEE Congr. Evol. Comput., Singapore, Sep. 2007, pp. 1934–1941.
[43] J. Zhang and A. C. Sanderson, “Adaptive differential evolution-A robust approach to multimodel problem optimization,” in Series of Adaptation, Learning, and Optimization, New York: Springer-Verlag, Aug. 2009.
[44] Z. Yang, J. Zhang, K. Tang, X. Yao, and A. C. Sanderson, “An adaptive coevolutionary differential evolution algorithm for large-scale optimization,” in Proc. IEEE Congr. Evol. Comput., May 2009, pp. 102–109.
[45] C. A. Coello Coello, “20 years of evolutionary multiobjective optimization: What has been done and what remains to be done,” in Computational Intelligence: Principles and Practice, Los Alamitos, CA: IEEE Computational Intelligence Society, 2006, pp. 73–88.
[46] J. Zhang and A. C. Sanderson, “Self-adaptive multiobjective differential evolution with direction information provided by archived inferior solutions,” in Proc. IEEE World Congr. Evol. Comput., Hong Kong, China, Jun. 2008, pp. 2801–2810.

更多推荐