遗传算法实践(一) 非线性函数寻优问题
2019-07-03 本文已影响1人
ChaoesLuol
前言
本系列总结一些遗传算法的应用案例，主要基于《网络模型与多目标遗传算法》（玄光男 林林 著）一书提到的一些应用，又根据个人兴趣对一些问题进行了修改。这里的遗传算法实现基于python的DEAP库，如果对该库不熟悉的朋友可以参见我之前的介绍文章。作者水平有限，错误难免，请方家不吝指正。
非线性函数寻优问题
为了验证遗传算法的寻优能力，这里用常用非线性函数Beale's Function进行测试。目标函数为：
该函数在搜索域内的最小值已知为。
函数在搜索域内的图像如下：
Beale's function
用遗传算法求解优化问题的一个关键在于如何对问题的搜索空间进行编码： 合理的编码方式能够大幅加快收敛速度，甚至编码的合理与否决定了算法能否寻找到最优解。
对于该问题，有两种编码思路：
- 最直接的就是对两个变量分别用实数进行编码, 这样也省去了解码的繁琐, 所见即所得；
- 另一种思路是用二进制编码对搜索空间进行离散化, 在给定需要的精度的情况下, 可以用有限位的二进制编码描述搜索空间。
接下来我们分别用这两种编码思路对问题进行求解。
实数编码
最简单直接的想法就是用实数对变量进行编码，一个变量对应一个实数。这里采用的遗传算法操作如下：
- 个体编码：使用实数编码
- 评价函数：直接使用原函数值进行评价
- 育种选择：binary锦标赛选择
- 变异算法：交叉-模拟二值交叉(SBX)，突变-有界多项式突变
- 环境选择：采用精英保存策略, 加速收敛
完整代码如下:
# ----------------------
# 问题定义
creator.create('FitnessMin', base.Fitness, weights=(-1.0,)) # 最小化问题
creator.create('Individual', list, fitness=creator.FitnessMin)
# 个体编码
low = [-4.5, -4.5] # 两个变量值的下界
up = [4.5, 4.5] # 两个变量值的上界
# 在下界和上界间用均匀分布生成实数变量
def genInd(low, up):
return [random.uniform(low[0], up[0]), random.uniform(low[1], up[1])]
toolbox = base.Toolbox()
toolbox.register('genInd', genInd, low, up)
toolbox.register('individual', tools.initIterate, creator.Individual, toolbox.genInd)
toolbox.register('population', tools.initRepeat, list, toolbox.individual)
# 评价函数
def evaluate(ind):
f = (1.5-ind[0]+ind[0]*ind[1])*(1.5-ind[0]+ind[0]*ind[1]) + \
(2.25-ind[0]+ind[0]*ind[1]*ind[1])*(2.25-ind[0]+ind[0]*ind[1]*ind[1]) + \
(2.625-ind[0]+ind[0]*ind[1]*ind[1]*ind[1])*(2.625-ind[0]+ind[0]*ind[1]*ind[1]*ind[1])
return f,
toolbox.register('eval', evaluate)
# 记录迭代数据
stats = tools.Statistics(key=lambda ind:ind.fitness.values)
stats.register('avg', np.mean)
stats.register('min', np.min)
stats.register('std', np.std)
logbook = tools.Logbook()
logbook.header = ['gen', 'nevals'] + (stats.fields)
# 注册遗传算法操作 - 选择,交叉,突变
toolbox.register('select', tools.selTournament, tournsize=2)
toolbox.register('mate', tools.cxSimulatedBinaryBounded, eta=20, low=low, up=up)
toolbox.register('mutate', tools.mutPolynomialBounded, eta=20, low=low, up=up, indpb=0.2)
# ----------------------
# 遗传算法主要部分
popSize = 100 # 族群规模
ngen = 200 # 迭代代数
cxpb = 0.8 # 交叉概率
mutpb = 0.1 # 突变概率
pop = toolbox.population(popSize) # 生成初始族群
# 评价初始族群
invalid_ind = [ind for ind in pop if not ind.fitness.valid]
fitnesses = toolbox.map(toolbox.eval, invalid_ind)
for ind, fitness in zip(invalid_ind, fitnesses):
ind.fitness.values = fitness
record = stats.compile(pop)
logbook.record(gen=0, nevals=len(invalid_ind), **record)
# 遗传算法迭代
for gen in range(1, ngen+1):
# 育种选择
offspring = toolbox.select(pop, popSize) # 子代规模与父代相同
offspring = [toolbox.clone(_) for _ in offspring]
# 变异操作
# 突变
for ind1, ind2 in zip(offspring[::2], offspring[1::2]):
if random.random() < cxpb:
toolbox.mate(ind1, ind2)
del ind1.fitness.values
del ind2.fitness.values
for ind in offspring:
if random.random() < mutpb:
toolbox.mutate(ind)
del ind.fitness.values
# 评价子代
invalid_ind = [ind for ind in offspring if not ind.fitness.valid]
fitnesses = toolbox.map(toolbox.eval, invalid_ind)
for ind, fitness in zip(invalid_ind, fitnesses):
ind.fitness.values = fitness
# 环境选择
combinedPop = pop + offspring # 采用精英策略,加速收敛
pop = tools.selBest(combinedPop, popSize)
# 记录数据
record = stats.compile(pop)
logbook.record(gen=gen, nevals=len(invalid_ind), **record)
print(logbook)
结果:
# 输出结果
bestInd = tools.selBest(pop,1)[0]
bestFit = bestInd.fitness.values[0]
print('最优解为: '+str(bestInd))
print('函数最小值为: '+str(bestFit))
## 结果为
## 最优解为: [3.0283366488777057, 0.5077404970344487]
## 函数最小值为: 0.00013976406990386456
对迭代过程可视化:
# 可视化迭代过程
minFit = logbook.select('min')
avgFit = logbook.select('avg')
plt.plot(maxFit, 'b-', label='Minimum Fitness')
plt.plot(avgFit, 'r-', label='Average Fitness')
plt.xlabel('# Gen')
plt.ylabel('Fitness')
plt.legend(loc='best')
# 保存计算结果
plt.tight_layout()
plt.savefig('Real number encoding.png')
Real number encoding
二进制编码
二进制编码方式
在给定了需求的精度，例如小数点后6位时，可以用一定长度的01编码来代表变量。在需求精度为小数点后6位时，我们需要将每单位区间分割为个数字，需要的编码长度为：
对于我们的问题来说，搜索域的下界为-4.5，上界为4.5。因此用24位编码代表一个变量就可以实现我们的精度要求()，一个完整的染色体包含两个变量，长度为48。
解码方式
在评价个体的时候，我们需要将二进制编码转换为落在搜索区间范围内的实数，才能根据函数公式求值。这里的二进制-十进制转换可以根据下式完成：
其中即二进制基因片段对应的十进制数，
为基因片段长度。
这里采用的遗传算法操作如下：
- 个体编码：二进制编码
- 评价函数：解码后根据函数求值
- 育种选择：binary锦标赛选择
- 变异算法：交叉-两点交叉，突变-位翻转突变
- 环境选择：采用精英保存策略, 加速收敛
完整代码如下:
# ----------------------
# 问题定义
creator.create('FitnessMin', base.Fitness, weights=(-1.0,)) # 最小化问题
creator.create('Individual', list, fitness=creator.FitnessMin)
# 个体编码 -- 二进制编码
geneLen = 48 # 基因长度
toolbox = base.Toolbox()
toolbox.register('Binary', random.randint, 0, 1) # 二进制编码
toolbox.register('individual', tools.initRepeat, creator.Individual, toolbox.Binary, geneLen)
toolbox.register('population', tools.initRepeat, list, toolbox.individual)
# 二进制编码的解码
low = [-4.5, -4.5] # 两个变量值的下界
up = [4.5, 4.5] # 两个变量值的上界
def decode(ind, low=low, up=up, geneLen=geneLen):
'''给定一条二进制编码,将其分割为两个变量,并且转换为对应的十进制数'''
# 将一条染色体分为两个变量,转换为10进制
subGeneLen = int(geneLen/2)
x1 = int(''.join(str(_) for _ in ind[:subGeneLen]), 2)
x2 = int(''.join(str(_) for _ in ind[subGeneLen:]), 2)
x1Rescaled = low[0] + x1 * (up[0] - low[0]) / (2**subGeneLen -1)
x2Rescaled = low[1] + x2 * (up[1] - low[1]) / (2**subGeneLen -1)
return x1Rescaled, x2Rescaled
# 评价函数
def evaluate(ind):
x1, x2 = decode(ind) # 先解码为10进制
f = (1.5-x1+x1*x2)*(1.5-x1+x1*x2) + \
(2.25-x1+x1*x2*x2)*(2.25-x1+x1*x2*x2) + \
(2.625-x1+x1*x2*x2*x2)*(2.625-x1+x1*x2*x2*x2)
return f,
toolbox.register('eval', evaluate)
# 迭代数据
stats = tools.Statistics(key=lambda ind:ind.fitness.values)
stats.register('avg', np.mean)
stats.register('min', np.min)
stats.register('std', np.std)
logbook = tools.Logbook()
logbook.header = ['gen', 'nevals'] + (stats.fields)
# 注册遗传算法操作 - 选择,交叉,突变
toolbox.register('select', tools.selTournament, tournsize=2)
toolbox.register('mate', tools.cxTwoPoint)
toolbox.register('mutate', tools.mutFlipBit, indpb=0.5)
# ----------------------
# 遗传算法主要部分
popSize = 100 # 族群规模
ngen = 200 # 迭代代数
cxpb = 0.8 # 交叉概率
mutpb = 0.1 # 突变概率
pop = toolbox.population(popSize) # 生成初始族群
# 评价初始族群
invalid_ind = [ind for ind in pop if not ind.fitness.valid]
fitnesses = toolbox.map(toolbox.eval, invalid_ind)
for ind, fitness in zip(invalid_ind, fitnesses):
ind.fitness.values = fitness
record = stats.compile(pop)
logbook.record(gen=0, nevals=len(invalid_ind), **record)
# 遗传算法迭代
for gen in range(1, ngen+1):
# 育种选择
offspring = toolbox.select(pop, popSize) # 子代规模与父代相同
offspring = [toolbox.clone(_) for _ in offspring]
# 变异操作
# 突变
for ind1, ind2 in zip(offspring[::2], offspring[1::2]):
if random.random() < cxpb:
toolbox.mate(ind1, ind2)
del ind1.fitness.values
del ind2.fitness.values
for ind in offspring:
if random.random() < mutpb:
toolbox.mutate(ind)
del ind.fitness.values
# 评价子代
invalid_ind = [ind for ind in offspring if not ind.fitness.valid]
fitnesses = toolbox.map(toolbox.eval, invalid_ind)
for ind, fitness in zip(invalid_ind, fitnesses):
ind.fitness.values = fitness
# 环境选择
combinedPop = pop + offspring # 采用精英策略,加速收敛
pop = tools.selBest(combinedPop, popSize)
# 记录数据
record = stats.compile(pop)
logbook.record(gen=gen, nevals=len(invalid_ind), **record)
print(logbook)
结果:
# 输出结果
bestInd = tools.selBest(pop,1)[0]
bestFit = bestInd.fitness.values[0]
print('最优解为: '+str(decode(bestInd)))
print('函数最小值为: '+str(bestFit))
## 结果
## 最优解为: (2.8837195267510136, 0.46988040029289735)
## 函数最小值为: 0.002469514927004176
对迭代过程可视化:
# 可视化迭代过程
minFit = logbook.select('min')
avgFit = logbook.select('avg')
plt.plot(maxFit, 'b-', label='Minimum Fitness')
plt.plot(avgFit, 'r-', label='Average Fitness')
plt.xlabel('# Gen')
plt.ylabel('Fitness')
plt.legend(loc='best')
# 保存计算结果
plt.tight_layout()
plt.savefig('Binary encoding.png')
Binary encoding
小结
- 对于这样简单的非线性函数寻优问题，实数编码和二进制编码都可以在一定精度上求得最优解；虽然在给出的答案中实数编码得到的答案似乎更接近全局最优，但考虑到遗传算法的随机性，多运行几次，二进制编码也能得到很好的结果；
- 对于同一个问题，遗传算法对个体的编码方式并不是唯一固定的。我们可以多思考多尝试，寻找更优的编码方式，例如在这个寻优问题中，二进制编码相比实数编码，不但需要更大的存储空间，而且还额外多了一个解码的过程，因此并非很优的选择。
