Pythonの進化計算ライブラリDeap

DEAP概要

Pythonで使える遺伝的アルゴリズムライブラリDeapを紹介したいと思います。
Pythonの遺伝的アルゴリズムライブラリは他にもPyevolveというのがあるのですが、Deapの方が開発が盛んらしいので、こちらを使ってみたいと思います。
以下がDeapで使用できる主なアルゴリズムおよび機能です。

• GA 遺伝的アルゴリズム
• GP 遺伝的プログラミング
• ES 進化戦略(CMA-ESなど)
• 多目的最適化(NSGA-II, SPEA-II)
• Co-evolution
• 並列化
• 個体の中の優等生の保持
• 定期的なチェックポイント
• ベンチマーク
• 進化の系図
• Particle Swarm Optimization, Differential Evolution, Estimation of Distribution Algorithm

Example

今回はGAのExampleを解説していこうと思います。
まず、モジュールのインポートです。

import random

from deap import base
from deap import creator
from deap import tools

__creator__はbaseのクラスを継承して新たなクラスを作成します。

creator.create("FitnessMax", base.Fitness, weights=(1.0,))
creator.create("Individual", list, fitness=creator.FitnessMax)

1行目は適合度を最大化することで最適化されるような適合度クラスを作成します。
weightsが配列になっていますが、単目的最適化ではサイズ1で良く(ただし","は必要みたいです)、多目的最適化では各目的の重みを配列にします。最小化で最適化を行う場合は(-1.0,)を設定します。

toolbox = base.Toolbox()
# Attribute generator
toolbox.register("attr_bool", random.randint, 0, 1)
# Structure initializers
toolbox.register("individual", tools.initRepeat, creator.Individual, 
    toolbox.attr_bool, 100)
toolbox.register("population", tools.initRepeat, list, toolbox.individual)

__toolbox__は関数を作成します。3行目のregisterではattr_boolという関数をrandom.randintに0,1の引数を与えて作成しています。attr_boolは0か1かをランダムで生成する関数になります。
__initRepeat__はToolboxにあらかじめ用意されている関数で、1つめの引数がデータを入れるコンテナの型、2つめが個々のデータを生成する関数、3つめがコンテナのサイズです。ここではindividual(個体)という関数とpopulation(集団)という関数をinitRepeatから作成しています。
次に評価関数、交差、突然変異、選択用の関数を作成します。(評価関数の返り値に","があるのにも注意してください。)

def evalOneMax(individual):
    return sum(individual),

toolbox.register("evaluate", evalOneMax)
toolbox.register("mate", tools.cxTwoPoints)
toolbox.register("mutate", tools.mutFlipBit, indpb=0.05)
toolbox.register("select", tools.selTournament, tournsize=3)

ここまででGAの設定が完了したので、実際の進化計算ルーチンを作成します。

def main():
    random.seed(64)
    # 初期の個体群を生成
    pop = toolbox.population(n=300)
    CXPB, MUTPB, NGEN = 0.5, 0.2, 40 # 交差確率、突然変異確率、進化計算のループ回数
    
    print("Start of evolution")
    
    # 初期の個体群の評価
    fitnesses = list(map(toolbox.evaluate, pop))
    for ind, fit in zip(pop, fitnesses):
        ind.fitness.values = fit
    
    print("  Evaluated %i individuals" % len(pop))
    
    # 進化計算開始
    for g in range(NGEN):
        print("-- Generation %i --" % g)
        
        # 次世代の個体群を選択
        offspring = toolbox.select(pop, len(pop))
        # 個体群のクローンを生成
        offspring = list(map(toolbox.clone, offspring))
    
        # 選択した個体群に交差と突然変異を適応する
        # 偶数番目と奇数番目の個体を取り出して交差
        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
    
        # 適合度が計算されていない個体を集めて適合度を計算
        invalid_ind = [ind for ind in offspring if not ind.fitness.valid]
        fitnesses = map(toolbox.evaluate, invalid_ind)
        for ind, fit in zip(invalid_ind, fitnesses):
            ind.fitness.values = fit
        
        print("  Evaluated %i individuals" % len(invalid_ind))
        
        # 次世代群をoffspringにする
        pop[:] = offspring
        
        # すべての個体の適合度を配列にする
        fits = [ind.fitness.values[0] for ind in pop]
        
        length = len(pop)
        mean = sum(fits) / length
        sum2 = sum(x*x for x in fits)
        std = abs(sum2 / length - mean**2)**0.5
        
        print("  Min %s" % min(fits))
        print("  Max %s" % max(fits))
        print("  Avg %s" % mean)
        print("  Std %s" % std)
    
    print("-- End of (successful) evolution --")
    
    best_ind = tools.selBest(pop, 1)[0]
    print("Best individual is %s, %s" % (best_ind, best_ind.fitness.values))

あとは、このmain関数をif __name__ == "__main__":main()のように呼び出せばOKです。

出力結果

Start of evolution
  Evaluated 300 individuals
-- Generation 0 --
  Evaluated 189 individuals
  Min 40.0
  Max 65.0
  Avg 54.7433333333
  Std 4.46289766358
-- Generation 1 --
  Evaluated 171 individuals
  Min 44.0
  Max 70.0
  Avg 58.48
  Std 3.98533980149

...(途中省略)...

-- Generation 38 --
  Evaluated 180 individuals
  Min 89.0
  Max 100.0
  Avg 97.7466666667
  Std 2.32719191779
-- Generation 39 --
  Evaluated 196 individuals
  Min 88.0
  Max 100.0
  Avg 98.1833333333
  Std 2.33589145486
-- End of (successful) evolution --
Best individual is [True, 1, 1, 1, 1, 1, True, 1, 1, 1, 1, True, 1, 1, 1, 1, True, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, True, True, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, True, 1, 1, 1, 1, True, True, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1], (100.0,)

交差や突然変異の関数の問題なのか配列の中に1とTrueが混じってますが、最適解が求められていることが分かります。

参考

DeapのExample:http://deap.gel.ulaval.ca/doc/default/examples/ga_onemax.html