python面向对象实现布谷鸟算法（CS）

Nixiawyc

于 2023-03-19 19:07:33 发布

阅读量380

点赞数 2

文章标签： python 算法机器学习

本文链接：https://blog.csdn.net/Nixiawyc/article/details/129654176

版权

以最小化一个m维变量的目标函数min f(X)为例，描述CS算法的实现步骤。

step1:定义目标函数f(X)， $X = (x_{1},...x_{m})^{T}$ ，问题维数m，寄生巢规模Popsize = N，淘汰概率Pa，最大迭代次数iter_max，随机生成N个鸟巢的初始位置 $X_{i},i=1,2,...,N$

step2:计算每个寄生巢的位置的目标函数 $f(X_{i})$ 的值，从中得到当前的最优函数值 $f(X^{*})$ 及最佳位置 $X^{*}$ 。保存 $X^{*}$ 和 $f(X^{*})$ 作为当前全局最佳目标值和目标解。

step3:莱维飞行：利用莱维飞行的公式对寄生巢位置进行更新。莱维飞行_百度百科 (baidu.com)

step4:计算更新位置后的目标函数值并与前期数值进行比较。若更小，则更新当前最优值。

step5:发现丢弃：针对每一个寄生巢，取随机数 $r\epsilon [0,1]$ 。若 r>Pa ，则随机更新该寄生巢位置，否则不变。

step6:若未达到最大迭代次数或最小误差要求，转step3，否则，输出全局最优结果，算法结束。

python面向对象实现布谷鸟算法(CS)

import numpy as np
import scipy.special as sc_special
import time
import matplotlib.pyplot as plt


def fit_func(nest):
    x = nest
    # return 3*(1-x)**2*np.e**(-x**2-(y+1)**2) - 10*(x/5-x**3-y**5)*np.e**(-x**2-y**2) - (np.e**(-(x+1)**2-y**2))/3
    return (1 + pow((1 + x[0] + x[1]), 2) * (19 - 14 * x[0] + 3 * x[0] * x[0] - 14 * x[1] + 6 * x[0] * x[1] +
                                             3 * x[1] * x[1])) * (
                   30 + pow((2 * x[0] - 3 * x[1]), 2) * (18 - 32 * x[0] + 12 * x[0] * x[0] +
                                                         48 * x[1] - 36 * x[0] * x[1] + 27 * x[1] * x[1]))


# 算法执行时间
def wrapper(func):
    def inner(*args, **kwargs):
        start_time = time.time()
        res = func(*args, **kwargs)
        end_time = time.time()
        result = end_time - start_time
        print('执行时间为 %.3fs' % result)
        return res

    return inner


# 生成莱维飞行的步长
def levy_flight(n, m, beta):
    """
    This function implements Levy's flight.
    ---------------------------------------------------
    Input parameters:
        n: Number of steps
        m: Number of dimensions
        beta: Power law index (note: 1 < beta < 2)
    Output:
        'n' levy steps in 'm' dimension
    """
    sigma_u = (sc_special.gamma(1 + beta) * np.sin(np.pi * beta / 2) / (
            sc_special.gamma((1 + beta) / 2) * beta * (2 ** ((beta - 1) / 2)))) ** (1 / beta)
    sigma_v = 1

    u = np.random.normal(0, sigma_u, (n, m))
    v = np.random.normal(0, sigma_v, (n, m))
    steps = u / ((np.abs(v)) ** (1 / beta))

    return steps


# 布谷鸟算法类
class Cuckoo(object):
    def __init__(self, n, m, fit_func, lower_boundary, upper_boundary):
        self.n = n  # 种群数量
        self.m = m  # 变量维数
        self.fit_func = fit_func  # 适应度函数
        self.lower_boundary = lower_boundary  # 所有变量的下界
        self.upper_boundary = upper_boundary  # 所有变量的上界
        self.iter_num = 1000  # 迭代次数， 默认为100
        self.pa = 0.25  # 寄生卵被宿主发现的概率，默认为0.25
        self.beta = 1.5  # 幂律指数，默认为1.5
        self.step_coefficient = 0.1  # 步长缩放因子，与问题的规模有关，默认为0.1
        self.nests = self.generate_nests()  # 所有的个体
        self.fitness = self.calc_fitness(self.nests)  # 所有个体对应的适应值
        # get the best nest and record it
        best_nest_index = np.argmin(self.fitness)  # 最优个体对应的索引
        self.best_fitness = self.fitness[best_nest_index]  # 最佳适应度值
        self.best_nest = self.nests[best_nest_index].copy()  # 最佳适应度值对应的个体
        self.min_record = [self.best_fitness]  # 记录每一代中的最优个体对应的适应度值
        self.best_record = [self.best_fitness]  # 记录全局的最优个体对应的适应度值

    # 产生所有鸟巢的位置
    def generate_nests(self):
        """
        Generate the nests' locations
        ---------------------------------------------------
        Output:
            generated nests' locations
        """
        lower_boundary = np.array(self.lower_boundary)
        upper_boundary = np.array(self.upper_boundary)
        nests = np.empty((self.n, self.m))

        for each_nest in range(self.n):
            nests[each_nest] = lower_boundary + np.array([np.random.rand() for _ in range(self.m)]) * (
                    upper_boundary - lower_boundary)

        return nests

    # 计算每个鸟巢的适应值（即函数值）
    def calc_fitness(self, nests):
        """
        calculate each nest's fitness
        ---------------------------------------------------
        Output:
            Every nest's fitness
        """

        n, m = nests.shape
        fitness = np.empty(n)

        for each_nest in range(n):
            fitness[each_nest] = self.fit_func(nests[each_nest])

        return fitness

    # 更新鸟巢的位置
    def update_nests(self):
        """
        This function is to get new nests' locations and use new better one to replace the old nest
        ---------------------------------------------------
        Output:
            Updated nests' locations
        """
        lower_boundary = np.array(self.lower_boundary)
        upper_boundary = np.array(self.upper_boundary)
        # generate steps using levy flight
        steps = levy_flight(self.n, self.m, self.beta)
        new_nests = self.nests.copy()

        for each_nest in range(self.n):
            # coefficient 0.01 is to avoid levy flights becoming too aggressive
            # and (nest[each_nest] - best_nest) could let the best nest be remained
            step_size = self.step_coefficient * steps[each_nest] * (self.nests[each_nest] - self.best_nest)
            step_direction = np.random.rand(self.m)
            new_nests[each_nest] += step_size * step_direction
            # apply boundary conditions
            new_nests[each_nest][new_nests[each_nest] < lower_boundary] = lower_boundary[
                new_nests[each_nest] < lower_boundary]
            new_nests[each_nest][new_nests[each_nest] > upper_boundary] = upper_boundary[
                new_nests[each_nest] > upper_boundary]

        new_fitness = self.calc_fitness(new_nests)
        self.nests[new_fitness > self.fitness] = new_nests[new_fitness > self.fitness]

    # 按照一定概率抛弃鸟蛋，换巢（局部搜索）
    def abandon_nests(self):
        """
        Some cuckoos' eggs are found by hosts, and are abandoned.So cuckoos need to find new nests.
        ---------------------------------------------------
        Output:
            Updated nests' locations
        """
        lower_boundary = np.array(self.lower_boundary)
        upper_boundary = np.array(self.upper_boundary)
        for each_nest in range(self.n):
            if np.random.rand() < self.pa:
                step_size = np.random.rand() * (
                        self.nests[np.random.randint(0, self.n)] - self.nests[np.random.randint(0, self.n)])
                self.nests[each_nest] += step_size
                # apply boundary conditions
                self.nests[each_nest][self.nests[each_nest] < lower_boundary] = lower_boundary[
                    self.nests[each_nest] < lower_boundary]
                self.nests[each_nest][self.nests[each_nest] > upper_boundary] = upper_boundary[
                    self.nests[each_nest] > upper_boundary]

    # 更新全局最优解
    def update_best(self):
        self.fitness = self.calc_fitness(self.nests)
        min_nest_index = np.argmin(self.fitness)
        min_fitness = self.fitness[min_nest_index]
        self.min_record.append(min_fitness)
        if min_fitness < self.best_fitness:
            self.best_fitness = min_fitness
            self.best_nest = self.nests[min_nest_index].copy()
        self.best_record.append(self.best_fitness)

    def cuckoo_search(self):
        """
        Cuckoo search function
        ---------------------------------------------------
        Output:
            The best solution and its value
        """
        for _ in range(self.iter_num):
            self.update_nests()
            # print('更新鸟巢位置')
            self.abandon_nests()
            # print('换巢')
            self.update_best()
            # print('更新全局最优值')


# 结果可视化
def plot_result(x):
    plt.plot(x)
    plt.show()


@wrapper
def main(n, m, fit_func, lower_boundary, upper_boundary):
    c = Cuckoo(n, m, fit_func, lower_boundary, upper_boundary)
    c.cuckoo_search()
    return c.best_nest, c.best_fitness, c.best_record, c.min_record


if __name__ == '__main__':
    n = 100
    m = 2
    lower_boundary = [-3, -3]
    upper_boundary = [3, 3]
    Result = main(n, m, fit_func, lower_boundary, upper_boundary)
    print("种群数量为{}，迭代次数为{}，最优值为{}，最优个体为{}".format(n, 1000, Result[1], Result[0]))
    print(fit_func(Result[0]))
    # print(Result[2])
    plot_result(Result[2])

文中仅介绍算法实现流程，代码根据参考[1]改写，具体更新公式见参考[2]