【遗传算法GA】--计算函数最值（Python）

文章目录

1.遗传算法概述
2.分步实现过程
3.完整Python代码
4.结果截图

1.遗传算法概述

对于遗传算法：遗传算法是一种用于解决最优化的搜索算法，也是进化算法的一种。取名遗传就是因为它借鉴了生物学中的一些概念，比如说遗传、变异、自然选择以及杂交等等。

对于个体：遗传算法中的个体可以抽象为染色体，然后使得种群向更好的方向进化。一般来说，染色体可以用一系列0和1的串来表示，通常使用的也就是二进制表示法。

对于进化：首先生成随机个体构成的种群，然后每个个体都有自己的适应度，基于适应度随机选择多个个体，然后通过自然选择以及突变产生新的生命种群，进行下一次的迭代。

GA参数(接下来代码中要使用的)：

参数名	含义
popsize	种群规模，我们使用100
binarylength	构成个体染色体的二进制串长度，我们使用10
pc	发生交叉的概率，我们使用0.6，60%的概率交叉
pm	发生变异的概率，我们使用0.001，0.1%的概率变异
pop	初始化种群
newpop	新一代种群

特点：

∙\bullet∙ 当适应度函数选择不当时，可能会收敛于局部最优。在这篇文章中我们的适应度函数就是计算的函数，函数值越大，适应性越好。

∙\bullet∙ 初始化种群的数量十分重要，初始种群多大，算法会占用大量资源；初始种群小时，可能不能求出最优解。

∙\bullet∙ 对于每个解，我们都要对其进行编码，这样才方便写下面的交叉、变异函数。

∙\bullet∙ 变异率也是一个重要的参数。

∙\bullet∙ 选择过程很重要，下面再具体的代码中讲解我们的选择思路。

∙\bullet∙ 遗传算法很快就能找到良好的解，即使是在很复杂的解空间中。

2.分步实现过程

任务:使用遗传算法求得函数y=10×sin(5x)+7∣x−5∣+10y = 10×sin(5x)+7\left |x-5 \right |+10y=10×sin(5x)+7∣x−5∣+10的最大值,定义域为[0,10]。函数图像如下所示：

∙\bullet∙ GA参数的设定（放在主函数中）：

popsize = 100 #种群规模
binarylength = 10 #二进制编码长度（DNA）
pc = 0.6 #交叉概率
pm = 0.001  #变异概率
pop = initpop(popsize,binarylength) #初始化种群

∙\bullet∙ 初始化种群函数initpop：两个参数，popsize值种群规模，binarylength值二进制编码长度。注意：我们认定最左边为最高位！

def initpop(popsize,binarylength):#生成popsize×binarylength的二维0、1序列pop = np.random.randint(0,2,(popsize,binarylength)) return pop

这里我们使用numpy的函数来生成0、1序列，给大家看一看结果：

∙\bullet∙ 将二进制编码转换为十进制数（指定定义域内）的函数bintodec：参数ypop为种群参数，因为计算过程中要修改二级制种群中的数值，所以我们拷贝了种群函数进行处理。

def bintodec(ypop):pop=ypop.copy() #拷贝种群[row,col] = pop.shape #得到pop种群的行数、列数#将二进制数组中每一个位置的值转换位对应的十进制数for i in range(col):pop[:,i]=2**(col-1-i)*pop[:,i]#每一行求和pop = np.sum(pop,axis=1)num=[]#因为二进制串为10位，所以我们除以1023将其限制再0-1之间，然后再乘以10得到定义域内的随机数。num=pop*10/1023return num

∙\bullet∙ 计算种群适应度函数cal_objval：适应度为其对应的函数值，函数值越大，适应度越高。

def cal_objval(pop):x = bintodec(pop)objval = 10*np.sin(5*x)+7*abs(x-5)+10return objval

∙\bullet∙ 选择函数selection，参数pop为种群，参数fitval为种群适应度，参数popsize为种群规模。根据其适应度选择个体，适应度高的个体我们多选择一些，来替换掉适应度低的个体（总个体数仍为popsize）。返回新生成的种群。

def selection(pop,fitval,popsize):idx = np.random.choice(np.arange(popsize),size=popsize,replace=True,p=fitval/fitval.sum())return pop[idx]

∙\bullet∙ 交叉函数crossover：参数pop为种群，参数pc为交叉概率。生成同样规格的newpop数组来存储新种群。

交叉规则：我们将其分成50对染色体，相邻的两个染色体为新染色体的父本和母本。每次随机生成一个数，如果生成数落在交叉概率pc外（大于pc）说明不需要进行交叉操作，那么将父本和母本原封不动的拷贝到两个子代个体中。如果生成数落在交叉概率pc内（小于pc）说明需要进行交叉操作，那么我们再随机生成一个0-9的数cpoint，我们将父本中[0:cpoint]和母本中[cpoint+1:py]组合形成第一个子代个体，将母本[0:cpoint]和父本[cpoint+1:py]组合形成子代第二个个体。对于50对染色体重复该操作，生成新的100个个体。

def crossover(pop,pc):[px,py] = pop.shapenewpop = np.ones((px,py))for i in range(0,px,2):if np.random.rand()<pc:cpoint = int(np.random.rand()*py*10//10)newpop[i,0:cpoint]=pop[i,0:cpoint]newpop[i,cpoint:py]=pop[i+1,cpoint:py]newpop[i+1,0:cpoint]=pop[i+1,0:cpoint]newpop[i+1,cpoint:py]=pop[i,cpoint:py]
#             newpop[i+1,:]=[pop[i+1,0:cpoint],pop[i,cpoint:py]]else:newpop[i,:]=pop[i,:]newpop[i+1,:]=pop[i+1,:]return newpop

∙\bullet∙ 变异函数mutation：参数pop为种群，参数pm为变异概率。

变异规则：对于每一个个体，我们首先生成一个随机数，如果随机数落在变异概率pm外（大于pm），那么就不发生变异。如果随机数落在变异概率pm内（小于pm）说明需要进行变异，我们再随机生成一个0-9的数mpoint，将该个体中mpoint位置的0、1互换，发生变异。

def mutation(pop,pm):[px,py] = pop.shapenewpop = np.ones((px,py))for i in range(px):if(np.random.rand()<pm):mpoint = int(np.random.rand()*py*10//10)if mpoint<=0:mpoint=1newpop[i,:]=pop[i,:]if newpop[i,mpoint]==0:newpop[i,mpoint]=1else:newpop[i,mpoint]=0else:newpop[i,:]=pop[i,:]return newpop

∙\bullet∙ 求得最优解函数best：参数pop为种群，参数fitvalue为适应度。选择种群中适应度最高的个体返回其二进制序列bestindividual以及适应度bestfit。

def best(pop,fitvalue):[px,py]=pop.shapebestindividual = pop[0,:]bestfit = fitvalue[0]for i in range(1,px):if fitvalue[i]>bestfit:bestindividual = pop[i,:]bestfit = fitvalue[i]return bestindividual,bestfit

∙\bullet∙ 主函数（包含绘图函数）：

if __name__=="__main__":popsize = 100 #种群规模binarylength = 10 #二进制编码长度（DNA）pc = 0.6 #交叉概率pm = 0.001  #变异概率pop = initpop(popsize,binarylength) #初始化种群#进行计算,迭代一百次，每十次画一张图for i in range(100):#计算当前种群适应度objval = cal_objval(pop)fitval = objval#选择操作newpop = selection(pop,fitval,popsize)#交叉操作newpop = crossover(newpop,pc);#变异操作newpop = mutation(newpop,pm);#更新种群pop = newpop;#寻找最优解并绘图[bestindividual,bestfit]=best(pop,fitval)x1 = bintodec(newpop)y1 = cal_objval(newpop)x = np.arange(0,10,0.1)y = 10*np.sin(5*x)+7*abs(x-5)+10if i%10==0:plt.figure()plt.rcParams["font.sans-serif"] = ["SimHei"]   #解决中文乱码问题plt.rcParams["axes.unicode_minus"] = False   #使一些符号正常显示plt.plot(x,y)plt.plot(x1,y1,'*')plt.title('迭代次数为:%d'%i)[n]=bestindividual.shapex2=0for i in range(n):x2+=2**(n-1-i)*bestindividual[i]print("The best X is ",x2*10/1023)print("The best Y is ",bestfit)

3.完整Python代码

import numpy as np
import matplotlib.pyplot as pltdef initpop(popsize,binarylength):pop = np.random.randint(0,2,(popsize,binarylength)) #生成popsize×binarylength的二维0、1序列return popdef bintodec(ypop):pop=ypop.copy()[row,col] = pop.shapefor i in range(col):pop[:,i]=2**(col-1-i)*pop[:,i]pop = np.sum(pop,axis=1)num=[]num=pop*10/1023return num#计算种群适应度
def cal_objval(pop):x = bintodec(pop)objval = 10*np.sin(5*x)+7*abs(x-5)+10return objvaldef selection(pop,fitval,popsize):idx = np.random.choice(np.arange(popsize),size=popsize,replace=True,p=fitval/fitval.sum())return pop[idx]def crossover(pop,pc):[px,py] = pop.shapenewpop = np.ones((px,py))for i in range(0,px,2):if np.random.rand()<pc:cpoint = int(np.random.rand()*py*10//10)newpop[i,0:cpoint]=pop[i,0:cpoint]newpop[i,cpoint:py]=pop[i+1,cpoint:py]newpop[i+1,0:cpoint]=pop[i+1,0:cpoint]newpop[i+1,cpoint:py]=pop[i,cpoint:py]
#             newpop[i+1,:]=[pop[i+1,0:cpoint],pop[i,cpoint:py]]else:newpop[i,:]=pop[i,:]newpop[i+1,:]=pop[i+1,:]return newpopdef mutation(pop,pm):[px,py] = pop.shapenewpop = np.ones((px,py))for i in range(px):if(np.random.rand()<pm):mpoint = int(np.random.rand()*py*10//10)if mpoint<=0:mpoint=1newpop[i,:]=pop[i,:]if newpop[i,mpoint]==0:newpop[i,mpoint]=1else:newpop[i,mpoint]=0else:newpop[i,:]=pop[i,:]return newpopdef best(pop,fitvalue):[px,py]=pop.shapebestindividual = pop[0,:]bestfit = fitvalue[0]for i in range(1,px):if fitvalue[i]>bestfit:bestindividual = pop[i,:]bestfit = fitvalue[i]return bestindividual,bestfitif __name__=="__main__":popsize = 100 #种群规模binarylength = 10 #二进制编码长度（DNA）pc = 0.6 #交叉概率pm = 0.001  #变异概率pop = initpop(popsize,binarylength) #初始化种群#进行计算for i in range(100):#计算当前种群适应度objval = cal_objval(pop)fitval = objval#选择操作newpop = selection(pop,fitval,popsize)#交叉操作newpop = crossover(newpop,pc);#变异操作newpop = mutation(newpop,pm);#更新种群pop = newpop;#寻找最优解并绘图[bestindividual,bestfit]=best(pop,fitval)x1 = bintodec(newpop)y1 = cal_objval(newpop)x = np.arange(0,10,0.1)y = 10*np.sin(5*x)+7*abs(x-5)+10if i%10==0:plt.figure()plt.rcParams["font.sans-serif"] = ["SimHei"]   #解决中文乱码问题plt.rcParams["axes.unicode_minus"] = False   #使一些符号正常显示plt.plot(x,y)plt.plot(x1,y1,'*')plt.title('迭代次数为:%d'%i)[n]=bestindividual.shapex2=0for i in range(n):x2+=2**(n-1-i)*bestindividual[i]print("The best X is ",x2*10/1023)print("The best Y is ",bestfit)

4.结果截图

数组结果为100次迭代后的最值横纵坐标：（X,Y）=（0.18572825024437928，52.8982183292736）。因为每次初始化数组会有不同，所以每次的结果都可能会不一样。这里选取四次迭代截图进行展示。