基于DEAP库的python进化算法-3.简单遗传算法的实现

文章目录

前言
一.问题描述与分析
二.问题的编码与解码
三.利用DEAP自带算法求解
- 1.DEAP自带的进化算法介绍
- 2.调用DEAP自带的进化算法
- 3.自行编写算法求解

前言

在上一篇中，我们已经介绍了如何在DEAP中实现进化算法的基本操作，这一篇中我们试图将各个操作组装起来，用进化算法解决一个简单的一元函数寻优问题

一.问题描述与分析

给定一个函数，求解其最大值

该函数的图像为：

该函数的最大值应该出现在28.309042309042308处，值为1657.4235763265594。

可以看到该函数有很多局部极值作为干扰项，如果进化算法过早收敛，很容易陷入某个局部最优。

二.问题的编码与解码

对于该问题，可以选取很多不同的编码方式。本文计划采用二进制编码，精确度要求到6位，那么首先应当确定编码的长度与解码方式。由于有：

所以我们需要的二进制编码长度为26

三.利用DEAP自带算法求解

1.DEAP自带的进化算法介绍

DEAP自带的算法都比较基础，通常可以用来测试问题描述、编码方式和交叉突变操作组合的有效性。需要比较复杂的进化算法时，可以通过在已有的算子上进行扩充

简单进化算法

deap.algorithms.eaSimple

DEAP中预置的简单进化算法流程描述如下：

1.根据工具箱中注册的toolbox.evaluate评价族群
2.根据工具箱中注册的toolbox.select选择与父代相同个数的育种个体
3.在族群中进行第一次循环，用工具箱中注册的toolbox.mate进行配种，并用生成的两个子代替换对应父代
4.在族群中进行第二次循环，用工具箱中注册的toolbox.mutate进行变异，用变异后的子代替换对应父代
5.从1开始重复循环，直到达到设定的迭代次数

需要注意的是在这个过程中，生成子代有四种情况：受到配种影响；受到变异影响；既受到配种也受到变异影响；既不受配种影响也不受变异影响
对应的伪代码可以表述为：

evaluate(population)
for g in range(ngen):population = select(population, len(population))offspring = varAnd(population, toolbox, cxpb, mutpb)evaluate(offspring)population = offspring

(μ + λ)进化算法

deap.algorithms.eaMuPlusLambda

该算法的流程如下：

1.根据工具箱中注册的toolbox.evaluate评价族群
2.在族群中进行循环，在每次循环中，随机选择crossover，mutation和reproduction三者之一：
如果选择到crossover，那么随机选择2个个体，用工具箱中注册的toolbox.mate进行配种，<将生成的第一个子代加入到后代列表中，第二个子代丢弃；
如果选择到mutation，用工具箱中注册的toolbox.mutate进行变异，将变异后的子代加入到后代列表中；
如果选择到reproduction，随机选择一个个体，将其复制加入到后代列表中
3.根据工具箱中注册的toolbox.select，在父代+子代中选择给定数量的个体作为子代
4.从1开始重复循环，直到达到设定的迭代次数

注意在这个子代生成的过程中，子代不会同时受到变异和配种影响。

对应的伪代码可以表述为：

evaluate(population)
for g in range(ngen):offspring = varOr(population, toolbox, lambda_, cxpb, mutpb)evaluate(offspring)population = select(population + offspring, mu)

(μ，λ)进化算法

deap.algorithms.eaMuCommaLambda

与(μ + λ)进化算法基本相同，唯一的区别在于生成子代族群时，只在产生的子代中选择，而丢弃所有父代。
对应伪代码可以表述为

evaluate(population)
for g in range(ngen):offspring = varOr(population, toolbox, lambda_, cxpb, mutpb)evaluate(offspring)population = select(offspring, mu)

2.调用DEAP自带的进化算法

在调用DEAP自带的算法时需要注意的是，由于内置算法调用的alias已经提前给定，因此我们在register的时候，需要按照给定名称注册。

例如toolbox.register('crossover', tools.cxUniform)就不能被内置算法识别，而应当按照
要求，命名为mate，并且显示给出交叉概率：
toolbox.register('mate', tools.cxUniform, indpb = 0.5)按照要求，使用预置的算法需要注册的工具有：
toolbox.evaluate：评价函数toolbox.select：育种选择toolbox.mate：交叉操作toolbox.mutate：突变操作

from deap import algorithms, creator, tools, base
from scipy.stats import bernoulli
import numpy as np
import randomrandom.seed(42)  # 保证结果可以复现
# 定义问题
creator.create('FitnessMax', base.Fitness, weights=(1.0,))  # 单目标优化，最大值问题
creator.create('Individual', list, fitness=creator.FitnessMax)# 生成个体
gene_size = 26  # 26位编码
toolbox = base.Toolbox()
toolbox.register('Binary', bernoulli.rvs, 0.5)
toolbox.register('Individual', tools.initRepeat, creator.Individual, toolbox.Binary, n=gene_size)# 解码 - 二进制转换为十进制
def decode(individual):# 解码到10进制num = int(''.join([str(_) for _ in individual]), 2)# 映射到-30到30区间x = -30 + num * 60 / (2 ** 26 - 1)return x# 评价函数-适应度
def eval(individual):# 转化到区间内实数x = decode(individual)return ((np.square(x) + x) * np.cos(2 * x) + np.square(x) + x),# 生成初始族群
pop_size = 100  # 族群中的个体数量
toolbox.register('Population', tools.initRepeat, list, toolbox.Individual)
pop = toolbox.Population(n=pop_size)# 在工具箱中注册遗传算法需要的工具
toolbox.register('evaluate', eval)
toolbox.register('select', tools.selTournament, tournsize=2)  # tournsize=2的锦标赛选择,参数少k
toolbox.register('mate', tools.cxUniform, indpb=0.5)  # 均匀交叉，参数少ind1,ind2
toolbox.register('mutate', tools.mutFlipBit, indpb=0.5)  # 位翻转变异# 注册计算过程中需要记录的数据
stats = tools.Statistics(key=lambda ind: ind.fitness.values)
stats.register('avg', np.mean)
stats.register('std', np.std)
stats.register('min', np.min)
stats.register('max', np.max)# 调用DEAP内置算法
resultPop, logbook = algorithms.eaSimple(pop, toolbox, cxpb=0.5, mutpb=0.2, ngen=50, stats=stats, verbose=False)# 输出计算过程
# nevals代表该迭代中调用evalute函数的次数
logbook.header = 'gen', 'nevals', 'avg', 'std', 'min', 'max'
print(logbook)

计算过程输出

gen  nevals  avg     std     min         max
0   100     279.668 376.539 -0.403678   1517.68
1   53      375.377 388.779 0.0422919   1641.3
2   59      335.908 371.105 0.0446177   1516.72
3   52      353.956 324.161 -0.280488   1505.81
4   59      445.318 406.245 0.841933    1528.44
5   62      519.898 415.008 1.39083     1606.41
6   71      499.191 430.979 0.184617    1417.72
7   62      511.865 399.538 0.0560637   1586.01
8   67      491.792 421.77  0.224034    1578.36
9   55      482.469 398.13  1.03128     1426.2
10  59      565.195 445.841 0.0425501   1427.1
11  50      655.547 444.508 0.127139    1547.92
12  49      670.214 466.569 0.0357719   1541.85
13  60      616.596 498.718 0.0158924   1544.06
14  56      573.002 484.627 -0.36737    1544.06
15  61      588.575 473.83  0.00965977  1533.09
16  57      671.817 476.531 0.0386515   1533.09
17  72      692.467 557.545 0.0437906   1639.39
18  66      700.388 527.192 0.0205807   1485.58
19  67      750.856 526.877 0.00740906  1542.36
20  64      719.56  535.471 -0.0687704  1387.21
21  58      691.961 546.773 0.00459343  1469.56
22  60      635.223 525.815 0.0344323   1602.17
23  49      759.276 519.709 0.0329074   1602.31
24  47      928.791 519.078 -0.395298   1628.81
25  57      948.308 530.236 0.0115315   1628.42
26  53      1032.62 496.975 0.0556756   1653.9
27  71      920.054 552.981 -0.24814    1656.97
28  52      1000.62 527.772 0.975663    1656.97
29  61      993.985 556.678 0.0756073   1653.9
30  56      988.46  553.916 2.33859     1653.91
31  68      938.748 598.801 -0.343211   1653.91
32  58      968.058 571.579 -0.374508   1653.91
33  61      1052.52 584.106 -0.395098   1657.17
34  60      1018.82 593.899 0.0344873   1657.17
35  50      1094.81 558.391 0.180522    1657.29
36  40      1149.37 577.618 -0.0701393  1657.29
37  51      1192.23 620.757 -0.0860435  1657.29
38  48      1284.57 593.595 0.434814    1657.29
39  55      1314.9  557.095 5.45569     1657.29
40  58      1279.23 632.164 0.00597406  1657.3
41  54      1339.26 578.89  -0.236449   1657.31
42  53      1329.61 591.209 0.0531222   1657.31
43  60      1319.97 607.591 0.834223    1657.31
44  65      1237    645.968 0.000313783 1657.31
45  71      1284.9  626.286 -0.364857   1657.31
46  51      1392.22 517.275 0.16591     1657.31
47  66      1310.41 592.914 0.304726    1657.31
48  66      1161.6  691.608 0.0934323   1657.31
49  50      1277.25 638.279 0.000197178 1657.31
50  61      1288.16 618.87  0.168657    1657.31

输出结果

# 输出最优解
index = np.argmax([ind.fitness for ind in resultPop])
x = decode(resultPop[index])
print('当前最优解:'+ str(x) + '\t对应的函数值为:'+ str(resultPop[index].fitness))
# 当前最优解:28.294433329916494  对应的函数值为:(1657.069165527869,)

可以看出遗传算法成功避开了局部最优，给出的结果非常接近全局最优解28.309

3.自行编写算法求解

自行编写通常来说会需要比内置函数更长的篇幅，但是也能获得更大的自由度。下面是一个用锦标赛交叉、位翻转突变求解同样问题的例子

import random
import numpy as np
from scipy.stats import bernoulli
from deap import creator,base,toolsrandom.seed(42)# 定义问题
creator.create('FitnessMax',base.Fitness,weights=(1.0,))
creator.create('Individual',list,fitness=creator.FitnessMax)# 生成个体-二进制编码
gene_size = 26
toolbox = base.Toolbox()
toolbox.register('Binary',bernoulli.rvs,0.5)
toolbox.register('Individual',tools.initRepeat,creator.Individual,toolbox.Binary,n=gene_size)# 生成初始种群
pop_size = 100      # 族群中个体的数量
toolbox.register('Population',tools.initRepeat,list,toolbox.Individual)
pop = toolbox.Population(n=pop_size)# 评价函数-二进制解码到十进制并计算适应度
def eval(individual):# 二进制解码到十进制num = int(''.join([str(_) for _ in individual]),2)# 转换为区间内的对应实数x = -30 + 60 * num /(2**26 - 1)return ((np.square(x) + x) * np.cos(2*x) + np.square(x) + x),
# 注册遗传算法工具-评价函数
toolbox.register('evaluate',eval)# 评价初始族群
fitnesses = map(toolbox.evaluate,pop)
for ind,fit in zip(pop,fitnesses):ind.fitness.values = fit# 进化迭代
N_GEN = 50      # 迭代次数
CXPB = 0.5      # 交叉概率
MUTPB = 0.2     # 突变概率# 注册进化过程中需要的工具：配种选择，交叉，变异
toolbox.register('select',tools.selTournament,tournsize = 2)       # tournsize = 2的锦标赛选择
toolbox.register('mate',tools.cxUniform)        # 均匀交叉
toolbox.register('mutate',tools.mutFlipBit)     # 位翻转变异# 用数据记录算法迭代过程
# 创建统计信息对象
stats = tools.Statistics(key=lambda ind : ind.fitness.values)
stats.register('avg',np.mean)
stats.register('std',np.std)
stats.register('min',np.min)
stats.register('max',np.max)
logbook = tools.Logbook()for gen in range(N_GEN):# 配种选择selectedTour = toolbox.select(pop,pop_size)     # 选择pop_size个个体# 复制个体，供交叉变异用selectedInd = list(map(toolbox.clone,selectedTour))# 对选出的育种族群两两进行交叉，对于被改变的个体，删除其适应度值# list[start:end:step]for child1,child2 in zip(selectedInd[::2],selectedInd[1::2]):if random.random() < CXPB:toolbox.mate(child1,child2,0.5)     # 均匀交叉del child1.fitness.valuesdel child2.fitness.values# 对选出的育种族群进行变异，对于被改变个体，删除适应度值for mutant in selectedInd:if random.random() < MUTPB:toolbox.mutate(mutant,0.5)      # 位翻转变异del mutant.fitness.values# 对于被改变的个体，重新评价其适应度值invalid_ind = [ind for ind in selectedInd if not ind.fitness.valid]fitnesses = map(toolbox.evaluate,invalid_ind)for ind,fit in zip(invalid_ind,fitnesses):ind.fitness.values = fit# 完全重插入pop[:] = selectedInd# 记录数据# 将每个注册功能应用于输入序列数据，并将结果作为字典返回record = stats.compile(pop)logbook.record(gen=gen,**record)# 输出计算过程
logbook.header = 'gen','avg','std','min','max'
print(logbook)

计算过程输出

gen  avg     std     min         max
0   466.229 468.083 -0.361262   1639.96
1   523.428 500.73  0.000387387 1639.96
2   591.327 499.476 0.761878    1639.96
3   511.309 484.249 -0.404202   1639.96
4   542.566 455.348 -0.367558   1549.18
5   511.607 477.01  0.0601817   1630.98
6   520.779 496.691 0.00812873  1630.98
7   588.768 520.461 -0.298492   1631.01
8   541.561 499.122 -0.202241   1631.01
9   730.499 538.137 0.0242079   1631.01
10  839.125 559.088 0.43618     1631.01
11  829.414 574.283 0.233167    1631.01
12  797.474 641.566 0.0437757   1654.98
13  813.062 648.566 0.999945    1657.27
14  890.102 685.853 -0.154962   1657.28
15  937.362 677.619 0.394182    1657.28
16  1107    656.154 -0.0508697  1657.28
17  1216.31 648.938 1.39068     1657.28
18  1272.64 625.1   -0.270614   1657.28
19  1218.76 612.227 0.0275138   1657.28
20  1082.98 678.12  0.0401458   1657.28
21  1100.7  630.675 0.0618022   1657.28
22  1244.65 595.943 7.8221e-05  1657.28
23  1235.14 641.172 -0.180385   1657.3
24  1219.11 660.604 0.051994    1657.3
25  1400.85 531.75  0.0735181   1657.3
26  1161.86 692.244 -0.22549    1657.3
27  1248.31 642.029 0.150307    1657.3
28  1303.98 619.694 1.40528     1657.3
29  1244.15 652.702 0.0203885   1657.3
30  1199.76 676.044 0.0233112   1657.3
31  1202.59 672.031 -0.136965   1657.3
32  1204.52 667.777 -0.0454565  1657.31
33  1248.02 630.396 0.283428    1657.31
34  1340.9  582.728 14.4446     1657.31
35  1359.91 597.02  -0.378353   1657.31
36  1286.13 640.359 -0.109894   1657.31
37  1310.39 618.029 0.898644    1657.31
38  1308.29 568.212 10.3651     1657.31
39  1283    635.098 0.0838442   1657.31
40  1338.32 584.495 -0.241494   1657.31
41  1334.45 591.224 0.0600735   1657.31
42  1327.41 606.735 1.83592     1657.31
43  1236.79 648.824 0.000148321 1657.31
44  1283.51 627.575 -0.384135   1657.36
45  1393.7  521.303 0.159144    1657.36
46  1310.79 590.335 0.529636    1657.36
47  1161.16 692.163 0.0592846   1657.36
48  1271.36 641.495 0.000130297 1657.36
49  1284.52 625.65  0.213304    1657.36

输出最优解

# 输出最优解
index = np.argmax([ind.fitness for ind in pop])
x = decode(pop[index])
print('当前最优解:'+str(x)+'\t对应的函数值为:'+str(pop[index].fitness))
# 当前最优解:28.308091138423848  对应的函数值为:(1657.4220750846391,)

结果可视化

import matplotlib.pyplot as plt# 用select方法从logbook中提取迭代次数,最大值、均值
gen = logbook.select('gen')
fit_maxs = logbook.select('max')
fit_average = logbook.select('avg')fig = plt.figure()
ax = fig.add_subplot(111)
ax.plot(gen,fit_maxs,'b-',linewidth=2.0,label='Max Fitness')
ax.plot(gen,fit_average,'r-',linewidth=2.0,label='Average Fitness')
ax.legend(loc='best')
ax.set_xlabel('Generation')
ax.set_ylabel('Fitness')plt.tight_layout()
plt.savefig('Generation_Fitness.png')