オペレーションズ・リサーチ

Pythonで遺伝的アルゴリズム

オペレーションズ・リサーチ

[PR]当サイトはアフィリエイト広告を利用しています。

はじめに

前回は遺伝的アルゴリズムについてザックリ説明しました。

今回は実際にプログラムを書きます。

とはいえ、いきなり全部のプログラムを載せて、ハイ終わりとするのはあまりにも丸投げなので、パーツごとにプログラムを作っていって、最後にまとめるようなアプローチにします。

前回の記事を読みながら取り組むと、より理解が深まると思います。

よろしければそちらもご確認ください。

また遺伝的アルゴリズムの活用例についてはこちらをご覧ください。

Pythonで作っていきます。

初期設定

まず初期設定からお話しします。

問題は前回と同様です。

まず変数の箱(遺伝子)を探索する個体分だけ用意します。

個体数は適当に決めます。本来は10個体くらいほしいですが、簡単のため3個にしています。

遺伝子の箱にランダムに数を入れます。

ただし

-5<=x1,x2,x3<=5 です。

ここまでのプログラムは次のようになります。

import random
agent = []
agentNum = 3
xNum = 3
rangeLeft = [-5,-5,-5]
rangeRight = [5,5,5]
for i in range(agentNum):
    tempX = [0,0,0]
    for j in range(xNum):
        tempX[j] = rangeLeft[j] + (rangeRight[j] - rangeLeft[j]) * random.random()
    agent.append(tempX)
print(agent)

agentが3つの個体を表します。

agentNumが個体数。

xNumが遺伝子の箱の数。変数の数でもあります。

rangeLeftとrangeRightがxの範囲です。

for文で遺伝子の箱にランダムに数を入れます。

次の図をご覧ください。

xの範囲の右側から左側を引くと、その範囲の大きさになります。

この範囲の大きさに0~1の乱数をかけると、その範囲の大きさをランダムに縮小できます。

その縮小した大きさを範囲の左側に足すと、範囲内にランダムな値を作ることができます。

それが

tempX[j] = rangeLeft[j] + (rangeRight[j] - rangeLeft[j]) * random.random()

の意味です。

agent.append(tempX)

で作成したランダムな値3つを個体の情報であるagentに追加します。

プログラムを実行した結果は次のようになります。

[[-1.7573163724126886, -4.844672388293249, -4.880300742212556], [-1.102002153547832, -3.3599328401193707, -2.366216548334669], [-2.6782797897269472, -3.7781370923601965, 1.0290968682799004]]

最後にこれを関数にして簡単に呼び出せるようにします。

def iniAgent(agent, agentNum, xNum, rangeLeft, rangeRight):
    for i in range(agentNum):
        tempX = [0,0,0]
        for j in range(xNum):
            tempX[j] = rangeLeft[j] + (rangeRight[j] - rangeLeft[j]) * random.random()
        agent.append(tempX)
    return

評価

各遺伝子の値を元に、yを計算します。

これをプログラムにすると次のようになります。

import random
agent = []
agentNum = 3
xNum = 3
rangeLeft = [-5,-5,-5]
rangeRight = [5,5,5]
def iniAgent(agent, agentNum, xNum, rangeLeft, rangeRight):…
iniAgent(agent, agentNum, xNum, rangeLeft, rangeRight)

agentY = []
for i in range(agentNum):
    agentY.append(0)
def calcY(agent,agentNum,agentY,xNum):
    for i in range(agentNum):
        sum = 0
        for j in range(xNum):
            sum += agent[i][j] * agent[i][j]
        agentY[i] = sum
    return
calcY(agent,agentNum,agentY,xNum)
print(agent)
print(agentY)

関数iniAgent()の部分は長いので「…」記号で省略表記です。

実際は上で述べた関数定義が入ります。

yの計算をするのがcalcY()です。

まあ要するに

y=x1×x1+x2×x2+x3×x3

を計算しています。

その結果をagentYに格納する関数です。

出力は次のようになります。

[[-2.6167458709119806, -2.734992315160948, 3.3386567533346856], [-0.3563569385825227, -0.8219666026397263, -4.393023150439718], [-4.633736527042172, 3.1739640619542246, 2.697647894447419]]

[25.474170833511646, 20.101271763830304, 38.82286623103841]

交叉

各遺伝子のyを計算して、次の「F:適応度」という値を求めます。

F(i)=(Imax-I(i))/(Imax-Imin)

そしてFを合計させFsを求めます。

Fs=sum(F(i))

そして0から1の乱数rを発生させます。

F(i)/Fs

を順番に足していって、rを超えたら、その時のiの個体を採用します。

まずF値を計算する関数を作ります。

def calcfitI(agentY, worst, best):
    if((worst - agentY) == 0 or (worst - best) == 0):
        return 0.01
    else:
        return (worst - agentY) / (worst - best) 

ここでworstはImax、bestはIminを表します。

Fの値にもし0が出たら、0.01などの小さい数値で置き換えます。

さらに遺伝子が収束してきて、すべてのIが同じ値になる場合も考えなければいけません。

そういう時は(Imax-Imin)が0になるのでFを0.01などの小さい数値で置き換えます。

これが

if((worst - agentY) == 0 or (worst - best) == 0):

の意味です。

それ以外は

F(i)=(Imax-I(i))/(Imax-Imin)

で計算します。

次にルーレットを行って交叉に使う2つのデータを決定する関数を作ります。

def roulette(fitIList, sumfit):
    r = random.random()
    s = 0
    index = 0
    for i in range(len(fitIList)):
        s += fitIList[i] / sumfit
        if(s > r):
            index = i
            break
    return index

ここで

fitIListは各個体のF値を格納したリスト、sumfitはFsを表します。

乱数を発生させて

F(i)/Fs

を順番に足していって、rを超えたら、その時のiの個体を採用します。

iの採用が

return index

の部分です。

実際に2つのデータを選択するために次のプログラムで、bestとworst、fitIListとsumfitを計算して、roulette()を二回実施してみます。

import random
agent = []
agentNum = 3
xNum = 3
rangeLeft = [-5,-5,-5]
rangeRight = [5,5,5]
def iniAgent(agent, agentNum, xNum, rangeLeft, rangeRight):…
iniAgent(agent, agentNum, xNum, rangeLeft, rangeRight)
agentY = []
for i in range(agentNum):
    agentY.append(0)
def calcY(agent,agentNum,agentY,xNum):…
calcY(agent,agentNum,agentY,xNum)
def calcfitI(agentY, worst, best):…
def roulette(fitIList, sumfit):…
bestIndex = agentY.index(min(agentY))
best = agentY[bestIndex]
worstIndex = agentY.index(max(agentY))
worst = agentY[worstIndex]

sumfit = 0
fitIList = []
for i in range(agentNum):
    fitIList.append(calcfitI(agentY[i], worst, best))
    sumfit += fitIList[i]
print(roulette(fitIList,sumfit))
print(roulette(fitIList,sumfit))

出力は

1

2

ちゃんと二つのデータが選択できました。

全部を結合してGAを走らせます

上記のプログラムの続きに次のようなプログラムを書きます。

iteration = 1000
crossoverR = 0.5
mutationR = 0.01
for k in range(iteration):
    sumfit = 0
    for i in range(agentNum):
        fitIList[i] = calcfitI(agentY[i], worst, best)
        sumfit += fitIList[i]
    for i in range(agentNum):
        agentTempX = []
        for p in range(xNum):
            agentTempX.append(0)
        parent1 = roulette(fitIList, sumfit)
        parent2 = roulette(fitIList, sumfit)
        for j in range(xNum):
            if(crossoverR < random.random()):
                agentTempX[j] = agent[parent1][j]
                if(mutationR > random.random()):
                    agentTempX[j] = rangeLeft[j] + (rangeRight[j] - rangeLeft[j]) * random.random()
            else:
                agentTempX[j] = agent[parent2][j]
                if(mutationR > random.random()):
                    agentTempX[j] = rangeLeft[j] + (rangeRight[j] - rangeLeft[j]) * random.random()
        agent[i] = agentTempX
    calcY(agent,agentNum,agentY,xNum)
    bestIndex = agentY.index(min(agentY))
    best = agentY[bestIndex]
    worstIndex = agentY.index(max(agentY))
    worst = agentY[worstIndex]
print(agent[bestIndex])
print(best)

解説していきます。

iteration:繰り返し回数

crossoverR:交叉の確率

mutationR:突然変異の確率

    sumfit = 0
    for i in range(agentNum):
        fitIList[i] = calcfitI(agentY[i], worst, best)
        sumfit += fitIList[i]

Fsを繰り返しの最初に計算します。

        agentTempX = []
        for p in range(xNum):
            agentTempX.append(0)

交叉した遺伝子を格納するための一時的な器です。

        parent1 = roulette(fitIList, sumfit)
        parent2 = roulette(fitIList, sumfit)

交叉する二つのデータを選びます。
そして実際に交叉します。

        for j in range(xNum):
            if(crossoverR < random.random()):
                agentTempX[j] = agent[parent1][j]
                if(mutationR > random.random()):
                    agentTempX[j] = rangeLeft[j] + (rangeRight[j] - rangeLeft[j]) * random.random()
            else:
                agentTempX[j] = agent[parent2][j]
                if(mutationR > random.random()):
                    agentTempX[j] = rangeLeft[j] + (rangeRight[j] - rangeLeft[j]) * random.random()
        agent[i] = agentTempX

crossoverR < random.random() で親1から、そうでないなら親2から遺伝子の箱をコピーします。

                if(mutationR > random.random()):
                    agentTempX[j] = rangeLeft[j] + (rangeRight[j] - rangeLeft[j]) * random.random()

突然変異で、たまにランダムな値を遺伝子の箱に入れます。

        agent[i] = agentTempX

作成した遺伝子の器を、個体の新たな遺伝子として登録します。

    calcY(agent,agentNum,agentY,xNum)
    bestIndex = agentY.index(min(agentY))
    best = agentY[bestIndex]
    worstIndex = agentY.index(max(agentY))
    worst = agentY[worstIndex]

yの値を計算して、bestとworstを登録します。

print(agent[bestIndex])
print(best)

反復が終わったら、最良の遺伝子とそのときのyを出力します。

プログラム全体

全体のコードを省略なしで載せると次のようになります。

import random
agent = []
agentNum = 10
xNum = 3
rangeLeft = [-5,-5,-5]
rangeRight = [5,5,5]
def iniAgent(agent, agentNum, xNum, rangeLeft, rangeRight):
    for i in range(agentNum):
        tempX = [0,0,0]
        for j in range(xNum):
            tempX[j] = rangeLeft[j] + (rangeRight[j] - rangeLeft[j]) * random.random()
        agent.append(tempX)
    return
iniAgent(agent, agentNum, xNum, rangeLeft, rangeRight)
agentY = []
for i in range(agentNum):
    agentY.append(0)
def calcY(agent,agentNum,agentY,xNum):
    for i in range(agentNum):
        sum = 0
        for j in range(xNum):
            sum += agent[i][j] * agent[i][j]
        agentY[i] = sum
    return
calcY(agent,agentNum,agentY,xNum)
def calcfitI(agentY, worst, best):
    if((worst - agentY) == 0 or (worst - best) == 0):
        return 0.01
    else:
        return (worst - agentY) / (worst - best)    
def roulette(fitIList, sumfit):
    r = random.random()
    s = 0
    index = 0
    for i in range(len(fitIList)):
        s += fitIList[i] / sumfit
        if(s > r):
            index = i
            break
    return index
bestIndex = agentY.index(min(agentY))
best = agentY[bestIndex]
worstIndex = agentY.index(max(agentY))
worst = agentY[worstIndex]
sumfit = 0
fitIList = []
for i in range(agentNum):
    fitIList.append(calcfitI(agentY[i], worst, best))
    sumfit += fitIList[i]

iteration = 1000
crossoverR = 0.5
mutationR = 0.01
for k in range(iteration):
    sumfit = 0
    for i in range(agentNum):
        fitIList[i] = calcfitI(agentY[i], worst, best)
        sumfit += fitIList[i]
    for i in range(agentNum):
        agentTempX = []
        for p in range(xNum):
            agentTempX.append(0)
        parent1 = roulette(fitIList, sumfit)
        parent2 = roulette(fitIList, sumfit)
        for j in range(xNum):
            if(crossoverR < random.random()):
                agentTempX[j] = agent[parent1][j]
                if(mutationR > random.random()):
                    agentTempX[j] = rangeLeft[j] + (rangeRight[j] - rangeLeft[j]) * random.random()
            else:
                agentTempX[j] = agent[parent2][j]
                if(mutationR > random.random()):
                    agentTempX[j] = rangeLeft[j] + (rangeRight[j] - rangeLeft[j]) * random.random()
        agent[i] = agentTempX
    calcY(agent,agentNum,agentY,xNum)
    bestIndex = agentY.index(min(agentY))
    best = agentY[bestIndex]
    worstIndex = agentY.index(max(agentY))
    worst = agentY[worstIndex]
print(agent[bestIndex])
print(best)

ここでagentNum = 10

としました。個体数3つだとさすがによい値が出ないので。

出力は

[0.016800221679470262, 0.014557315436035445, 0.0013207670648451852]

0.0004959073068231582

0に近い値が出ましたね。

まとめ

今回は遺伝的アルゴリズム(GA)のプログラミングをテーマに書きました。

特定の手法で各個体の値を更新することで、最小値が求まっていく過程を感じられたら幸いです。

遺伝的アルゴリズムのことをもっと知りたい、他の最適化手法が知りたいというときは次の書籍がおすすめです。ただしある程度プログラミングができる前提で書かれています。