基于DEAP库的Python进化算法从入门到入土 --(四)遗传算法的改进

前言

前面一节我们尝试了用GA求解TSP问题,简单遗传算法总是不能很好收敛到一个较优的解,在用时和求解精度上都被贪心算法吊打。在末尾我们总结了三个可能的改进方向,这次我们想要沿着这三个方向试着改进简单遗传算法,以提升其性能。

精英保存策略

简介

在简单遗传算法中,我们总是用产生的育种后代族群完全取代父代族群。在真实的生物进化中,总是会有一部分后代由于不能适应环境而夭折,留下来的后代往往是相对来说更加优秀的个体。仿照这个原理,我们可以在SGA中加入自然选择,用精英保存策略(Eliteness preservation strategy)来筛选产生的育种后代并保留其中较为优秀的个体,这样可以使得算法更快收敛。

精英保存策略的步骤如下:

  • 在选择育种后代时,选择出比原始族群数量更大的一个育种族群,即n_{offspring}>n_{pop}
  • 对该育种族群进行变异操作 -- 交叉和突变;
  • 对变异后的育种族群进行自然选择,只取其中n_{pop}个个体以保持种群规模。

代码实现

为了获取更大的自由度,此处不再调用DEAP中内置的eaSimple函数,而是自行用注册的工具箱实现。代码如下:

def GA_improved(cities, nCities, npop = 100, cxpb = 0.5, mutpb = 0.2, ngen = 200):
    global cityDist
    ## 问题定义
    creator.create('FitnessMin', base.Fitness, weights=(-1.0,)) # 最小化问题
    creator.create('Individual', list, fitness=creator.FitnessMin)

    ## 定义个体编码
    toolbox = base.Toolbox()
    toolbox.register('indices', random.sample, range(nCities), nCities) # 创建序列
    toolbox.register('individual', tools.initIterate, creator.Individual, toolbox.indices)

    ## 生成族群
    toolbox.register('population', tools.initRepeat, list, toolbox.individual)
    pop = toolbox.population(npop)

    ## 注册所需工具
    cityDist = cityDistance(cities)
    toolbox.register('evaluate', routeDistance)
    toolbox.register('select', tools.selTournament, tournsize = 2)
    toolbox.register('mate', tools.cxOrdered)
    toolbox.register('mutate', tools.mutShuffleIndexes, indpb = 0.2)

    ## 数据记录
    stats = tools.Statistics(key=lambda ind: ind.fitness.values)
    stats.register('avg', np.mean)
    stats.register('min', np.min)
    logbook = tools.Logbook()
    logbook.header = ['gen', 'nevals'] + (stats.fields)

    ## 实现遗传算法
    # 评价族群
    invalid_ind = [ind for ind in pop if not ind.fitness.valid]
    fitnesses = toolbox.map(toolbox.evaluate, invalid_ind)
    for ind, fit in zip(invalid_ind, fitnesses):
        ind.fitness.values = fit
    # 记录数据
    record = stats.compile(pop)
    logbook.record(gen=0, nevals=len(invalid_ind), **record)
    
    for gen in range(1, ngen+1): # 方便输出数据好看
        # 配种选择
        offspring = toolbox.select(pop, 2*npop)
        offspring = [toolbox.clone(_) for _ in offspring] #一定要复制,否则在交叉和突变这样的原位操作中,会改变所有select出来的同个体副本
        # 变异操作 - 交叉
        for child1, child2 in zip(offspring[::2], offspring[1::2]):
            if random.random() < cxpb:
                toolbox.mate(child1, child2)
                del child1.fitness.values
                del child2.fitness.values
        # 变异操作 - 突变
        for mutant in offspring:
            if random.random() < mutpb:
                toolbox.mutate(mutant)
                del mutant.fitness.values
        # 评价当前没有fitness的个体,确保整个族群中的个体都有对应的适应度
        invalid_ind = [ind for ind in offspring if not ind.fitness.valid]
        fitnesses = toolbox.map(toolbox.evaluate, invalid_ind)
        for ind, fit in zip(invalid_ind, fitnesses):
            ind.fitness.values = fit
        # 环境选择 - 保留精英
        pop = tools.selBest(offspring, npop, fit_attr='fitness') # 选择精英,保持种群规模
#         pop[:] = offspring
        # 记录数据
        record = stats.compile(pop)
        logbook.record(gen=gen, nevals=len(invalid_ind), **record)
    return pop, logbook

与上节相同,在随机生成的规模为30的TSP问题上进行测试:

nCities = 30
cities = genCity(nCities) # 随机生成nCities个城市坐标
resultPopGA_improved, logbookGA_improved = GA_improved(cities, nCities)
plt.plot(logbookGA_improved.select('min'), 'r-', label='Minimum Fitness')
plt.plot(logbookGA_improved.select('avg'), 'b-', label='Average Fitness')
plt.ylabel('Fitness')
plt.xlabel('# Iteration')
plt.legend(loc='best')
plt.tight_layout()
plt.title(f'GA with eliteness preservation strategy iterations, Problem size:{nCities}')

改进后的GA迭代过程如下:

GA_eliteness

相比于SGA的迭代过程:

SGA

可以看到在加入了精英保存策略之后,遗传算法的收敛速度大大加快了。

在规模为30的TSP问题上运行10次,族群规模100,迭代代数200,对比改进前后所得解的情况:

# SGA GA_Eliteness preservation
1 6035.558845495109 5470.286031491102
2 5292.642425419087 5068.347131947512
3 5805.933713727669 5817.778732730951
4 5767.204287839886 4815.938824190787
5 5270.255845888708 5543.729805228135
6 6137.969866274479 4887.10488671049
7 5300.2776396096715 4684.56020813833
8 5405.8737178780775 4790.2302076385295
9 6278.724906181772 5613.217834196187
10 6360.790872210711 5138.992421217108
summary 5765.52321\pm405.67208 5183.01861\pm379.57487

可以看到从统计意义上,用了精英保存策略之后,在同样代数下,遗传算法的寻优能力得到了提升。并且有4次计算结果已经接近或超过了nnTSP+2-OPT(参考上一节中的数据)。虽然与repNNTSP+2-OPT仍然有一定差距,但是已经获得了不小的提升。

结合局部搜索算法

遗传算法依靠较大的种群数和变异操作,能够在很大的空间内寻找最优解,它的一个重要优点时不需要问题的具体信息。但是正由于没有利用问题的具体信息,导致GA的局部搜索能力差,在解空间较大的优化问题中,往往需要“人海战术”,产生很大的族群来尝试寻找最优解。

在上一篇文章中我们已经知道2-OPT这种局部搜索方法能够通过uncross路径来改善解的质量,我们可以在变异操作之后,用2-OPT改善族群中的优秀个体,提高搜索效率。

代码实现如下:

def opt(route, k=2):
    # 用2-opt算法优化路径
    # 输入:cityDist -- [n,n]矩阵,记录了城市间的距离
    # route -- sequence,记录路径
    # 输出: 优化后的路径optimizedRoute及其路径长度
    nCities = len(route) # 城市数
    optimizedRoute = route # 最优路径
    minDistance = routeDistance(route) # 最优路径长度
    for i in range(1,nCities-2):
        for j in range(i+k, nCities):
            if j-i == 1:
                continue
            reversedRoute = route[:i]+route[i:j][::-1]+route[j:]# 翻转后的路径
            reversedRouteDist = routeDistance(reversedRoute)
            # 如果翻转后路径更优,则更新最优解
            if  reversedRouteDist < minDistance:
                minDistance = reversedRouteDist
                optimizedRoute = reversedRoute
    return optimizedRoute, minDistance

def GA_eliteness_2Opt(cities, nCities, npop = 100, cxpb = 0.5, mutpb = 0.2, ngen = 200):
    global cityDist
    ## 问题定义
    creator.create('FitnessMin', base.Fitness, weights=(-1.0,)) # 最小化问题
    creator.create('Individual', list, fitness=creator.FitnessMin)

    ## 定义个体编码
    toolbox = base.Toolbox()
    toolbox.register('indices', random.sample, range(nCities), nCities) # 创建序列
    toolbox.register('individual', tools.initIterate, creator.Individual, toolbox.indices)

    ## 生成族群
    toolbox.register('population', tools.initRepeat, list, toolbox.individual)
    pop = toolbox.population(npop)

    ## 注册所需工具
    cityDist = cityDistance(cities)
    toolbox.register('evaluate', routeDistance)
    toolbox.register('select', tools.selTournament, tournsize = 2)
    toolbox.register('mate', tools.cxOrdered)
    toolbox.register('mutate', tools.mutShuffleIndexes, indpb = 0.2)
    toolbox.register('localOpt', opt) # 注册2-Opt

    ## 数据记录
    stats = tools.Statistics(key=lambda ind: ind.fitness.values)
    stats.register('avg', np.mean)
    stats.register('min', np.min)
    logbook = tools.Logbook()
    logbook.header = ['gen', 'nevals'] + (stats.fields)

    ## 实现遗传算法
    # 评价族群
    invalid_ind = [ind for ind in pop if not ind.fitness.valid]
    fitnesses = toolbox.map(toolbox.evaluate, invalid_ind)
    for ind, fit in zip(invalid_ind, fitnesses):
        ind.fitness.values = fit
    # 记录数据
    record = stats.compile(pop)
    logbook.record(gen=0, nevals=len(invalid_ind), **record)
    
    for gen in range(1, ngen+1): # 方便输出数据好看
        # 配种选择
        offspring = toolbox.select(pop, 2*npop)
        offspring = [toolbox.clone(_) for _ in offspring] #一定要复制,否则在交叉和突变这样的原位操作中,会改变所有select出来的同个体副本
        # 变异操作 - 交叉
        for child1, child2 in zip(offspring[::2], offspring[1::2]):
            if random.random() < cxpb:
                toolbox.mate(child1, child2)
                del child1.fitness.values
                del child2.fitness.values
        # 变异操作 - 突变
        for mutant in offspring:
            if random.random() < mutpb:
                toolbox.mutate(mutant)
                del mutant.fitness.values
        # 评价当前没有fitness的个体,确保整个族群中的个体都有对应的适应度
        invalid_ind = [ind for ind in offspring if not ind.fitness.valid]
        fitnesses = toolbox.map(toolbox.evaluate, invalid_ind)
        for ind, fit in zip(invalid_ind, fitnesses):
            ind.fitness.values = fit
        # 环境选择 - 保留精英
        pop = tools.selBest(offspring, npop, fit_attr='fitness') # 选择精英,保持种群规模
        
        # 对族群中的精英进行优化
        nOpt = int(npop/10)
        toBeOpt = tools.selBest(pop, nOpt)
        for ind in toBeOpt:
            ind[:], ind.fitness.values = toolbox.localOpt(ind)
#             print(ind.fitness.values)

        # 记录数据
        record = stats.compile(pop)
        logbook.record(gen=gen, nevals=len(invalid_ind), **record)
    return pop, logbook

训练过程如图:

GA_eliteness_2Opt

给出的最优路线:

最优路径为:[1, 29, 22, 3, 28, 19, 16, 2, 14, 7, 0, 15, 17, 8, 25, 27, 21, 20, 6, 12, 5, 26, 4, 9, 10, 23, 18, 24, 13, 11, 1]
最优路径距离为:(4400.1228564415105,)
tourGA_eliteness_2Opt_30cities_tour

可以看到,在加入局部搜索之后,获得的解质量大幅提升,已经超越了repNNTSP+2-OPT的水平,另外迭代收敛所需的代数也大大减小了。

超参数调参

遗传算法中的超参数及其一般取值

一般遗传算法中的超参数有:

参数 参数过小 参数过大 常用值
族群规模 搜索速度过慢,可能较难找到最优解 占用计算资源 50-100
迭代步数 在收敛之前就结束搜索 浪费计算资源 100-500
交叉概率 族群不能得到有效更新,寻优能力受限 可能破坏最优解 0.5-0.99
突变概率 容易陷入局部最优 可能破坏最优解 0.0001-0.2

通常来说这些超参数是通过经验和反复试错获得的,并且超参数的选择和问题的类型紧密相关。

除了这些一般超参数以外,各项操作--育种选择、变异和环境选择的方法选择也可以视为一些额外的超参数。

超参数寻优

迭代步数 -- 停止准则

GA的停止准则通常包含以下几种:

  • 当适应度函数收敛到某个特定值:一般在对问题有深刻了解的情况下或者问题有特定约束时才能使用;
  • 当最佳适应度在较长时间内不再改变:这是一个常用的准则,但是如果以此作为准则,有在局部最优提前结束的可能;
  • 当族群没有进化潜力时:这也是较为普遍采用的准则,一般来说可以通过计算族群的variance来确定族群的多样性是否还能支持进一步进化。当variance小于某个值时,可以判定已经收敛。
  • 当迭代步数达到预设值:这是为了防止计算时间无限拖延的一个安全设定。迭代步数达到预设值并不能保证遗传算法收敛到最优,因此通常只是作为一个最后保险。

此处采用第三种准则,当族群适应度的标准差小于1\times10^{-9}时,结束计算。结果如下:

StoppingCriterion_probsize_100

可以看到对于规模为100的TSP问题,只用了不到100代就提前停止了计算。这就比迭代到200步再停止省下了一半多的计算时间。

Tour_StoppingCriterion_probsize_100

可以从获得的路线肉眼判定,它应该已经非常接近最优路径了。计算的提前结束并没有降低解的质量。

族群规模

各路论文中推荐的族群规模千差万别,从几十到几百,不一而足。这篇文章认为一种"rule of thumb"是将族群规模选定为问题维度的10倍。但是该规则并不能scalable的,如果面对一个规模300的TSP问题,难道将族群规模设置为3000?这显然是一个非常过剩的数字,并且有研究者认为单纯增加族群规模未必有用,参见这篇文章

在笔者的尝试下,对于一个规模100的TSP问题,将族群规模从100提升到1000,计算时间增加了13倍,但是解的质量没有明显提升。下面列出单纯改变族群规模对计算的影响,TSP问题规模为100:

族群规模 迭代停止 最优路径长度 计算时间
50 106 7616.677 52.9 s
100 96 7522.976 1 min 35 s
200 84 7795.912 2 min 44 s
300 100 7738.037 4 min 44 s
1000 96 7319.198 13 min 47 s

考虑到用时和精度,100 应该是一个合理的取值。

变异概率 -- 自适应遗传算法

在通常的遗传算法中,交叉概率和突变概率是固定的值,与个体以及迭代过程都没有关系。在M.Srinivas 和 L .M. Patnaik 1994年的一篇文章 Adaptive probabilities of crossover and mutation in genetic algorithms中,他们提出了一些不同的思考:

  • 不同适应度的个体不应该具有相同的交叉概率p_c与突变概率p_m -- 对于优秀的个体,应当降低变异概率使其能最大限度保留;而对于不良个体,则应当增大变异概率,使其能脱出不良状态。
  • 作为族群整体来说,在求解的前期需要较大变异概率以便将个体在解空间中扩散出去,帮助寻优;而在求解的后期则需要较小的变异概率,加速收敛。

在这样思想的指导下,他们设计了一套自适应的交叉概率p_c与突变概率p_m计算规则:
p_c=k_1(f_{max}-f^\prime)/(f_{max}-\overline{f}),\qquad f^\prime\ge\overline{f} \\ p_c=k_3, \qquad f^\prime<\overline{f} \\ p_m=k_2(f_{max}-f)/(f_{max}-\overline{f}),\qquad f\ge\overline{f}\\ p_m = k_4, f>\overline{f}
其中f^\prime代表交叉操作中选择的两个个体中较大的适应度,\overline{f}代表族群适应度的平均值,f_{max}代表族群中最大的适应度。k_1,k_2,k_3,k_4是常数。在这篇文章中,他们推荐的值为1.0, 0.5, 1.0, 0.5。

用这套参数实测结果如下:

AGAiterations_probsize_100

可以看到,由于突变概率过大,算法在这套参数下完全无法收敛。转而尝试论文Optimal reactive power dispatch using an adaptive genetic algorithm中设定的一套参数:0.85, 0.5, 1.0, 0.05。降低了突变概率之后使得能够正常收敛,但是迭代步数相比未使用适应性策略会普遍偏大:

AGAiterations_probsize_100_paramset2
Tour_AGA_probsize_100_paramset2

优化后的算法测试

多核运行

遗传算法由于其较大的计算量,对于大规模问题通常来说需要很长的计算时间,为了减少等待时间,这里采用多核计算。对于local parallelization,可以用python的multiprocessing模块或者concurrent.futures来实现;对于cluster parallelization,可以用SCOOP来实现。

import multiprocessing

nbCpu = multiprocessing.cpu_count()
pool = multiprocessing.Pool(processes = nbCpu) # 调用所有CPU
print('Conducting calculation with %d workers'%nbCpu)
toolbox.register('map', pool.map)
# ....#
pool.close() # 运行完毕后记得关闭计算池

注意通过以上方式DEAP只会对评价进行并行计算,而其他操作如交叉和突变仍然按常规执行。如果评价函数非常简单,那么并行计算并不会节约太多时间。

如果是用控制台跑单独的python file,需要将multiprocess放在main中,参考这里

算法测试

在问题规模30,100的TSP问题中测试改进后的算法,以之前测试中表现较好的repNNTSP+2Opt作为参照:

问题规模30:

# Run 收敛代数 最优路径长度 用时
1 37 4247.60993 1.55 s
2 32 4294.63696 1.35 s
3 29 4269.27221 1.32 s
4 36 4190.33008 1.63 s
5 31 4190.33008 1.36 s
总结 4238.43585\pm42.00462
repNNTSP+2Opt 4479.93341 83.41 ms

问题规模100:

# Run 收敛代数 最优路径长度 用时
1 156 7533.05931 2 min 26 s
2 NA 7543.14155 3 min 12 s
3 144 7730.58752 2 min 19 s
4 170 7431.41310 2 min 41 s
5 159 7589.11515 2 min 27 s
总结 7565.46333\pm97.30689
repNNTSP+2Opt 8532.54589 353 ms

可以看到,我们改进之后的遗传算法能够在计算精度上远远超过之前的benchmark repNNTSP+2OPT,哪怕是最差的一次运算精度也较之更优。

但是当问题规模扩张到300以上时,计算用时非常之长,哪怕使用多核并行计算,一次运算也需要半小时以上。

小结

在本文中,我们从三个方面探索了GA的改进,使之在求解TSP时精度大幅提高:

  • 精英保存策略能使GA收敛速度大大加快,但是也有早熟,使得算法陷入局部最优的风险;
  • 结合局部搜索之后,GA从盲目搜索到能够利用问题的性质,解的质量得到大幅改善;
  • 结合了自适应算法之后,GA脱出局部最优的能力增强,但是迭代收敛需要的步数也增加了。

在进行了这三方面改善之后,GA求解TSP的能力大幅提高,从上一节的被各种吊打,已经超越了benchmark repNNTSP+2opt。

在问题规模继续扩大到300以上时,当前很难在较短时间内完成计算,需要尝试将选择、突变、评价操作进行并行计算,提高计算速度。

©著作权归作者所有,转载或内容合作请联系作者
  • 序言:七十年代末,一起剥皮案震惊了整个滨河市,随后出现的几起案子,更是在滨河造成了极大的恐慌,老刑警刘岩,带你破解...
    沈念sama阅读 212,542评论 6 493
  • 序言:滨河连续发生了三起死亡事件,死亡现场离奇诡异,居然都是意外死亡,警方通过查阅死者的电脑和手机,发现死者居然都...
    沈念sama阅读 90,596评论 3 385
  • 文/潘晓璐 我一进店门,熙熙楼的掌柜王于贵愁眉苦脸地迎上来,“玉大人,你说我怎么就摊上这事。” “怎么了?”我有些...
    开封第一讲书人阅读 158,021评论 0 348
  • 文/不坏的土叔 我叫张陵,是天一观的道长。 经常有香客问我,道长,这世上最难降的妖魔是什么? 我笑而不...
    开封第一讲书人阅读 56,682评论 1 284
  • 正文 为了忘掉前任,我火速办了婚礼,结果婚礼上,老公的妹妹穿的比我还像新娘。我一直安慰自己,他们只是感情好,可当我...
    茶点故事阅读 65,792评论 6 386
  • 文/花漫 我一把揭开白布。 她就那样静静地躺着,像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上,一...
    开封第一讲书人阅读 49,985评论 1 291
  • 那天,我揣着相机与录音,去河边找鬼。 笑死,一个胖子当着我的面吹牛,可吹牛的内容都是我干的。 我是一名探鬼主播,决...
    沈念sama阅读 39,107评论 3 410
  • 文/苍兰香墨 我猛地睁开眼,长吁一口气:“原来是场噩梦啊……” “哼!你这毒妇竟也来了?” 一声冷哼从身侧响起,我...
    开封第一讲书人阅读 37,845评论 0 268
  • 序言:老挝万荣一对情侣失踪,失踪者是张志新(化名)和其女友刘颖,没想到半个月后,有当地人在树林里发现了一具尸体,经...
    沈念sama阅读 44,299评论 1 303
  • 正文 独居荒郊野岭守林人离奇死亡,尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 36,612评论 2 327
  • 正文 我和宋清朗相恋三年,在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 38,747评论 1 341
  • 序言:一个原本活蹦乱跳的男人离奇死亡,死状恐怖,灵堂内的尸体忽然破棺而出,到底是诈尸还是另有隐情,我是刑警宁泽,带...
    沈念sama阅读 34,441评论 4 333
  • 正文 年R本政府宣布,位于F岛的核电站,受9级特大地震影响,放射性物质发生泄漏。R本人自食恶果不足惜,却给世界环境...
    茶点故事阅读 40,072评论 3 317
  • 文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹,春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 30,828评论 0 21
  • 文/苍兰香墨 我抬头看了看天上的太阳。三九已至,却和暖如春,着一层夹袄步出监牢的瞬间,已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 32,069评论 1 267
  • 我被黑心中介骗来泰国打工, 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留,地道东北人。 一个月前我还...
    沈念sama阅读 46,545评论 2 362
  • 正文 我出身青楼,却偏偏与公主长得像,于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子,可洞房花烛夜当晚...
    茶点故事阅读 43,658评论 2 350

推荐阅读更多精彩内容