見出し画像

Processing でグラフを描く⑲　遺伝的アルゴリズム（後編）

creativival

Processing でグラフを描く　第１９回目です。

前回に引き続き、遺伝的アルゴリズムを考えていきます。遺伝的アルゴリズムは生物の進化の過程をシミュレートして、複雑な問題の解を求めるプログラムのことです。今回のミッションは、このアルゴリズムを適用して、ドット絵の頂点座標から「頂点をつなぐ最短経路」を求めることです。

前回の振り返り

tsp random

前回は「ランダムに経路を求めて、ベストルート（その時点での最短ルート）を更新する」方法（ランダム更新と呼ぶ）を考えましたが、なかなか求める最短経路には近づけないところまでやりました。

今回はベストルートを一部変更して、ベストルートを更新する方法を実装します。

ベストルートを一部入れ替える

import random

BACKGROUND_COLOR = color(0, 44, 77)
SIZE = 1000
STEP = 200
# 都市の座標リスト
cities = [
    (23, 1),
    (39, 1),
    (21, 3),
    (23, 3),
    # (25, 4),
    # (27, 4),
    (39, 3),
    (41, 3),
    # (25, 6),
    # (27, 6),
    (31, 12),
    (41, 12),
    (31, 14),
    (37, 14),
    (29, 16),
    (37, 16),
    (1, 16),
    (3, 16),
    (19, 16),
    (21, 16),
    (16, 18),
    (19, 18),
    (3, 20),
    (5, 20),
    (13, 20),
    (16, 20),
    (29, 20),
    (33, 20),
    (29, 22),
    (31, 22),
    (5, 22),
    (7, 22),
    (11, 22),
    (13, 22),
    (31, 24),
    (33, 24),
    (7, 24),
    (11, 24),
    (1, 28),
    (3, 28),
    (27, 29),
    (29, 29),
    (3, 30),
    (5, 30),
    (5, 32),
    (7, 32),
    (25, 32),
    (27, 32),
    (7, 34),
    (9, 34),
    (23, 34),
    (25, 34),
    (9, 36),
    (11, 36),
    (17, 36),
    (19, 36),
    (15, 38),
    (17, 38),
    (19, 38),
    (21, 38),
    (13, 40),
    (15, 40),
    (13, 42),
    (15, 42),
    (23, 42),
    (25, 42),
    (11, 44),
    (15, 44),
    (21, 44),
    (25, 44)
]
city_num = len(cities)
best_route = None
shortest_distance = 0
random_update, replace_update = 0, 0


def setup():
    global best_route, shortest_distance
    size(SIZE, SIZE)
    # 初期ルート
    best_route = Route()
    shortest_distance = best_route.distance


def draw():
    global best_route, shortest_distance, random_update, replace_update
    background(BACKGROUND_COLOR)
    noFill()
    for i in range(int(SIZE / STEP)):
        stroke(0, 255, 0, 50)
        line(-2000, i * STEP, 2000, i * STEP)
        line(i * STEP, -2000, i * STEP, 2000)
        for j in range(10):
            stroke(0, 255, 0, 25)
            line(-2000, i * STEP + j * STEP / 10.0, 2000, i * STEP + j * STEP / 10.0)
            line(i * STEP + j * STEP / 10.0, -2000, i * STEP + j * STEP / 10.0, 2000)

    # ルートを表示
    best_route.display()

    # 新しいルートをランダムに作成
    for _ in range(5):
        new_route = Route()
        if new_route.distance < best_route.distance:
            best_route = new_route
            shortest_distance = new_route.distance
            random_update += 1

    # ベストルートを入れ替える
    for i in range(5):
        new_route = best_route.replace(i + 1)
        new_route.calc()
        if new_route.distance < best_route.distance:
            best_route = new_route
            shortest_distance = new_route.distance
            replace_update += 1
    # テキスト表示
    fill(255)
    textSize(20)
    text('number of trials: ' + str(frameCount), 10, 30)
    text('shortest distance: ' + str(shortest_distance), 10, 60)
    text('random update: ' + str(random_update), 10, 90)
    text('replace update: ' + str(replace_update), 10, 120)


class Route:
    def __init__(self, a=[]):
        if len(a) == 0:
            self.city_nums = random.sample(
                list(range(city_num)),
                city_num
            )
        else:
            self.city_nums = a
        self.distance = 0
        self.distance_list = []
        self.calc()

    def display(self):
        stroke(127)
        beginShape()
        for i in self.city_nums:
            x = cities[i][0] * STEP / 10.0
            y = cities[i][1] * STEP / 10.0
            vertex(x, y)
            ellipse(x, y, 4, 4)
            text(str(i), x + 10, y)
        endShape(CLOSE)

    def calc(self):
        self.distance_list = []
        for i in range(city_num):
            start_city_num = self.city_nums[i]
            end_city_num = self.city_nums[(i + 1) % city_num]
            x0 = cities[start_city_num][0]
            y0 = cities[start_city_num][1]
            x1 = cities[end_city_num][0]
            y1 = cities[end_city_num][1]
            distance = dist(x0, y0, x1, y1)
            self.distance_list.append(distance)
        self.distance = sum(self.distance_list)

    def replace(self, num):
        indexes = random.sample(list(range(city_num)), num)
        replace_route = Route(self.city_nums[:])
        for i in range(num - 1):
            replace_num1 = indexes[i]
            replace_num2 = indexes[i + 1]
            replace_route.city_nums[replace_num1], replace_route.city_nums[replace_num2] = \
                replace_route.city_nums[replace_num2], replace_route.city_nums[replace_num1]
        return replace_route
tsp replace

Routeクラスに replaceメソッドを追加しました。replaceメソッドはベストルートからいくつかの都市を選び、それを入れ替えたルートを返します。そして入れ替えたルートと元のベストルートを比較して、経路がより短くなっていたらベストルートを更新します（入れ替え更新と呼ぶ）。

１０００回試行した結果をアニメーションにしました。ランダム更新に比べて入れ替え更新は、圧倒的に効率よく最短経路に近づいていきます。ランダム更新は２回に対して、入れ替え更新はなんと８９回も行われています。ベストルートの経路は 500近くまで短くすることができました。

ベストルートの長い都市間を切断する

import random

...省略...
random_update, replace_update, long_cut_update = 0, 0, 0  # 最後の変数を追加

def setup():
    ...省略...

def draw():
    global best_route, shortest_distance, random_update, replace_update, \
        long_cut_update  # 追加
    ...省略...

    # 経路の長い部分を切断（追加）
    for i, distance in enumerate(best_route.distance_list):
        if distance > 16:  # 距離が 16以上の都市を選ぶ
            new_route = best_route.cut(i)
            new_route.calc()
            if new_route.distance < best_route.distance:
                best_route = new_route
                shortest_distance = new_route.distance
                long_cut_update += 1
                
    # テキスト表示
    fill(255)
    textSize(20)
    text('number of trials: ' + str(frameCount), 10, 30)
    text('shortest distance: ' + str(shortest_distance), 10, 60)
    text('random update: ' + str(random_update), 10, 90)
    text('replace update: ' + str(replace_update), 10, 120)
    text('long cut update: ' + str(long_cut_update), 10, 150)  # 追加


class Route:
   ...省略...
   # cutメソッドを追加
    def cut(self, num):
        cut_route = Route(self.city_nums[:])
        cut_num1 = num
        cut_num2 = (num + random.randint(2, city_num - 1)) % city_num
        cut_route.city_nums[cut_num1], cut_route.city_nums[cut_num2] = \
            cut_route.city_nums[cut_num2], cut_route.city_nums[cut_num1]
        return cut_route
tsp long cut

さらに別のアプローチを考えました。経路の変化を見ていると、斜めに長く伸びている線がなかなか切断されないことがわかりました。この長い経路を優先的に切断して、ベストルートを更新させてみます（切断更新と呼ぶ）。

Routeクラスに cutメソッドを追加しました。cutメソッドは特定の都市間の経路を切断して、新しいルートを返します。ベストルートのインスタンス変数distance_list から都市間の経路が 16 以上の都市を選んで、cutメソッドをつかって新しいルートを作成します。そして新しいルートとベストルートを比較して、より短くなっていたらベストルートを更新します。

１０００回試行した結果のアニメーションを見ると、切断更新は最高の効率を示しました。ランダム更新２回、入れ替え更新28回に対して、切断更新は106回も実施されています。経路も400台まで短くすることができました。

最高性能を示した切断更新ですが、ある程度短くなると更新が止まってしまいます。そこで出てくるのが、今回の主題である「遺伝的アルゴリズム」です。

（コラム）「都市間の経路が 16 以上の都市」を選んだ理由が気になった方がおられるかもしれません（するどい！）。ドット絵（恐竜）で、頂点間の距離が一番長いものを探すと、16でした（頂点0 と頂点1をつないだ線）。つまり長さ 16 以上の経路はすべて更新の対象とすることができるからです。

遺伝的アルゴリズム

import random
from operator import attrgetter

BACKGROUND_COLOR = color(0, 44, 77)
SIZE = 1000
STEP = 200
# 都市の座標リスト
cities = [
    (23, 1),
    (39, 1),
    (21, 3),
    (23, 3),
    # (25, 4),
    # (27, 4),
    (39, 3),
    (41, 3),
    # (25, 6),
    # (27, 6),
    (31, 12),
    (41, 12),
    (31, 14),
    (37, 14),
    (29, 16),
    (37, 16),
    (1, 16),
    (3, 16),
    (19, 16),
    (21, 16),
    (16, 18),
    (19, 18),
    (3, 20),
    (5, 20),
    (13, 20),
    (16, 20),
    (29, 20),
    (33, 20),
    (29, 22),
    (31, 22),
    (5, 22),
    (7, 22),
    (11, 22),
    (13, 22),
    (31, 24),
    (33, 24),
    (7, 24),
    (11, 24),
    (1, 28),
    (3, 28),
    (27, 29),
    (29, 29),
    (3, 30),
    (5, 30),
    (5, 32),
    (7, 32),
    (25, 32),
    (27, 32),
    (7, 34),
    (9, 34),
    (23, 34),
    (25, 34),
    (9, 36),
    (11, 36),
    (17, 36),
    (19, 36),
    (15, 38),
    (17, 38),
    (19, 38),
    (21, 38),
    (13, 40),
    (15, 40),
    (13, 42),
    (15, 42),
    (23, 42),
    (25, 42),
    (11, 44),
    (15, 44),
    (21, 44),
    (25, 44)
]

city_num = len(cities)
best_route = None
shortest_distance = 0
routes = []
route_num = 1000
random_update, replace_update, long_cut_update, genetic_update = 0, 0, 0, 0


def setup():
    global best_route, shortest_distance
    size(SIZE, SIZE)
    # 初期ルート
    best_route = Route()
    shortest_distance = best_route.distance
    # ルートグループを作る（遺伝的アルゴリズム）
    for i in range(route_num):
        routes.append(Route())


def draw():
    global best_route, shortest_distance, random_update, replace_update, \
        long_cut_update, routes, genetic_update
    background(BACKGROUND_COLOR)
    noFill()
    for i in range(int(SIZE / STEP)):
        stroke(0, 255, 0, 50)
        line(-2000, i * STEP, 2000, i * STEP)
        line(i * STEP, -2000, i * STEP, 2000)
        for j in range(10):
            stroke(0, 255, 0, 25)
            line(-2000, i * STEP + j * STEP / 10.0, 2000, i * STEP + j * STEP / 10.0)
            line(i * STEP + j * STEP / 10.0, -2000, i * STEP + j * STEP / 10.0, 2000)
    pushMatrix()

    # # ルートを表示
    best_route.display()

    # 新しいルートをランダムに作成
    for _ in range(5):
        new_route = Route()
        if new_route.distance < best_route.distance:
            best_route = new_route
            shortest_distance = new_route.distance
            random_update += 1
            # ベストルートを加える
            routes.append(best_route)

    # ベストルートを入れ替える
    for i in range(5):
        new_route = best_route.replace(i + 1)
        new_route.calc()
        if new_route.distance < best_route.distance:
            best_route = new_route
            shortest_distance = new_route.distance
            replace_update += 1
            # ベストルートを加える
            routes.append(best_route)

    # 経路の長い部分を切断
    for i, distance in enumerate(best_route.distance_list):
       6 if distance > 1:  # 距離が 16以上の都市を選ぶ
            new_route = best_route.cut(i)
            new_route.calc()
            if new_route.distance < best_route.distance:
                best_route = new_route
                shortest_distance = new_route.distance
                long_cut_update += 1
                # ベストルートを加える
                routes.append(best_route)

    # 遺伝的アルゴリズム
    # 経路の短い順に並べる
    routes.sort(key=attrgetter('distance'))
    # 上位 1％の精鋭部隊を作る（次の世代にそのまま移行できる）
    expert_routes = routes[:int(len(routes) / 100)]
    # 上位 10％の上流階級を作る（次の世代を作ることができる）
    better_routes = routes[:int(len(routes) / 10)]
    # 交配のルールで上流階級の子供たちを作成
    new_routes = []
    while len(new_routes) < route_num:
        partners = random.sample(better_routes, 2)
        new_routes.append(partners[0].cross_breed(partners[1]))
    for new_route in new_routes:
        if new_route.distance < best_route.distance:
            # print('update: ', new_route.distance, best_route.distance)
            best_route = new_route
            shortest_distance = best_route.distance
            genetic_update += 1
    # 次の世代（精鋭部隊 + 上流階級の子供たち）
    routes = expert_routes + new_routes
    popMatrix()

    fill(255)
    textSize(20)
    text('number of trials: ' + str(frameCount), 10, 30)
    text('shortest distance: ' + str(shortest_distance), 10, 60)
    text('random update: ' + str(random_update), 10, 90)
    text('replace update: ' + str(replace_update), 10, 120)
    text('long cut update: ' + str(long_cut_update), 10, 150)
    text('genetic update: ' + str(genetic_update), 10, 180)


class Route:
    def __init__(self, a=[]):
        if len(a) == 0:
            self.city_nums = random.sample(
                list(range(city_num)),
                city_num
            )
        else:
            self.city_nums = a
        self.distance = 0
        self.distance_list = []
        self.calc()

    def display(self):
        stroke(127)
        beginShape()
        for i in self.city_nums:
            x = cities[i][0] * STEP / 10.0
            y = cities[i][1] * STEP / 10.0
            vertex(x, y)
            ellipse(x, y, 4, 4)
            text(str(i), x + 10, y)
        endShape(CLOSE)

    def calc(self):
        self.distance_list = []
        for i in range(city_num):
            start_city_num = self.city_nums[i]
            end_city_num = self.city_nums[(i + 1) % city_num]
            x0 = cities[start_city_num][0]
            y0 = cities[start_city_num][1]
            x1 = cities[end_city_num][0]
            y1 = cities[end_city_num][1]
            distance = dist(x0, y0, x1, y1)
            self.distance_list.append(distance)
        self.distance = sum(self.distance_list)

    def replace(self, num):
        indexes = random.sample(list(range(city_num)), num)
        replace_route = Route(self.city_nums[:])
        for i in range(num - 1):
            replace_num1 = indexes[i]
            replace_num2 = indexes[i + 1]
            replace_route.city_nums[replace_num1], replace_route.city_nums[replace_num2] = \
                replace_route.city_nums[replace_num2], replace_route.city_nums[replace_num1]
        return replace_route

    def cut(self, num):
        cut_route = Route(self.city_nums[:])
        cut_num1 = num
        cut_num2 = (num + random.randint(2, city_num - 1)) % city_num
        cut_route.city_nums[cut_num1], cut_route.city_nums[cut_num2] = \
            cut_route.city_nums[cut_num2], cut_route.city_nums[cut_num1]
        return cut_route

    def cross_breed(self, partner):
        # 並びの切断個所をランダムに選ぶ
        cut_id = random.randint(1, city_num)
        part1 = self.city_nums[:][:cut_id]
        # 50% の確率で並びを反転する
        if random.random() > 0.5:
            part1.reverse()
        # 1% の確率でパートナーを変更
        if random.random() > 0.99:
            partner = Route()
        # パートナーから不足している id を抜き出す
        part2 = [c for c in partner.city_nums[:] if c not in part1]
        # 交配
        return Route(part1 + part2)
tsp genetic

上記のコードが遺伝的アルゴリズム（交配更新と呼ぶ）を取り入れた完成形です。比較のため、ランダム更新、入れ替え更新、切断更新も同時に実行しています。左上の数値を見れば、各更新の効率がわかります。１０００回試行までの数値は、ランダム更新（0）、入れ替え更新（6）、切断更新（50）、交配更新（115）となっています。たった1000回でほぼ恐竜を再現できましたから、かなり性能が上がりました。

遺伝的アルゴリズムの実装方法を説明します。

  1. 1000個のランダムなルートを作成（ルートグループ）

  2. ルートグループを経路の短さによって並べる

  3. 上位 10% の上流階級を作成（子供を作ることができる）

  4. 上流階級の子供とベストルートを選んで、経路が短いときはベストルートを更新

  5. 次の世代（ルートグループ）として、上流階級の子供を選ぶ

  6. 2 ～ 6 を繰り返す

交配のルール

    def cross_breed(self, partner):
        # 並びの切断個所をランダムに選ぶ
        cut_id = random.randint(1, city_num)
        part1 = self.city_nums[:][:cut_id]
        # 50% の確率で並びを反転する
        if random.random() > 0.5:
            part1.reverse()
        # 1% の確率でパートナーを変更
        if random.random() > 0.99:
            partner = Route()
        # パートナーから不足している id を抜き出す
        part2 = [c for c in partner.city_nums[:] if c not in part1]
        # 交配
        return Route(part1 + part2)

Routeクラスに実装した cross_breedメソッドについて説明します。

あるルート（パートナー１)と別のルート（パートナー２）から新しいルート（子供）を作成します。ルートの並びを遺伝子と考えて、優秀な遺伝子同士を掛け合わせて、さらに優秀な遺伝子を作っていくことを目的としています。

  1. パートナー１の遺伝子を切断する箇所をランダムに選ぶ

  2. パートナー１の遺伝子を切断後、50% の確率で反転させる

  3. 1% の確率でパートナー２をランダム遺伝子を持つパートナーに変更（更新が止まってしまうのを防ぐため）

  4. パートナー２の遺伝子から不足している遺伝子を抜き出す（並びはそのまま）

  5. パートナー１の遺伝子とパートナー２の遺伝を持つ子供を作成する

まとめ

遺伝的アルゴリズムの応用として、ドット絵（恐竜）の頂点をつなぐ経路を選ぶプログラムを作成してきました。ランダム更新、入れ替え更新、切断更新、交配更新の４つの考え方で、恐竜の形をほぼ再現できることがわかりました。

遺伝的アルゴリズムは強力です。この考え方を身に着ければ、コンピューターが何万、何億、何京年かかる難題でも、有限の時間で答えを導き出すことができます。

つまり私たちは遺伝的アルゴリズムの考え方を使って、世界最速コンピューター「富岳」にも勝てるということです。数学って、すごくないですか！

前の記事
Processing でグラフを描く⑱　遺伝的アルゴリズム（前編）
次の記事
Processing でグラフを描く⑳　ライフゲーム

その他のタイトルはこちら


この記事が気に入ったらサポートをしてみませんか？