粒子群算法
大概思路
代码
# -*- coding: utf-8 -*- """ Created on Wed Jul 20 17:14:56 2022 @author: 83989 """ import math import random import numpy as np import matplotlib.pyplot as plt import pylab as mpl mpl.rcParams['font.sans-serif'] = ['SimHei'] class PSO: def __init__(self, dimension, time, size, low, up, v_low, v_high): # 初始化 self.dimension = dimension # 变量个数 self.time = time # 迭代的代数 self.size = size # 种群大小 self.bound = [] # 变量的约束范围 self.bound.append(low) self.bound.append(up) self.v_low = v_low self.v_high = v_high self.x = np.zeros((self.size, self.dimension)) # 所有粒子的位置 self.v = np.zeros((self.size, self.dimension)) # 所有粒子的速度 self.p_best = np.zeros((self.size, self.dimension)) # 每个粒子最优的位置 self.g_best = np.zeros((1, self.dimension))[0] # 全局最优的位置 # 初始化第0代初始全局最优解 temp = -1000000 for i in range(self.size): for j in range(self.dimension): self.x[i][j] = random.uniform(self.bound[0][j], self.bound[1][j]) self.v[i][j] = random.uniform(self.v_low, self.v_high) self.p_best[i] = self.x[i] # 储存最优的个体 fit = self.fitness(self.p_best[i]) # 做出修改 if fit > temp: self.g_best = self.p_best[i] temp = fit def fitness(self, x): """ 个体适应值计算 """ x1 = x[0] x2 = x[1] x3 = x[2] x4 = x[3] x5 = x[4] y = math.floor((x2 * np.exp(x1) + x3 * np.sin(x2) + x4 + x5) * 100) / 100 # print(y) return y def update(self, size): c1 = 2.0 # 学习因子 c2 = 2.0 w = 0.8 # 自身权重因子 for i in range(size): # 更新速度(核心公式) self.v[i] = w * self.v[i] + c1 * random.uniform(0, 1) * ( self.p_best[i] - self.x[i]) + c2 * random.uniform(0, 1) * (self.g_best - self.x[i]) # 速度限制 for j in range(self.dimension): if self.v[i][j] < self.v_low: self.v[i][j] = self.v_low if self.v[i][j] > self.v_high: self.v[i][j] = self.v_high # 更新位置 self.x[i] = self.x[i] + self.v[i] # 位置限制 for j in range(self.dimension): if self.x[i][j] < self.bound[0][j]: self.x[i][j] = self.bound[0][j] if self.x[i][j] > self.bound[1][j]: self.x[i][j] = self.bound[1][j] # 更新p_best和g_best if self.fitness(self.x[i]) > self.fitness(self.p_best[i]): self.p_best[i] = self.x[i] if self.fitness(self.x[i]) > self.fitness(self.g_best): self.g_best = self.x[i] def pso(self): best = [] self.final_best = np.array([1, 2, 3, 4, 5]) for gen in range(self.time): self.update(self.size) if self.fitness(self.g_best) > self.fitness(self.final_best): self.final_best = self.g_best.copy() print('当前最佳位置:{}'.format(self.final_best)) temp = self.fitness(self.final_best) print('当前的最佳适应度:{}'.format(temp)) best.append(temp) t = [i for i in range(self.time)] plt.figure() plt.grid(axis='both') plt.plot(t, best, color='red', marker='.', ms=10) plt.rcParams['axes.unicode_minus'] = False plt.margins(0) plt.xlabel(u"迭代次数") # X轴标签 plt.ylabel(u"适应度") # Y轴标签 plt.title(u"迭代过程") # 标题 plt.show() if __name__ == '__main__': time = 50 size = 100 dimension = 5 v_low = -1 v_high = 1 low = [1, 1, 1, 1, 1] up = [25, 25, 25, 25, 25] pso = PSO(dimension, time, size, low, up, v_low, v_high) pso.pso()
特点
收敛速度快但容易陷入局部最优解
遗传算法
大概思路
代码
from __future__ import division import numpy as np import matplotlib.pyplot as plt def f(x): return 5*np.sin(6*x) + 6*np.cos(5*x) class GeneticAlgorithm(object): """ 遗传算法求多极小函数最优值 """ def __init__(self, cross_rate, mutation_rate, n_population, n_iterations, DNA_size): self.cross_rate = cross_rate #交配的可能性大小 self.mutate_rate = mutation_rate #基因突变率 self.n_population = n_population #种群大小 self.n_iterations = n_iterations #迭代次数 self.DNA_size = 6 # DNA的长度 self.x_bounder = [-3, 6] #搜索范围 def init_population(self): """ 初始化种群 :return: """ population = np.random.randint(low=0, high=2, size=(self.n_population, self.DNA_size)).astype(np.int8) return population def transformDNA(self, population): """ 编码:十进制转化为二进制 :param population: :return: """ population_decimal = ( (population.dot(np.power(2, np.arange(self.DNA_size)[::-1])) / np.power(2, self.DNA_size) - 0.5) * (self.x_bounder[1] - self.x_bounder[0]) + 0.5 * (self.x_bounder[0] + self.x_bounder[1]) ) return population_decimal def fitness(self, population): """ 适应性函数 :param population: 基因个体 :return: 与当前最优的解 """ transform_population = self.transformDNA(population) fitness_score = f(transform_population) return fitness_score - fitness_score.min() def select(self, population, fitness_score): """ 根据适应值进行选择 :param population: :param fitness_score: :return: """ fitness_score = fitness_score + 1e-4 fitness_score = np.array([1/i for i in fitness_score]) #求最小值,需要反一下,值越小,越适应 idx = np.random.choice(np.arange(self.n_population), size=self.n_population, replace=True, p=fitness_score/fitness_score.sum()) return population[idx] def create_child(self, parent, pop): """ 交叉 :param parent: :param pop: :return: """ if np.random.rand() < self.cross_rate: index = np.random.randint(0, self.n_population, size=1) cross_points = np.random.randint(0, 2, self.DNA_size).astype(np.bool) parent[cross_points] = pop[index, cross_points] return parent def mutate_child(self, child): """ 变异 :param child: :return: """ for i in range(self.DNA_size): if np.random.rand() < self.mutate_rate: child[i] = 1 else: child[i] = 0 return child # 进化 def evolution(self): population = self.init_population() NUM = [] for i in range(self.n_iterations): fitness_score = self.fitness(population) best_person = population[np.argmin(fitness_score)] if (i+1)%100 == 0: NUM.append(f(self.transformDNA(best_person))) print(u'第%-4d次进化后, 最优结果在: %s处取得, 对应的最小值为: %f' % (i, self.transformDNA( best_person), f(self.transformDNA( best_person)))) population = self.select(population, fitness_score) population_copy = population.copy() for parent in population: child = self.create_child(parent, population_copy) child = self.mutate_child(child) parent[:] = child population = population x_num = [i for i in range(len(NUM))] plt.rcParams['font.sans-serif'] = ['SimHei'] plt.rcParams['axes.unicode_minus'] = False plt.figure(1) plt.title('遗传算法') plt.ylabel('适应度(最小值)') plt.xlabel('进化次数(*100)') plt.plot(x_num, NUM, 'g') plt.show() self.best_person = best_person return self.transformDNA(best_person) def main(): ga = GeneticAlgorithm(cross_rate=0.9, mutation_rate=0.1, n_population=200, n_iterations=2001, DNA_size=8) best = ga.evolution() # 绘图 plt.figure(2) x = np.linspace(start=ga.x_bounder[0], stop=ga.x_bounder[1], num=200) plt.plot(x, f(x)) plt.scatter(ga.transformDNA(ga.best_person), f(ga.transformDNA(ga.best_person)), s=200, lw=0, c='red', alpha=0.5) ax = plt.gca() ax.spines['right'].set_color('none') ax.spines['top'].set_color('none') ax.xaxis.set_ticks_position('bottom') ax.yaxis.set_ticks_position('left') plt.show() return best if __name__ == '__main__': Bests = [] for _ in range(1): Bests.append(f(main())) Avg = np.mean(Bests) Best = np.min(Bests) Worst = np.max(Bests) Var = np.var(Bests) print("结果:", Bests) print("\nAVG: ", Avg, "\nBest: ", Best, "\nWorst: ", Worst, "\nVar: ", Var)
特点
算法能在较短的时间当中得到最优值,且种群的概念使得算法有并行的模式,多的个体同时搜索最优值。但算法的实现比较麻烦,设计到问题的编码和解码,且参数较多,不同的参数选择选择也没有理论的好坏,只能通过多次尝试。同时由于遗传算法的迭代眼里,使得初始解对最终解的影响比较大。
模拟退火算法
大概思路
代码
import matplotlib.pyplot as plt import numpy as np import math import time np.random.seed(int(time.time())) num = 300 K = 0.1 R = 0.9 #控制温度降低快慢 T_max = 30 #初始温度 T_min = 0.1 #下限温度 def Func(x): return 5*np.sin(6*x) + 6*np.cos(5*x) def main(): x = np.random.uniform(0,2*math.pi) Best_A =Func(x) Best_array = [] T_array = [] T = T_max while(T >T_min): for i in range(num): x_temp = np.random.uniform(0,2*math.pi) current = Func(x_temp) dE = Best_A - current if dE>=0: Best_A = current x = x_temp else: if math.exp(dE/T) > np.random.uniform(0,1): Best_A = current x = x_temp T = R*T T_array.append(T) Best_array.append(Best_A) Plot(Best_array,T_array) # print("最小值 : ",Best_A) return Best_A,x def Plot(num,T_array): plt.figure(1) x_num = [i for i in range(len(num))] plt.rcParams['font.sans-serif'] = ['SimHei'] plt.rcParams['axes.unicode_minus'] = False plt.subplot(211) plt.title('模拟退火算法') plt.ylabel('目标函数') plt.plot(x_num,num,'g') plt.subplot(212) plt.ylabel('温度') plt.plot(x_num,T_array,'r') plt.show() if __name__ == '__main__': N = 1 Bests = [] xs = [] for _ in range(N): y,x = main() Bests.append(y) xs.append(x) Avg = np.mean(Bests) Best = np.min(Bests) Worst = np.max(Bests) Var = np.var(Bests) print("结果:",Bests,"坐标:",xs) print("\nAVG: ",Avg, "\nBest: ",Best, "\nWorst: ",Worst, "\nVar: ",Var) plt.figure(2) x = np.linspace(start=0, stop=2*math.pi, num=200) plt.plot(x, Func(x)) plt.scatter(xs, Bests, s=200, lw=0, c='red', alpha=0.5) ax = plt.gca() ax.spines['right'].set_color('none') ax.spines['top'].set_color('none') ax.xaxis.set_ticks_position('bottom') ax.yaxis.set_ticks_position('left') plt.show()
特点
模拟退火算法计算过程比较简单,且鲁棒性强,多次实验都能达到较好的结果,且结果方差小,统计结果都有较好的表现。同时也发现了模拟退火算法的一些缺点,如果温度下降较快,很容易得到的不是全局最优解。为了得到全局最优解,就得使温度下降的足够慢,这样导致执行时间比较长。