启发式算法回顾(自学使用,侵权秒删

一、粒子群算法

1.算法定义

        粒子群算法,也称粒子群优化算法或鸟群觅食算法(Particle Swarm Optimization),缩写为 PSO, 是近年来由J. Kennedy和R. C. Eberhart等开发的一种新的进化算法(Evolutionary Algorithm - EA)。PSO 算法属于进化算法的一种,和模拟退火算法相似,它也是从随机解出发,通过迭代寻找最优解,它也是通过适应度来评价解的品质,但它比遗传算法规则更为简单,它没有遗传算法的“交叉”(Crossover) 和“变异”(Mutation) 操作,它通过追随当前搜索到的最优值来寻找全局最优。

2.算法优势

        其迭代过程受自身经验、种群经验、速度惯性影响,这种算法以其实现容易、精度高、收敛快等优点引起了学术界的重视,并且在解决实际问题中展示了其优越性,在高维特征的收敛上有较好的表现。

代码实现(随机权重改进):

# 第一步,绘制函数图像
import numpy as np
import matplotlib.pyplot as plt
import mpl_toolkits.mplot3d

def func(x,y):
    return x**2+y**2-x*y-10*x-4*y+60

x0 = np.linspace(-15,15,100)
y0 = np.linspace(-15,15,100)
x0,y0 = np.meshgrid(x0,y0)
z0 = func(x0,y0)
fig = plt.figure(constrained_layout=True)
ax = fig.add_subplot(projection='3d')
ax.plot_surface(x0,y0,z0,cmap=plt.cm.viridis,alpha=0.7)
#ax.plot_wireframe(x0,y0,z0) # 另一种绘图方式
ax.set_title('$y = x_1^2+x_2^2-x_1x_2-10x_1-4x_2+60$')

# 第二步,设置粒子群算法的参数
n = 30 # 粒子数量
narvs = 2 # 变量个数
c1 = 2 # 个体学习因子
c2 = 2 # 社会学习因子
w_max = 0.9 # 惯性权重
w_min = 0.4
K = 40 # 迭代次数
vxmax = np.array([(15-(-15))*0.2,(15-(-15))*0.2]) # 粒子在x方向的最大速度
x_lb = np.array([-15,-15]) # x和y的下界
x_ub = np.array([15,15]) # x和y的上界

# 第三步,初始化粒子   
x = x_lb + (x_ub-x_lb)*np.random.rand(n,narvs)
v = -vxmax + 2*vxmax*np.random.rand(n,narvs)

# 第四步,计算适应度
fit = func(x[:,0],x[:,1]) # 计算每个粒子的适应度
pbest = x # 初始化这n个例子迄今为止找到的最佳位置
ind = np.argmax(fit) # 找到适应度最大的那个粒子的下标
gbest = x[ind,:]
gbest_total = np.zeros(K)

# 第五步,更新粒子速度和位置
for j in range(K): # 外层循环,共K次
    w = np.random.rand() # 随机数作为惯性权重
    for p in range(n):
        v[p,:] = w*v[p,:] + c1*np.random.rand(narvs)*(pbest[p,:]-x[p,:]) + c2*np.random.rand(narvs)*(gbest-x[p,:])
    loc_v = np.where(v<-vxmax)
    v[loc_v] = -vxmax[loc_v[1]] # 速度小于-vmax的元素赋值为-vmax
    loc_v = np.where(v>vxmax)
    v[loc_v] = vxmax[loc_v[1]] # 速度大于vmax的元素赋值为vmax
    x = x + v # 更新第i个粒子的位置
    loc_x = np.where(x<x_lb)
    x[loc_x] = x_lb[loc_x[1]]
    loc_x = np.where(x>x_ub)
    x[loc_x] = x_ub[loc_x[1]]
    
    # 第六步,重新计算适应度并找到最优粒子
    fit = func(x[:,0],x[:,1]) # 重新计算n个粒子的适应度
    for k in range(n): # 更新第k个粒子迄今为止找到的最佳位置
        if fit[k]<func(pbest[k,0],pbest[k,1]):
            pbest[k,:] = x[k,:]
    if np.min(fit)<func(gbest[0],gbest[1]): # 更新所有粒子迄今找到的最佳位置
        gbest = x[np.argmin(fit),:]
    gbest_total[j] = func(gbest[0],gbest[1])

ax.scatter(x[:,0],x[:,1],fit,c='r',marker='x')
fig,ax = plt.subplots()
ax.plot(np.arange(K),gbest_total)
plt.rcParams['font.sans-serif']=['SimHei'] #用来正常显示中文标签
plt.rcParams['axes.unicode_minus'] = False #用来显示负号
ax.set_title('算法找到的最小值随优化次数的变化曲线')
ax.set_xlabel('循环次数')
ax.autoscale(axis='x',tight=True)
print('找到的最优解为:',func(gbest[0],gbest[1]))

 原文链接:改进的粒子群算法(PSO)算法Python代码——随机惯性权重_随即权重的pso-CSDN博客

二、遗传算法

1. 算法定义

遗传算法是一种模拟自然选择和遗传机制的优化算法,通过模拟基因的变异、交叉和选择等操作,逐代演化产生新的解,最终找到全局最优解。

2. 遗传算法的基本步骤

遗传算法的基本步骤包括:

  1. 初始化种群: 随机生成初始解的种群。
  2. 适应度评估: 计算每个个体的适应度,即解的优劣程度。
  3. 选择操作: 根据适应度选择个体,将适应度高的个体更有可能被选中。
  4. 交叉操作: 选中的个体进行基因交叉,产生新的个体。
  5. 变异操作: 对新个体进行基因变异,引入新的基因信息。
  6. 更新种群: 根据选择、交叉和变异等操作更新种群。
  7. 重复迭代: 重复进行选择、交叉、变异等操作,直到满足停止条件。
3. 个体的编码方法

在遗传算法中,个体的编码方式通常包括二进制编码、实数编码、排列编码等。选择适当的编码方式取决于具体问题的特点。

4.代码实现

import matplotlib.pyplot as plt
import numpy as np
import random
 
 
# (-1, 2)
# 初始化原始种群
def ori_popular(num):
    popular = []
    for i in range(num):
        x = random.uniform(-1, 2)  # 在此范围内生成一个随机浮点数
        popular.append(x)
    return popular
 
 
# 编码,也就是由表现型到基因型,性征到染色体
def encode(popular):  # popular应该是float类型的列表
    popular_gene = []
    for i in range(0, len(popular)):
        data = int((popular[i]-(-1)) / 3 * 2**18)  # 染色体序列为18bit
        bin_data = bin(data)  # 整形转换成二进制是以字符串的形式存在的
        for j in range(len(bin_data)-2, 18):  # 序列长度不足补0
            bin_data = bin_data[0:2] + '0' + bin_data[2:]
        popular_gene.append(bin_data)
    return popular_gene
 
 
# 解码,即适应度函数。通过基因,即染色体得到个体的适应度值
def decode(popular_gene):
    fitness = []
    for i in range(len(popular_gene)):
        x = (int(popular_gene[i], 2) / 2**18) * 3 - 1
        value = x * np.sin(10 * np.pi * x) + 2
        fitness.append(value)
    return fitness
 
 
# 选择and交叉。选择用轮牌赌,交叉概率为0.66
def choice_ex(popular_gene):
    fitness = decode(popular_gene)
    sum_fit_value = 0
    for i in range(len(fitness)):
        sum_fit_value += fitness[i]
    # 各个个体被选择的概率
    probability = []
    for i in range(len(fitness)):
        probability.append(fitness[i]/sum_fit_value)
    # 概率分布
    probability_sum = []
    for i in range(len(fitness)):
        if i == 0:
            probability_sum.append(probability[i])
        else:
            probability_sum.append(probability_sum[i-1] + probability[i])
 
    # 选择
    popular_new = []
    for i in range(int(len(fitness)/2)):
        temp = []
        for j in range(2):
            rand = random.uniform(0, 1)  # 在0-1之间随机一个浮点数
            for k in range(len(fitness)):
                if k == 0:
                    if rand < probability_sum[k]:
                        temp.append(popular_gene[k])
                else:
                    if (rand > probability_sum[k-1]) and (rand < probability_sum[k]):
                        temp.append(popular_gene[k])
 
        # 交叉,交叉率为0.66。
        is_change = random.randint(0, 2)
        if is_change:
            temp_s = temp[0][9:15]
            temp[0] = temp[0][0:9] + temp[1][9:15] + temp[0][15:]
            temp[1] = temp[1][0:9] + temp_s + temp[1][15:]
 
        popular_new.append(temp[0])
        popular_new.append(temp[1])
    return popular_new
 
 
# 变异.概率为0.05
def variation(popular_new):
    for i in range(len(popular_new)):
        is_variation = random.uniform(0, 1)
        # print([len(k) for k in popular_new])
        if is_variation < 0.02:
            rand = random.randint(2, 19)
            if popular_new[i][rand] == '0':
                popular_new[i] = popular_new[i][0:rand] + '1' + popular_new[i][rand+1:]
            else:
                popular_new[i] = popular_new[i][0:rand] + '0' + popular_new[i][rand+1:]
    return popular_new
 
 
if __name__ == '__main__':  # alt+enter
    # 初始化原始种群, 一百个个体
    num = 100
    ori_popular = ori_popular(num)
    # 得到原始种群的基因
    ori_popular_gene = encode(ori_popular)  # 18位基因
    new_popular_gene = ori_popular_gene
    y = []
    for i in range(1000):  # 迭代次数。繁殖1000代
        new_popular_gene = choice_ex(new_popular_gene)  # 选择和交叉
        new_popular_gene = variation(new_popular_gene)  # 变异
        # 取当代所有个体适应度平均值
        new_fitness = decode(new_popular_gene)
        sum_new_fitness = 0
        for j in new_fitness:
            sum_new_fitness += j
        y.append(sum_new_fitness/len(new_fitness))
    
    # 画图
    x = np.linspace(0, 1000, 1000)
    fig = plt.figure()  # 相当于一个画板
    axis = fig.add_subplot(111)  # 坐标轴
    axis.plot(x, y)
    plt.show() 

三、模拟退火算法 

1.算法定义

        用固体退火模拟组合优化问题,将内能E模拟为目标函数值f,温度T演化成控制参数t,即得到解组合优化问题的模拟退火算法:由初始解i和控制参数初值t开始,对当前解重复“产生新解→计算目标函数差→接受或舍弃”的迭代,并逐步衰减t值,算法终止时的当前解即为所得近似最优解。

        退火过程由冷却进度表(Cooling Schedule)控制,包括控制参数的初值t及其衰减因子Δt、每个t值时的迭代次数L和停止条件Tf。而温度的作用就是来计算转移概率P的。当温度每次下降后,转移概率也发生变化,因此在所有温度下迭代L次的结果也都是不相同的。在每个温度下迭代L次来寻找当前温度下的最优解,然后降低温度继续寻找,直到到达终止温度,即转移概率P接近于0。

        

2.退火过程中参数控制
(1)初始的温度T(0)应选的足够高,使的所有转移状态都被接受。初始温度越高,获得高质量的解的概率越大,耗费的时间越长。

(2)退火速率,即温度下降,最简单的下降方式是指数式下降:
T(n) = α \alphaαT(n) ,n =1,2,3,…
其中α \alphaα是小于1的正数,一般取值为0.8到0.99之间。使的对每一温度,有足够的转移尝试,指数式下降的收敛速度比较慢。

(3)终止温度
如果温度下降到终止温度或者达到用户设定的阈值,则退火完成。

import math
from random import random
import matplotlib.pyplot as plt

def func(x, y):                  #函数优化问题
    res= 4*x**2-2.1*x**4+x**6/3+x*y-4*y**2+4*y**4
    return res
#x为公式里的x1,y为公式里面的x2
class SA:
    def __init__(self, func, iter=100, T0=100, Tf=0.01, alpha=0.99):
        self.func = func
        self.iter = iter         #内循环迭代次数,即为L =100
        self.alpha = alpha       #降温系数,alpha=0.99
        self.T0 = T0             #初始温度T0为100
        self.Tf = Tf             #温度终值Tf为0.01
        self.T = T0              #当前温度
        self.x = [random() * 11 -5  for i in range(iter)] #随机生成100个x的值
        self.y = [random() * 11 -5  for i in range(iter)] #随机生成100个y的值
        self.most_best =[]
        """
        random()这个函数取0到1之间的小数
        如果你要取0-10之间的整数(包括0和10)就写成 (int)random()*11就可以了,11乘以零点多的数最大是10点多,最小是0点多
        该实例中x1和x2的绝对值不超过5(包含整数5和-5),(random() * 11 -5)的结果是-6到6之间的任意值(不包括-6和6)
        (random() * 10 -5)的结果是-5到5之间的任意值(不包括-5和5),所有先乘以11,取-6到6之间的值,产生新解过程中,用一个if条件语句把-5到5之间(包括整数5和-5)的筛选出来。
        """
        self.history = {'f': [], 'T': []}

    def generate_new(self, x, y):   #扰动产生新解的过程
        while True:
            x_new = x + self.T * (random() - random())
            y_new = y + self.T * (random() - random())
            if (-5 <= x_new <= 5) & (-5 <= y_new <= 5):  
                break                                  #重复得到新解,直到产生的新解满足约束条件
        return x_new, y_new 

    def Metrospolis(self, f, f_new):   #Metropolis准则
        if f_new <= f:
            return 1
        else:
            p = math.exp((f - f_new) / self.T)
            if random() < p:
                return 1
            else:
                return 0

    def best(self):    #获取最优目标函数值
        f_list = []    #f_list数组保存每次迭代之后的值
        for i in range(self.iter):
            f = self.func(self.x[i], self.y[i])
            f_list.append(f)
        f_best = min(f_list)
        
        idx = f_list.index(f_best)
        return f_best, idx    #f_best,idx分别为在该温度下,迭代L次之后目标函数的最优解和最优解的下标

    def run(self):
        count = 0
        #外循环迭代,当前温度小于终止温度的阈值
        while self.T > self.Tf:       
           
            #内循环迭代100次
            for i in range(self.iter): 
                f = self.func(self.x[i], self.y[i])                    #f为迭代一次后的值
                x_new, y_new = self.generate_new(self.x[i], self.y[i]) #产生新解
                f_new = self.func(x_new, y_new)                        #产生新值
                if self.Metrospolis(f, f_new):                         #判断是否接受新值
                    self.x[i] = x_new             #如果接受新值,则把新值的x,y存入x数组和y数组
                    self.y[i] = y_new
            # 迭代L次记录在该温度下最优解
            ft, _ = self.best()
            self.history['f'].append(ft)
            self.history['T'].append(self.T)
            #温度按照一定的比例下降(冷却)
            self.T = self.T * self.alpha        
            count += 1
            
            # 得到最优解
        f_best, idx = self.best()
        print(f"F={f_best}, x={self.x[idx]}, y={self.y[idx]}")

sa = SA(func)
sa.run()

plt.plot(sa.history['T'], sa.history['f'])
plt.title('SA')
plt.xlabel('T')
plt.ylabel('f')
plt.gca().invert_xaxis()
plt.show()

评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值