基于python的缎蓝园丁鸟优化算法Satin Bowerbird Optimizer(SBO)

# ！usr/bin/env python # -*- coding: utf-8 -*- # Time : 2021/12/3 17:36 # @Author : LucXiong # @Project : Model # @File : SBO.py """ Ref:https://blog.csdn.net/u011835903/article/details/107857884 """ import random # random Function import numpy as np # numpy operations import matplotlib.pyplot as plt import math import test_function class SBO(): def __init__(self, pop_size=50, n_dim=2, alpha=0.94, lb=-1e5, ub=1e5, max_iter=1000, func=None): self.pop = pop_size self.n_dim = n_dim self.alpha = alpha # 步长的最大阈值 self.func = func self.max_iter = max_iter # max iter self.z = 0.02 # z是缩放比例因子 self.r_mutate = 0.05 # 变异概率 self.lb, self.ub = np.array(lb) * np.ones(self.n_dim), np.array(ub) * np.ones(self.n_dim) assert self.n_dim == len(self.lb) == len(self.ub), 'dim == len(lb) == len(ub) is not True' assert np.all(self.ub > self.lb), 'upper-bound must be greater than lower-bound' self.X = np.random.uniform(low=self.lb, high=self.ub, size=(self.pop, self.n_dim)) self.Y = [self.func(self.X[i]) for i in range(len(self.X))] # y = f(x) for all particles self.pbest_x = self.X.copy() # personal best location of every particle in history self.pbest_y = [np.inf for i in range(self.pop)] # best image of every particle in history self.fit = [1 / (1 + self.Y[i]) if self.Y[i] > 0 else 1 - self.Y[i] for i in range(self.pop)] self.prob = [self.fit[i] / sum(self.fit) for i in range(self.pop)] self.gbest_x = self.pbest_x.mean(axis=0).reshape(1, -1) # global best location for all particles self.gbest_y = np.inf # global best y for all particles self.gbest_y_hist = [] # gbest_y of every iteration self.update_gbest() def update_pbest(self): ''' personal best :return: ''' for i in range(len(self.Y)): if self.pbest_y[i] > self.Y[i]: self.pbest_x[i] = self.X[i] self.pbest_y[i] = self.Y[i] def update_gbest(self): ''' global best :return: ''' idx_min = self.pbest_y.index(min(self.pbest_y)) if self.gbest_y > self.pbest_y[idx_min]: self.gbest_x = self.X[idx_min, :].copy() self.gbest_y = self.pbest_y[idx_min] def cal_prob(self): self.fit = [1 / (1 + self.Y[i]) if self.Y[i] > 0 else 1 - self.Y[i] for i in range(self.pop)] self.prob = [self.fit[i] / sum(self.fit) for i in range(self.pop)] def update(self, iter_num): idx_min = self.Y.index(min(self.Y)) for i in range(self.pop): # roulette wheel for k in range(self.n_dim): select_list = [] while len(select_list) < 1: select_list = [] r = np.random.rand(1) for kk in range(len(self.prob)): if self.prob[kk] > (r[0]): select_list.append(kk) j = random.choice(select_list) lemta = self.alpha / (1 + self.prob[j]) try: self.X[i, k] += lemta * 0.5 * (self.X[j, k] + self.gbest_x[0][k] - 2 * self.X[i, k]) except: self.X[i, k] += lemta * 0.5 * (self.X[j, k] + self.gbest_x[k] - 2 * self.X[i, k]) # 变异 if np.random.rand(1) < self.r_mutate and i != idx_min: for j in range(self.n_dim): # 正态分布 Satin bowerbird optimizer: A new optimization algorithm to optimize # ANFIS for software development effort estimation.pdf self.X[i][j] += self.z * np.random.normal() * (self.ub[j] - self.lb[j]) self.X = np.clip(self.X, self.lb, self.ub) self.Y = [self.func(self.X[i]) for i in range(len(self.X))] # def mutate(self): # # 没有用上这个函数 # idx_min = self.Y.index(min(self.Y)) # for i in range(self.pop): # if np.random.rand(1)[0] < self.r_mutate and i != idx_min: # for j in range(self.n_dim): # self.X[i][j] += self.z * np.random.normal() * (self.ub[j] - self.lb[j]) # self.X = np.clip(self.X, self.lb, self.ub) # self.Y = [self.func(self.X[i]) for i in range(len(self.X))] def run(self): for i in range(self.max_iter): print(i) self.update(i) self.cal_prob() self.update_pbest() self.update_gbest() self.gbest_y_hist.append(self.gbest_y) self.best_x, self.best_y = self.gbest_x, self.gbest_y return self.best_x, self.best_y if __name__ == '__main__': # 寻优效果不错，但是计算时间较长，有一部分原因是我频繁整体计算适应度值 n_dim = 30 lb = [-100 for i in range(n_dim)] ub = [100 for i in range(n_dim)] demo_func = test_function.fu5 pop_size = 100 max_iter = 100 sbo = SBO(n_dim=n_dim, pop_size=pop_size, max_iter=max_iter, lb=lb, ub=ub, func=demo_func) best_x, bext_y = sbo.run() print(f'{demo_func(sbo.gbest_x)}\t{sbo.gbest_x}') plt.plot(sbo.gbest_y_hist) plt.show()