遺伝的アルゴリズムで関数の最適値を求める（その１）

進化アルゴリズムとは

進化的アルゴリズム（evolutionary algorithm, EA）とは、最適化アルゴリズムの一つで、自然界の淘汰、遺伝子組み換え、突然変異などをモデルとして組み込んだ最適化手法である。

進化的アルゴリズムのは下記の4種類に大分される。

• 遺伝的アルゴリズム
• 遺伝的プログラミング
• 進化戦略
• 進化的プログラミング

今回は遺伝的アルゴリズム（genetic algorithm, GA）を扱う。

使用したライブラリ

下記のDEAPというpythonライブラリを用いる。
https://deap.readthedocs.io/en/master/

コード

目的関数

下記に示すBird Bunctionの最小値を探索する。
$$f(x, y) = sin(x)e^{(1-cos(y))^2}+cos(y)e^{(1-sin(x))^2}+(x-y)^2$$
Bird Functionは多峰性の関数で、(x,y)=(4.70104, 3.15294),(-1.58214, -3.13024)に大域的最適解-106.764537、その周辺に多数の局所解を持つ。

sample_GA.py

import numpy as np
def bird(x):
    x1 = x[0]
    x2 = x[1]
    t = np.sin(x1)*np.exp((1-np.cos(x2))**2) + np.cos(x2)*np.exp((1-np.sin(x1))**2) + (x1-x2)**2 
    return t,

# meshgrid上で計算
c1 = c2 = np.arange(-10, 10, 0.1)
xx, yy = np.meshgrid(c1, c2)
length = np.size(c1)  
z = np.empty((length, length))
for i in range(length):
    for j in range(length):
        xin=np.array([xx[i,j], yy[i,j]])
        z[i,j], = bird(xin)

# 描画
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

fig0 = plt.figure(figsize=(12,12), dpi=80).add_subplot(111, projection='3d')
fig0.plot_surface(xx, yy, z, alpha=1)
plt.xlim(10, -10)            
plt.ylim(-10, 10)

ライブラリのインポート

必要なライブラリをインポートしておく。
DEAPはrandomライブラリを使用するため、1試行ごとに結果が変化する。
これを避けるためには、random.seed(1)とすること。

sample_GA.py

import pandas as pd
import random
from deap import algorithms
from deap import base,creator,tools
# random.seed(1)

初期個体群の生成

初期個体を生成する関数を作っておく。
ここでは、x, yともに-10〜10の範囲で、10×10の格子状に初期個体を生成した。

sample_GA.py

def initPopulation(pcls, ind_init):
    x = np.linspace(-10, 10, 10)
    y = np.linspace(-10, 10, 10)
    inix, iniy = np.meshgrid(x,y)
    contents = np.concatenate([inix.reshape([-1,1]), iniy.reshape([-1,1])], axis=1)
    pcls = [ind_init(c) for c in contents]
    return pcls

toolboxに登録

評価（"evaluate"）、個体群（"population"）、交叉（"mate"）、突然変異（"mutate"）、選択（"select"）等の手法をtoolboxに登録しておく。

例えば、下記のtoolbox.register("mate", tools.cxBlend, alpha=0.2)は、
tool.cxblendという関数を交叉手法として登録します、という意味（alpha=0.2はtools.cxBlendの引数）。

交叉方法にブレンド交叉、突然変異にガウス変異、選択形式は10個体のトーナメントで実施する。
初期個体群を指定したい場合、"population_guess"に登録しておく。

sample_GA.py

creator.create("FitnessMax", base.Fitness, weights=(-1.0,))
creator.create("Individual", np.ndarray, fitness=creator.FitnessMax)

toolbox = base.Toolbox()
toolbox.register("evaluate", bird)
toolbox.register("population_guess", initPopulation, list, creator.Individual)
toolbox.register("mate", tools.cxBlend, alpha=0.2)
toolbox.register("mutate", tools.mutGaussian, mu=[0.0, 0.0], sigma=[2.0, 2.0], indpb=1)
toolbox.register("select", tools.selTournament, tournsize=10)

実行

ここまで設定すれば、あとは簡単なコマンドで動く。

下記の例では各世代の評価個体数、平均値、標準偏差、最大値、最小値をlogbookに記録していく。
hofとは、Holl of Fameの略で、評価が高い順に指定した数だけ保存する事ができる。

algorithms.eaSimpleの引数cxpbは交叉率、mutpbは突然変異率、ngenは世代数を示す。
1世代あたりの個体数はpopで定義したものを引き継ぎ10×10=100個体となる。
初期サンプリングと世代交代時の個体数を別々に指定することも可能。

sample_GA.py

# 初期個体群を生成
pop = toolbox.population_guess()

# 設定
hof = tools.HallOfFame(5, similar=np.array_equal)
stats = tools.Statistics(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)

# 実行
pop, logbook= algorithms.eaSimple(pop, toolbox,
                        cxpb=0.9, mutpb=0.1, ngen=20,
                        stats=stats, halloffame=hof, verbose=True)

# 上位3個体を表示
print("holl of fame.1 : %s" % (hof[0]))
print("holl of fame.2 : %s" % (hof[1]))
print("holl of fame.3 : %s" % (hof[2]))

terminal

gen	nevals	avg    	std    	min     	max   
0  	100   	80.0665	99.3005	-83.8352	414.98
1  	74    	7.49092	62.2371	-103.255	301.724
2  	89    	11.499 	108.449	-105.693	665.594
3  	84    	-60.283	53.2721	-106.71 	161.645
4  	87    	-95.3896	32.299 	-106.723	40.9594
5  	75    	-99.3516	28.4967	-106.758	23.0967
6  	73    	-104.068	15.7185	-106.764	4.33984
7  	76    	-103.476	18.6955	-106.764	10.0824
8  	80    	-101.665	22.5172	-106.764	16.8155
9  	92    	-102.631	20.6472	-106.764	26.921 
10 	77    	-102.882	19.2791	-106.764	6.34586
11 	90    	-99.4555	30.4939	-106.764	56.3788
12 	89    	-100.566	27.1489	-106.764	30.2934
13 	79    	-100.978	25.2596	-106.764	51.745 
14 	78    	-98.4071	32.1796	-106.764	85.5625
15 	76    	-105.728	10.3096	-106.764	-3.14868
16 	89    	-95.2292	38.3427	-106.764	80.6272 
17 	91    	-102.44 	25.6436	-106.764	96.6545 
18 	80    	-105.432	11.2501	-106.764	4.33866 
19 	83    	-102.271	23.504 	-106.764	67.6912 
20 	79    	-103.856	16.5553	-106.764	-3.86946
holl of fame.1 : [4.70021671 3.1529326 ]
holl of fame.2 : [4.70021671 3.15293287]
holl of fame.3 : [4.70021671 3.15293358]

大域的最適解-106.765付近に収束した。
各個体の微分情報を使わないアルゴリズムなので、微小な誤差は出てしまう。

まとめと所感

• 遺伝的アルゴリズムDEAPで関数の最適化を行った
• 乱数を使うため、試行ごとに収束するかどうかが怪しい
• 離散値で多変数の最適化向けか？