0 引言

生物地理学优化 (Biogeography-Based Optimization, BBO) 算法是Simon受到生物地理学理论的启发，在对生物物种迁移数学模型的研究基础上，于2008年提出的一种新的智能优化算法^[1]。该算法的基本思想来源于Wallace^[2]在19世纪提出的生物地理学理论，该理论研究了生物物种栖息地的分布、迁移和灭绝规律。BBO算法自提出之后，已经得到了较广泛的研究，BBO算法的理论分析、改进以及应用都取得了较好的研究成果:文献[3]提出了混合迁移策略对BBO算法进行改进，提高了算法的收敛精度; 文献[4]提出了动态迁移策略，利用自适应约束方法来约束处理约束条件，更充分地利用了个体有效信息; 文献[5-6]利用马尔可夫链模型对迁移模型进行分析，验证不同的迁移率模型对算法性能产生重要影响，并得出符合自然规律的复杂迁移率模型的性能优于简单的线性迁移率模型的性能的结论; 文献[7]提出一种新的变异算法来改进BBO算法，提高了算法的全局搜索与局部搜索能力; 文献[8-9]使用差分变异算法与BBO算法结合，利用差分变异算子的开采能力，有效改善BBO算法变异能力差的情况，避免陷入局部最优，从而增强了算法的全局搜索能力。为了进一步改进BBO算法的探测与搜索能力，本文将人工萤火虫群优化 (Glowworm Swarm Optimization, GSO) 算法的局部决策域策略引入BBO算法的迁移操作；并且，将改进的BBO算法与原始BBO算法及文献[7-8]的改进算法作比较。实验结果证明，本文提出的改进算法无论是探测能力还是搜索能力都能取得更好的结果。

1 生物地理学优化算法

在生态系统中，每一个栖息地是否适合物种生存可以用栖息适宜指数 (Habitat Suitability Index, HSI) 来描述。一个栖息地的适宜指数可以由适宜指数变量 (Suitability Index Variables, SIV) 来表示，该变量用来描述与适宜指数相关的特征，例如栖息地的降水量、湿度、植被多样性等。HSI与栖息地中物种的多样性成正比。对于某一块栖息地，假设它所能容纳的最大物种数量为S_max，物种的迁入率为λ，其最大值为I，迁出率为μ，其最大值为E，并且假设E=I。当该栖息地物种数为0时，显然，此时该栖息地物种的迁入率为最大值I。随着物种的迁入，栖息地中物种数量增加，迁入率随之减小而迁出率随之增大，当栖息地中的物种数量达到所能容纳的上限S_max时，由于栖息地已经达到饱和状态，因此，迁入率下降为0，而此时迁出率达到最大值E。在图 (1) 中，S₁ < S₂，S₁处的HSI小于S₂处的HSI。显然，当栖息地物种数量为S₁时的迁入率λ(S₁) 大于S₂处的迁入率λ(S₂)。类似地，S₁时的迁出率μ(S₁) 小于S₂处的迁入率μ(S₂)。

有了对迁入率与迁出率的描述之后，可以计算某一时间段内栖息地中所存在的物种数量的概率P_s(t+Δt)。它表示t+Δt时刻栖息地包含s个物种的概率。计算公式如式 (1) 所示:

$ {P_s}(t + \Delta t) = {P_s}(t)(1-{\lambda _s}\Delta t-{\mu _s}\Delta t) + {P_{s-1}}{\lambda _{s - 1}}\Delta t + {P_{s + 1}}{\mu _{s + 1}}\Delta t $

(1)

其中：λ_s与μ_s分别表示当栖息地物种数量为s时的迁入率与迁出率。

根据图迁入率与迁出率的线性关系，可得到物种数量为k时的迁入率与迁出率的计算方法，如式 (2) 与式 (3) 所示:

$ {\lambda _k} = I(1-k/n) $

(2)

$ {\mu _k} = E(k/n) $

(3)

其中，n为栖息地可能存在的最大数目物种时的物种数量。由于假设栖息地物种的迁入率最大值与迁出率最大值相等，因此可得到λ_k+μ_k=E，则当Δt→0时，可以将概率公式转化为式 (4)。如果一个栖息地存在物种数量概率较低，则该栖息地容易发生突变行为。也就是说，突变概率函数与该栖息地的存在物种数量概率成反比，因此可以得到相应的突变概率如式 (5) 所示。

$ {P_s} = \left\{ \begin{array}{l} -({\lambda _s} + {\mu _s}){P_s} + {\mu _{s + 1}}{P_{s + 1}}, \;\;S = 0\\ -({\lambda _s} + {\mu _s}){P_s} + {\lambda _{s-1}}{P_{s - 1}} + {\mu _{s + 1}}{P_{s + 1}}, \;1 \le S \le {S_{\max }} - 1\\ - ({\lambda _s} + {\mu _s}){P_s} + {\lambda _{s - 1}}{P_{s - 1}}, \;\;S = {S_{\max }} \end{array} \right. $

(4)

$ {m_s} = {m_{\max }}(1-{P_s}/{P_{\max }}) $

(5)

令随机初始种群用H_i^(g)表示，综上可得到BBO算法的迁移算子与突变算子如算法1与算法2所示。

算法1  迁移算子。

input：initial population

for i=1 to N do

  Select (H_i^(g))    //根据probability∝λ_i选择H_i^(g)

  if rand (0, 1) < λ_i then

    Select (H_j^(g))    //根据probability∝μ_j, j∈[1, N]选择H_j^(g)

    if rand (0, 1) < μ_j then

     H_j^(g)(SIV)←H_i^(g)(SIV)

     end if

  end if

end for

output：solution of immigration ()

算法2  突变算子。

input：initial population

for i=1 to N do

  P_s(λ_i, μ_i)  //根据式 (1) 计算存在物种数量概率P_s

  m_s(P_s)       //根据式 (5) 计算物种突变概率m_s

   Select (H_i^(g))    //根据probability∝m_s选择H_i^(g)

   if H_i^(g) is selected then

    H_i^(g)(SIV′)←H_i^(g)(SIV)    //在H_i^(g)中随机选择一个SIV用SIV′替换

  end if

end for

output：solution of mutation ()

BBO是通过迁移操作实现可行解之间的信息共享的，栖息地中物种的数量正比于HSI。对于不同的栖息地之间的迁移运动而言，一块栖息地发生迁入或迁出行为的概率与其迁入率或迁出率成正比。因此，将栖息地的迁入率的大小作为其是否发生迁入行为的判断标准，进一步通过对其他栖息地的迁出率考察选择从哪块栖息地进行迁出操作。突变操作能够使得低物种数量存在概率 (existing probability of the number of species) 的解具有一定的几率突变为较好的解，防止算法陷入局部最优，可以增加解的多样性。为了防止最优解被执行突变操作，可以采用精英策略，将每一代的最优解都保存到精英表中，参与到下一代的算法执行过程中。如何执行突变操作应该根据具体解决的问题进行设定。

2 人工萤火虫优化算法

人工萤火虫群优化 (Glowworm Swarm Optimization, GSO) 算法是由Krishnanad等^[10]提出的一种仿生新型群智能优化算法。在GSO算法中，群智能体 (swarm of Agents) 随机分布于搜索空间中。这些Agent根据自然界中萤火虫群的活动规律进行建模。每一个Agent都具有一个被称为虫荧光素 (luciferin) 的属性，这一属性用来表示每个智能体所代表的萤火虫的发光量 (luminescent quantity)。每个智能体的虫荧光素的值与其所发出的光的强度成正比，并且它们在变邻域 (variable neighborhood) 中相互影响。该邻域被定义为局部决策域 (local-decision domain)，其范围r_dⁱ由径向传感器范围r_s(0 < r_dⁱ≤r_s) 界定。在搜索空间中，假设存在Agent i与Agent j。如果Agent j在Agent i的邻域范围内，并且Agent j的虫荧光素水平高于Agent i, 则称Agent i将Agent j视为邻居。决策域使得搜索空间中的Agent能够选择可以相互作用的邻居并同时形成不相交的子群。每一个Agent都能够被其虫荧光素水平更高的邻居吸引。

在BBO算法中，迁移操作仅考虑栖息地各自的迁入率与迁出率大小来判断如何执行迁移操作，并没有考虑到栖息地之间的相互影响。所以，本文将GSO算法的局部决策域思想引入BBO算法的迁移操作，能够增强算法的探测能力。GSO算法中的局部决策域策略能够实现决策邻域内个体间吸引关系的定量分析，通过计算栖息地之间的吸引度来判断它们之间的吸引关系，能够在迁入率与迁出率的基础上更精确地执行迁移操作。

3 基于局部决策域策略的改进BBO算法

局部决策域策略的特点适合作为BBO中栖息地之间消息传递的辅助工具。在各个栖息地之间，待迁入的栖息地能够吸引待迁出栖息地的物种。下面首先给出应用于BBO的局部决策域策略相关定义。

定义1 栖息地迁移矩阵。设栖息地数量为n，计算任意两个栖息地之间的相似度S(H_i, H_j)(相似度可以是对称的也可以是不对称的)，表示栖息地H_j吸引物种从栖息地H_i迁出并且迁入到栖息地H_j的可能性，并将H_i称为待迁出栖息地，H_j称为待迁入栖息地。栖息地之间的相似度组成的n×n的相似度矩阵S称为栖息地迁移矩阵。S(H_i, H_j) 按照式 (6) 计算：

$ S({H_i}, {H_j}) =-\left\| {SIV_i^m-SIV_j^m} \right\| $

(6)

定义2 栖息地迁入参考度 (Habit Immigration Reference, HIR)。相似度S(H_k, H_k) 表示栖息地H_k能否成为待迁入栖息地的参考度，并将其记为HIR(k)(k=1, 2，…, n)。HIR(k) 的取值能够影响待迁入栖息地与待迁出栖息地数量的比例，决定了LBBO的迁移操作能否有效地搜索到优化问题的全局最优解。

根据文献[7]对参考度的定义可以得到，对于所有的栖息地，待迁入栖息地的数量正比于HIR(k)。对于任意栖息地H_k，根据式 (1) 可知，HIR(k) 与迁入率λ_k成正比，与迁出率μ_k成反比。在本文中，令S(H_k, H_k) 的取值为栖息地H_k迁入率λ_k与迁出率μ_k的比值λ_i/μ_k，并且在后面的仿真实验中采用HIR(k) 取值为相似度S的中值、均值、中值的一半以及均值的一半等取值方法作比较，讨论如何选择更优的HIR(k) 取值方法。

定义3 栖息地吸引度 (Habit Attraction Reference, HAR)。将GSO算法思想引入BBO算法，GSO中Agent的虫荧光素属性可以视为BBO算法中任一栖息地对其他相邻栖息地物种的吸引强度。令l_Hi(t) 表示任意栖息地H_i在t时刻的HAR取值。在算法初始化时，令每一块栖息地都具有相等的HAR值l₀。HAR的更新公式如式 (7) 所示:

$ {l_{{H_i}}}(t + 1) = (1-\rho ){l_{{H_i}}}(t) + \gamma {S_{t + 1}}({H_i}, {H_j}) $

(7)

其中: l_Hi(t) 表示栖息地H_i在t时刻的吸引度值; ρ(0 < ρ≤1) 为吸引度衰变常数; γ为吸引度增强常数; S_t(H_i, H_j) 表示t时刻栖息地H_i的HIR值，由式 (1) 计算。

栖息地迁入参考度HIR反映一块栖息地适合其他物种迁入的程度，HIR值越大则表明该栖息地适合其他物种迁入的程度越强，也就是说对其他栖息地的吸引度越强。对于任意栖息地H_i与H_j，其吸引度的更新需要考虑当前吸引度与迁入参考度HIR的双重影响。

根据栖息地的HAR能够确定物种的迁移方向，每一块栖息地的物种都会向该栖息地邻域范围内其他吸引度更高的栖息地迁移。对于任意栖息地i，其物种向邻域内栖息地j的迁移概率如式 (8) 所示:

$ {P_{{H_i}{H_j}}}(t) = \frac{{{l_{{H_j}}}(t)-{l_{{H_i}}}(t)}}{{\sum\limits_{k \in {N_i}(t)} {{l_{{H_k}}}(t)-{l_{{H_i}}}(t)} }} $

(8)

其中：j∈N_i(t)，N_i(t)={j:S_t(H_i, H_j) < r_d^H_i(t); l_i(t) < l_j(t)}表示栖息地H_i在t时刻的邻居集合; r_d^H_i(t) 表示t时刻栖息地H_i的邻域半径。令H_i中物种按照概率P_HiHj(t) 迁入栖息地H_j∈N_Hi(t)。

根据式 (8) 可以表明，物种选择待迁入栖息地的概率是与其吸引度成正比的。

定义4 栖息地决策域半径 (Habit Local-decision Domain Range, HLDR)。设栖息地数量为n，任一栖息地H_i的邻居栖息地数量为N_j(t) 且n=j+1。每一块栖息地近在其邻域范围内传递信息。r_d^H_i(t) 表示t时刻栖息地H_i的邻域半径，其更新公式如式 (9) 所示:

$ r_d^{{H_i}}(t + 1) = \min \{ {r_s}, \max \{ 0, r_d^{{H_i}} + \beta ({n_t}-{N_{{H_i}}}(t))\} \} $

(9)

其中: β为常参数，n_t为控制邻居数量的参数。根据局部决策域策略为改进的BBO算法作出定义。

定义5 栖息地迁移策略Ω(l, r)。Ω(l, r):Hⁿ→H用来控制LBBO迁移操作。通过计算l与r求得栖息地吸引度HAR值来判断物种何时从待迁出栖息地H_i迁出并且迁入到待迁入栖息地H_j。

LBBO的栖息地迁移策略执行的基本流程如下所示。

1) 对于任意栖息地，分别根据式 (2)、式 (3) 与式 (6) 计算得到迁入率λ、迁出率μ以及栖息地相似度矩阵S，并根据λ与μ求得栖息地迁入参考度HIR(k)，根据式 (7)、式 (8) 与式 (9) 求得l、r与迁移概率。

2) 根据式 (7) 与式 (8) 以及步骤1) 所求得的HIR(k)、l、r与迁移概率计算并更新每一块每块栖息地的栖息地吸引度HAR。

3) 将HAR值排序，选择取值最大时所对应的栖息地H_i与栖息地H_j。

4) 在H_j的SIV中随机选择SIV_m。

5) 用SIV_m随机替换H_i的一个SIV。

定义6  栖息地突变策略M(λ, μ)。M(λ, μ):H→H通过栖息地物种数量存在概率P_s来确定随机改变栖息地适宜度指数SIV_m的栖息地突变操作，通过式 (4) 来计算。栖息地突变策略基本流程见BBO突变操作。

定义7  LBBO。LBBO={Φ, m, n, λ, μ, Ω, T}:Hⁿ→Hⁿ 是一个八元组，能够将已经初始化的栖息地进行迭代优化。其中Φ=∅→{Hⁿ, HSIⁿ}是初始化一组栖息地，并且计算栖息地的HSI值的函数。T=Hⁿ(True, False) 是算法终止的判定标准。

LBBO的伪代码如算法3所示。

算法3  LBBO算法。

input：initial population

for g=1 to MAX_GEN do

  for i=1 to N do

    Calculate λ_i, μ_i, HIR, HAR

  end for

  for i=1 to N do     //迁移算子

     Select (H_i^(g), H_j^(g))    //根据定义5的策略选择H_i^(g)与H_j^(g)

    H_j^(g)(SIV)←H_i^(g)(SIV)

  end for

  mutation ()    //根据算法2执行突变算子

  for i=1 to N do

    评价H_j^(g)(SIV) 种群中的新个体

  if H_j^(g)(SIV) 优于H_i^(g)(SIV)

     H_i^(g)(SIV)=H_i^(g)(SIV)

       end if

   end for

end for

output：best solutions

4 实验研究 4.1 测试函数

在本文中，通过benchmark函数优化问题来比较LBBO与基本BBO性能。为了提高实验结果的可靠性，选取了在函数极值个数、极值点分布方面有代表性的7个单峰函数Sphere function、Schwefel’s Problem 1.2、Schwefel’s Problem 2.21、Schwefel’s Problem 2.22、Rosenbrock function、Step function、Quartic function与5个多峰函数Schwefel’s Problem 2.26、Rastrigin function、Ackley function、Griewank function、Fletcher-Powell function作为测试函数^[11]。为了描述方便将其标记为F1~F12，如表 1所示。其中，F1与F7~F8为分段函数，F1、F2、F5~F7以及F10与F12为正则函数。

4.2 LBBO算法性能分析 4.2.1 实验参数设置

为了比较LBBO算法与基本BBO算法及两种改进BBO算法的性能，种群规模及执行代数均按照文献[1]进行设置，令栖息地最大物种数N=50;每个函数运行代数MAX_GEN=50。

LBBO、BBO、改进BBO (Improved BBO，IBBO)^[7]及基于差分进化的BBO (Differential Evolution/BBO, DE/BBO)^[8]的参数根据文献[1]设置为：栖息地迁移概率P_mod=1，突变概率m=0，精英栖息地数量K=2，最大迁移率设置为E=I=1，最大突变率m_max，初始HAR l₀=5，初始决策域半径r_d^H_i(t₀)=1;每个函数运行30次取平均值。

4.2.2 实验结果与分析

实验结果如表 2所示。对于单峰函数F1~F7，通过对实验结果的观察可以发现，各项指标LBBO均优于BBO算法及其他两种改进算法。对于函数F5，是非凸的病态函数，其全局最优点位于一个平滑、狭长的抛物线形山谷内，其最优点较难求得。而LBBO采用了局部决策域策略，能够通过不同解之间的相互吸引关系较为准确的判断全局最优解。并且除了函数F3外，LBBO运算30次所得到结果的标准差也小于BBO，表明LBBO在求解函数F1、F2以及F4~F7时更加稳定。对于函数F3，通过对实验结果观察可以发现，LBBO在求解过程中与BBO算法及其他两种算法相比，大部分最优解的取值与均值差距较大。

对于多峰函数F8~F12，优化结果LBBO算法优于BBO算法与其他两种比较算法。其中，函数F9峰形呈高低起伏不定跳跃性的出现, 因此，LBBO中的吸引子策略能够通过迁移信息在解之间的传播而比较有效地发现函数最优解。函数F10十分复杂，优化过程中易落入局部最优的陷阱。LBBO采用局部决策域策略，能够在算法陷入局部最优时使得算法在一定概率下接受较差解，使得算法具有更强的跳出局部最优的能力。因此LBBO能够较准确地得到函数F10优化问题的最优解。函数F11具有许多局部极小值，较难求得全局最优点，LBBO从全局搜索与局部搜索两个方面对该复杂函数进行优化，具有更强的全局寻优能力与跳出局部最优的能力，所以LBBO能够非常有效地求解函数F11这样的复杂函数。

4.3 LBBO的HIR参数取值分析实验 4.3.1 实验设置

LBBO的吸引子策略HIR取值为栖息地迁入率与栖息地迁出率的比值，为了比较HIR不同取值情况下算法的效率设置对比实验。分别令HIR取值为栖息地迁入率与迁出率的比值、栖息地相似度矩阵S的中值 (median)、栖息地相似度矩阵S中值的一半、栖息地相似度矩阵S的均值以及栖息地相似度矩阵均值的一半，采用4.1节所示12种测试函数，通过30次独立测试结果的均值与标准差来比较算法的效率。五种HIR取值方法的LBBO的参数按照4.2.1节进行设置，令栖息地迁移概率P_mod，突变概率m=0，精英栖息地数量K=2，最大迁移率设置为E=I=1，最大突变率m_max，初始HAR l₀=5，初始决策域半径=1;栖息地最大物种数N=50;每个函数运行代数MAX_GEN=50。

4.3.2 实验结果与分析

实验结果如表 3与表 4所示，通过观察表 3数据可以发现在12个函数求解中LBBO (HIR(k)=λ_k/μ_k) 效果要优于其他4种情况。主要原因是栖息地相似度矩阵S仅考虑了单一栖息地对迁移操作的影响，而迁入率与迁出率则综合考虑了两个栖息地之间的相互关系，因此HIR的值设置为栖息地迁入率与迁出率的比值能够更好地利用栖息地物种数量信息，在栖息地之间传递吸引度与归属度信息时能够更加准确地指导栖息地间物种迁移行为。

4.4 LBBO初始种群与迭代次数参数分析实验 4.4.1 实验设置

为了分析算法的初始最大物种数以及算法执行代数对算法效果的影响，在实验设置中均采用较大数值来分析较大取值的物种数及执行代数是否能对算法效果带来明显的影响。采用4.1节所示F1、F5、F12等3个测试函数，通过30次独立测试结果的均值与标准差来比较算法的效率。栖息地最大物种数N=500;每个函数运行代数MAX_GEN=1 000，其他参数设置与4.2.1节相同。

4.4.2 实验结果与分析

对于简单的单峰函数F1与复杂的多峰函数F5与F12增大物种数与执行代数后能够取得最优值。对于测试函数F1，表 4与图 1(a)中显示在算法执行600代之前随着种群规模的增大算法的精确度会略微增加，但是600代之后基本上没有差异。对于F5由图 1(b)可知在算法迭代800代之前初始种群数量对算法的结果有着显著影响。图 1(c)表明相对F5更加复杂的函数F12在迭代1 000代初始种群的规模才开始不明显。根据这一实验结果可以得出，测试函数的复杂程度与其对初始种群规模的敏感正度是成正比的，对于简单函数而言更大的初始种群没有很明显的影响。

对于执行代数而言，从图 1中观察可知在一定数量的代数之后精度几乎不再提升，其主要原因是，算法能否求得最优值主要取决于算法跳出局部最优的能力。算法迭代次数的增加无疑会对算法执行时间造成很大影响，然而，如表 5所示，执行代数相同时初始种群的规模对算法执行时间的影响更加明显。这是因为初始种群规模增加会使得迁入迁出率、HIR等参数的计算次数激增而显著地增加算法执行时间。因此对于算法在实际问题的应用应该通过大量实验考察初始种群规模与迭代次数在何时无法对算法的精度带来提升，进而选择恰当的取值来平衡算法精度与算法效率的关系。

5 结语

本文提出了一种基于局部决策域策略的LBBO算法，采用萤火虫算法局部决策域策略来改进算法的局部搜索能力。通过对LBBO算法在12个基准函数的测试结果分析，该算法能够较好地求解函数极值优化问题。通过与BBO算法的比较，可以发现LBBO的性能明显优于原始算法。当然，BBO算法的研究与其他进化算法 (Evolutionary Algorithm, EA)(例如遗传算法 (Genetic Algorithm, GA)、粒子群优化 (Particle Swarm Optimization, PSO) 算法等) 相比仍处于较初级的研究阶段。在今后的工作中，将尝试利用其他的算法框架模型 (例如文化算法框架模型) 来改进BBO的性能。