蝙蝠算法学习

2018-12-17 本文已影响0人星辰研创

序

因需要研究一下这个蝙蝠算法，Bat Algorithm，遂有了这篇学习博客，也作为我的第一篇博客吧，谈谈这方面的知识，练练手。首先，了解到，蝙蝠算法是由Yang教授于2010年提出的高效生物启发式算法，一种搜索全局最优解的算法。该算法是基于迭代的优化技术，初始化为一组随机解，然后通过迭代搜寻最优解，且在最优解周围通过随机飞行产生局部新解，加强了局部搜索。与其他算法相比，BA 在准确性和有效性方面远优于其他算法，且没有许多参数要进行调整。（以上部分摘自百度百科）

什么是蝙蝠算法

描述

简单地说，蝙蝠算法就是模拟蝙蝠回声发射与检测这样的一个机制。一般来说，仿生算法要抽象出来，需要做一些必要的假设和简化：

所有的蝙蝠都使用回声定位来感知距离，并且可以判断出是food还是障碍物
在一定位置以一定速度随机飞行，且可以自动调整发射脉冲的频率或波长，并依据距离调整脉冲发射率
假设响度从一个正值变化到最小值
估计时延和三维地形时不使用射线追踪
频率f在[ $f_{min}$ , $f_{max}$ ]范围内，对应的波长λ在[ $\lambda_{min}$ , $\lambda_{max}$ ]范围内
Tips：对于给定的问题，应该考虑改变波长λ或频率f，λ和f是相关的，λf是一个定值。
根据2010年Yang的文章，有一些公式：

我们需要在每一时间步长t内模拟蝙蝠的位置和速度更新，于是有了如下公式：

$f_{i}=f_{min}+(f_{max}-f_{min})\beta$
$v_{i}^{t+1}=v_{i}^{t}+(x_{i}^{t}-x_{*})f_{i}$
$x_{i}^{t+1}=x_{i}^{t}+v_{i}^{t+1}$
式中， $\beta\in[0,1]$ 是服从均匀分布的随机变量， $x_{*}$ 是全局最优解。

确定一个解之后，使用随机游走产生一个新解：

$x_{new}=x_{old}+\epsilon*A^{(t)}$
式中， $\epsilon\in[-1,1]$ 是一个随机数， $A^{(t)}$ 是当前时步内所有蝙蝠的平均响度。在实现时，提供一个缩放参数来控制步长：
$x_{new}=x_{old}+\sigma\epsilon_{t}*A^{(t)}$
式中， $\epsilon_{t}$ 服从高斯正态分布 $N[0,1]$ ， $\sigma$ 是缩放因子。

此外，响度和脉冲发射率也需要更新：

$A_{i}^{t+1}=\alpha A_{i}^{t}$
$r_{i}^{t+1}=r_{i}^{0}[1-\exp(-\gamma t)]$
式中， $\alpha$ 和 $\gamma$ 是常数，对于任意的 $0<\alpha<1$ 、 $\gamma>0$ ，有：
$A_{i}^{t} \rightarrow 0,r_{i}^{t} \rightarrow r_{t}^{0}, t \rightarrow \infty$

伪代码

初始化种群 $x_{i}$ 和 $v_{i}$
初始化频率 $f_{i}$ 、脉冲发射率 $r_{i}$ 及响度 $A_{i}$
For 1:MAX_ITER
通过调整频率产生新解
根据上述公式更新速度与位置
if rand> $r_{i}$
从最佳解中选择一个并产生局部解
end if
随机飞行产生新解
if $(rand<A_{i})$ && $(f(x_{i})<f(x_{*}))$
接受新解并根据公式更新响度和脉冲发射率
end if
找出当前最佳解 $x_{*}$
End For

如何用编程来实现

这里使用Python语言来描述蝙蝠算法，在程序中面向对象的思想将算法部分编为一个类，并利用算法优化一个简单的例子。

首先是算法类，

import numpy as np

'''
蝙蝠算法-Bat Algorithm
'''


class BA(object):
    def __init__(self, d, N_p, N_gen, Qmin, Qmax, lower_bound, upper_bound, func):
        self.d = d  # 搜索维度
        self.N_p = N_p  # 个体数
        self.N_gen = N_gen  # 迭代次数
        self.A = 1 + np.random.random(self.N_p)  # 响度
        self.r = np.random.random(self.N_p)  # 脉冲发射率
        self.Qmin = Qmin  # 最小频率
        self.Qmax = Qmax  # 最大频率
        self.lower_bound = lower_bound  # 搜索区间下限
        self.upper_bound = upper_bound  # 搜索区间上限
        self.func = func  # 目标函数
        self.alpha = 0.85
        self.gamma = 0.9
        self.r0 = self.r

        self.Lb = self.lower_bound * np.ones(self.d)
        self.Ub = self.upper_bound * np.ones(self.d)
        self.Q = np.zeros(self.N_p)  # 频率
        self.v = np.zeros((self.N_p, self.d))  # 速度

        self.sol = np.zeros((self.N_p, self.d))  # 种群
        self.fitness = np.zeros(self.N_p)  # 个体适应度
        self.best = np.zeros(self.d)  # 最好的solution

        self.fmin = 0.0  # 最小fitness

    # 初始化蝙蝠种群
    def init_bat(self):
        for i1 in range(self.N_p):
            self.sol[i1, :] = self.Lb + (self.Ub - self.Lb) * np.random.uniform(0, 1, self.d)
            self.fitness[i1] = self.func(self.sol[i1, :])
        self.fmin = np.min(self.fitness)
        fmin_arg = np.argmin(self.fitness)
        self.best = self.sol[fmin_arg, :]

    # 越界检查
    def simplebounds(self, s, lb, ub):
        for j1 in range(self.d):
            if s[j1] < lb[j1]:
                s[j1] = lb[j1]
            if s[j1] > ub[j1]:
                s[j1] = ub[j1]
        return s

    # 迭代部分
    def start_iter(self):
        S = np.zeros((self.N_p, self.d))
        self.init_bat()
        for step in range(self.N_gen):
            for i2 in range(self.N_p):
                self.Q[i2] = self.Qmin + (self.Qmin - self.Qmax) * np.random.uniform(0, 1)
                self.v[i2, :] = self.v[i2, :] + (self.sol[i2, :] - self.best) * self.Q[i2]
                S[i2, :] = self.sol[i2, :] + self.v[i2, :]
                S[i2, :] = self.simplebounds(S[i2, :], self.Lb, self.Ub)  # 越界检查

                if np.random.random() > self.r[i2]:
                    S[i2, :] = self.best + 0.001 * np.random.randn(self.d)  # 此处没有实现乘以响度平均值
                    S[i2, :] = self.simplebounds(S[i2, :], self.Lb, self.Ub)  # 越界检查

                Fnew = self.func(S[i2, :])
                if (Fnew <= self.fitness[i2]) and (np.random.random() < self.A[i2]):
                    self.sol[i2, :] = S[i2, :]
                    self.fitness[i2] = Fnew
                    self.A[i2] = self.alpha * self.A[i2]  # 响度更新
                    self.r[i2] = self.r0[i2] * (1 - np.exp(-1 * self.gamma * step))  # 脉冲发射率更新

                if Fnew <= self.fmin:
                    self.best = S[i2, :]
                    self.fmin = Fnew
            print(step, ':', '\n', 'BEST=', self.best, '\n', 'min of fitness=', self.fmin)
        return self.best, self.fmin

这里用一个简单的Griewan函数来测试蝙蝠算法，此函数在 $(0,0...0)$ 有全局极小值0。

import numpy as np
from ba import BA


# 测试用例,Griewan函数，x=(0, 0...,0)处有全局极小值
def func(x):
    y1 = 1 / 4000 * sum(np.power(x, 2))
    y2 = 1
    for h in range(x.size):
        y2 = y2 * np.cos(x[h] / np.sqrt(h + 1))
    y = y1 - y2 + 1
    return y


if __name__ == '__main__':
    ba = BA(10, 20, 100, 0, 2, -2, 2, func)
    best, fmin = ba.start_iter()
    print('=============================================')
    print('BEST=', best, '\n', 'min of fitness=', fmin)

蝙蝠算法学习

序

什么是蝙蝠算法

描述

伪代码

如何用编程来实现

猜你喜欢

热点阅读