遗传算法(Genetic Algorithm,GA)起源于对生物系统所进行的计算机模拟研究,用于解决寻找最优解的问题。
它是模仿自然界生物进化机制发展起来的随机全局搜索和优化方法,借鉴了达尔文的进化论以及孟德尔的遗传学说。
遗传算法的一般步骤为:
1.随机产生种群
2.确定个体适应度函数,判断个体适应度,是否符合优化准则,若符合,输出最佳个体及最优解(即最优基因型),否则,进行下一步
3.依据适应度选择父母,适应度高的个体被选中的概率高,适应度低的个体被淘汰
4.用父母的染色体按照一定方法进行交叉,生成子代
5.子代染色体按照一定方式变异
6.由交叉变异后的形成的新一代种群,返回步骤2,直至最优解产生
下面通过解决实际问题来进一步深化对于遗传算法的认识。
问题1:求函数y=10sin5x+7cos4x,x∈[0,10]的最大值
此函数在区间上的图像如图:
下面通过遗传算法来近似找出最大值。
首先是初始化种群。我们构建一个种群个体数为200,DNA数为20的种群,并且定义一些参数,代码如下:
DNA_SIZE=20
POP_SIZE=200
CROSSOVER_RATE=0.8
MUTATION_RATE=0.005
N_GENERATIONS=100
X_BOUND=[0,10]
下面是各个过程的解释以及代码。
解码过程。之所以有这个过程是因为我们编码使用的二进制会产生一个位数很大的数据,这让我们处理起来不大方便,我们希望找到一个一一映射,将二进制数与我们规定区间内的一些实数一一对应起来,下面的代码可以帮助我们实现这个过程:
def translateDNA(pop):
x=pop.dot(2**np.arange(DNA_SIZE)[::-1])/float(2**DNA_SIZE-1)*(X_BOUND[1]-X_BOUND[0])+X_BOUND[0]
return x
判断适应度:
def f(x):
return 10*np.sin(5*x)+7*np.cos(4*x)
def get_fitness(pop):
x=translateDNA(pop)
pred=f(x)
return (pred-np.min(pred))+1e-3
以上函数是为了让所有适应度均为大于零的正数,方便计算。
让种群产生交叉变异:
def crossover_and_mutation(pop,CROSSOVER_RATE=0.8):
new_pop=[]
for father in pop:
child = father
if np.random.rand() < CROSSOVER_RATE:
mother = pop[np.random.randint(POP_SIZE)]
cross_points=np.random.randint(low=0,high=DNA_SIZE)
child[cross_points:]=mother[cross_points:]
mutation(child)
new_pop.append(child)
return new_pop
def mutation(child,MUTATION_RATE=0.005):
if np.random.rand() < MUTATION_RATE:
mutate_point=np.random.randint(0,DNA_SIZE)
child[mutate_point]=child[mutate_point]^1
自然选择,自然选择过程我们应该以一定概率选择各个适应度的个体,numpy.random.choice函数可以帮助我们实现这一功能:
def select(pop,fitness):
idx=np.random.choice(np.arange(POP_SIZE),size=POP_SIZE,replace=True,p=(fitness)/float(fitness.sum()))
return pop[idx]
接下来就是将以上过程进行迭代,迭代一定次数后我们便可选择出相对较优的个体。
完整代码如下:
import numpy as np
DNA_SIZE=20
POP_SIZE=200
CROSSOVER_RATE=0.8
MUTATION_RATE=0.005
N_GENERATIONS=100
X_BOUND=[0,10]
def f(x):
return 10*np.sin(5*x)+7*np.cos(4*x)
def translateDNA(pop):
x=pop.dot(2**np.arange(DNA_SIZE)[::-1])/float(2**DNA_SIZE-1)*(X_BOUND[1]-X_BOUND[0])+X_BOUND[0]
return x
def get_fitness(pop):
x=translateDNA(pop)
pred=f(x)
return (pred-np.min(pred))+1e-3
def crossover_and_mutation(pop,CROSSOVER_RATE=0.8):
new_pop=[]
for father in pop:
child = father
if np.random.rand() < CROSSOVER_RATE:
mother = pop[np.random.randint(POP_SIZE)]
cross_points=np.random.randint(low=0,high=DNA_SIZE)
child[cross_points:]=mother[cross_points:]
mutation(child)
new_pop.append(child)
return new_pop
def mutation(child,MUTATION_RATE=0.005):
if np.random.rand() < MUTATION_RATE:
mutate_point=np.random.randint(0,DNA_SIZE)
child[mutate_point]=child[mutate_point]^1
def select(pop,fitness):
idx=np.random.choice(np.arange(POP_SIZE),size=POP_SIZE,replace=True,p=(fitness)/float(fitness.sum()))
return pop[idx]
def print_info(pop):
fitness=get_fitness(pop)
max_fitness_index=np.argmax(fitness)
print("max_fitness:",fitness[max_fitness_index])
x=translateDNA(pop)
print("最优基因型:",pop[max_fitness_index])
print("x:",x[max_fitness_index])
def main():
pop = np.random.randint(2,size=(POP_SIZE,DNA_SIZE))
for i in range(N_GENERATIONS):
x=translateDNA(pop)
pop=np.array(crossover_and_mutation(pop,CROSSOVER_RATE))
fitness=get_fitness(pop)
pop=select(pop,fitness)
print_info(pop)
main()
最终的运行结果为:
我们通过之前的图像可以看到 ,该段函数在x=1.5与7.8左右有两个峰值,可见我们的程序跑出了一个相对令人满意的结果。如果增加DNA数目、种群个体数以及迭代次数,结果会更加精确。
问题2:在x,y∈【-3,3】内求如下二元函数的最大值:
z=sin(x**2+y)+cos(x+y**2)
完整代码以及注释如下:
#应用遗传算法计算二元函数最大值
import numpy as np
from matplotlib import pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
from matplotlib import cm
DNA_SIZE=24
POP_SIZE=200#种群个数
CROSSOVER_RATE=0.8#交叉概率
MUTATION_RATE=0.005#变异概率
N_GENERATIONS=100#迭代代数
X_BOUND=[-3,3]
Y_BOUND=[-3,3]
#定义二元函数:
def F(x,y):
return np.sin(x**2+y)+np.cos(x+y**2)
def plot_3d(ax):
X = np.linspace(*X_BOUND, 100)
Y = np.linspace(*Y_BOUND, 100)
X,Y = np.meshgrid(X, Y)
Z = F(X, Y)
ax.plot_surface(X,Y,Z,rstride=1,cstride=1,cmap=cm.coolwarm)
ax.set_zlim(-10,10)
ax.set_xlabel('x')
ax.set_ylabel('y')
ax.set_zlabel('z')
plt.pause(3)
plt.show()
#解码过程:
def translateDNA(pop):#pop表示种群矩阵,一行表示一个二进制编码的DNA,行数为种群数目
x_pop=pop[:,1::2]#奇数列表示X
y_pop=pop[:,::2] #偶数列表示Y
#pop:(POP_SIZE,DNA_SIZE)*(DNA_SIZE,1) --> (POP_SIZE,1)
x = x_pop.dot(2**np.arange(DNA_SIZE)[::-1])/float(2**DNA_SIZE-1)*(X_BOUND[1]-X_BOUND[0])+X_BOUND[0]
y = y_pop.dot(2**np.arange(DNA_SIZE)[::-1])/float(2**DNA_SIZE-1)*(Y_BOUND[1]-Y_BOUND[0])+Y_BOUND[0]
return x,y
#求适应度:
def get_fitness(pop):
x,y=translateDNA(pop)
pred=F(x,y)
return (pred-np.min(pred))+1e-3
#交叉变异过程:
def crossover_and_mutation(pop,CROSSOVER_RATE=0.8):
new_pop=[]
for father in pop:
child = father
if np.random.rand() < CROSSOVER_RATE:
mother = pop[np.random.randint(POP_SIZE)]
cross_points=np.random.randint(low=0,high=DNA_SIZE*2)
child[cross_points:] = mother[cross_points:]
#变异:
mutation(child)
new_pop.append(child)
return new_pop
#变异过程:
def mutation(child,MUTATION_RATE=0.005):
if np.random.rand() < MUTATION_RATE:
mutate_point=np.random.randint(0,DNA_SIZE*2)
child[mutate_point]=child[mutate_point]^1
#自然选择过程:
def select(pop,fitness):
idx=np.random.choice(np.arange(POP_SIZE),size=POP_SIZE,replace=True,p=(fitness)/(fitness.sum()))
return pop[idx]
#输出:
def print_info(pop):
fitness=get_fitness(pop)
max_fitness_index=np.argmax(fitness)
print("max_fitness:",fitness[max_fitness_index])
x,y=translateDNA(pop)
print("最优基因型:",pop[max_fitness_index])
print("(x,y):",(x[max_fitness_index],y[max_fitness_index]))
def main():
fig = plt.figure()
ax = Axes3D(fig)
plt.ion()#将画图模式改为交互模式,程序遇到plt.show不会暂停,而是继续执行
plot_3d(ax)
pop = np.random.randint(2, size=(POP_SIZE, DNA_SIZE*2)) #matrix (POP_SIZE, DNA_SIZE)
for _ in range(N_GENERATIONS):#迭代N代
x,y = translateDNA(pop)
if 'sca' in locals():
sca.remove()
sca = ax.scatter(x, y, F(x,y), c='black', marker='o');plt.show();plt.pause(0.1)
pop = np.array(crossover_and_mutation(pop, CROSSOVER_RATE))
#F_values = F(translateDNA(pop)[0], translateDNA(pop)[1])#x, y --> Z matrix
fitness = get_fitness(pop)
pop = select(pop, fitness) #选择生成新的种群
print_info(pop)
plt.ioff()
plot_3d(ax)
main()
我们这次加上了显示图像的代码,用于更加明显的看到结果:
从图像中我们看到遗传算法给出的(x,y)数对值在图像中对应中部高峰,基本解决了求最大值的问题。
以上多数内容从下面这篇文章学习,感谢博主。