遺伝的アルゴリズムを使って「離れていても心は一つ」を自動化してみた in Python

先日バズった岐阜新聞の広告「離れていても、心はひとつ」をご存知でしょうか。(yahoo! newsなどで紹介されています)

一見すると黒地に白い水玉模様なのですが、離れて（あるいはメガネを外して）見ると白地に書かれた黒い文字が浮かび上がってきます。

同じように「黒い背景に」「白い円を配置し」「結果として黒い字が読める」画像を自動で生成できたら面白いなと思ってやってみました。¹²

結論から

遺伝的アルゴリズムを使って、それなりの画像を生成できたと思います。

近くで見ると	離れると

実装は[github](https://github.com/tancematrix/figure-ground-ga)

本記事の内容

• どうやったら「離れたら読める画像」の生成ができるか、遺伝的アルゴリズムを採用するに至るまでの思考の経緯
• 遺伝的アルゴリズム実装例略解
• 結果と改良の過程

どうやったら「離れたら読める画像」の生成ができるか?

準備: プログラムの流れ

プログラムの概要はこんな感じです。

1. 目標画像を用意する。
2. 目標画像と同じサイズの真っ黒の画像を用意する。
3. パラメータとして、[円の位置(xy座標), 円の半径] を円の数だけ保持する。（ランダムに初期化する）
4. パラメータを元に、真っ黒の画像に白い円を描画する(生成画像)。
5. 生成画像の「離れると読める」度の評価値を計算する。
6. 評価値を最適化するようなパラメータを求める。

1,2,3,4は問題なさそうです。ちなみに4については、各円の座標・半径から円形のmask画像を生成し, 全ての円に関して積を取ることで実現しています。以下のような感じです。（都合上、白地に黒い円を描く感じになっています）

def circle_mask(shape, cx, cy, r):
    # cx, cyは円の中心の座標, rは半径
    x, y = np.ogrid[:shape[0], :shape[1]]
    r2 = (x-cx)*(x-cx) + (y-cy)*(y-cy) # 各点 に対して(cx, cy)からの距離 ** 2　を計算
    mask = r2 > r ** 2 # 半径 ** 2 以下なら円内
    return mask

def encode(shape, params):
    # paramsは[x座標, y座標, 半径] * 円の数 であるnumpy 2d-array, shapeは生成画像のサイズ
    _circle_mask = lambda arr: circle_mask(shape, *arr) # 部分適用的な
    mask = np.multiply.reduce(np.apply_along_axis(_circle_mask, 1, params))
    return np.ones(shape) * mask

問題は5と6です。以下詳細です。

「離れたら読める」 はどう評価するか

そもそもなぜあの広告は「離れると読める」 のか?

基本的にはあの広告は有名な錯視「ルビンの壺」と同じ原理です。人は白と黒（白と黒でなくてもいいですが）からなる画像を見ると、どちらかを「図」（意味的なまとまり）どちらかを「地」（背景） と認識します。

「図」と「地」というのはゲシュタルト心理学(wikipedia)における用語で、デザイン分野なんかでは馴染みのある言葉のようです。
問題の広告においては、近くで見ると円が「図」になり、離れて見ると文字が「図」になります。
これは、人間が「まとまりがある方」を「図」として認識する、という法則（プレグナンツの法則）によります。近くでみたとき、円(白)は円というまとまりで認識できる一方で文字(黒)はガタガタしていて、まとまりがありません。
また、一度円が目に入ってしまうと、文字の一部しか目に入らないことになるのでなかなか文字は認識できません。

しかし、離れてみると画像はぼやけて解像度が下がります。ぼやけることによって、文字(黒)のガタガタは平滑化され、気にならなくなります。こうなると、文字を「まとまり」として認識するため、文字が見えるという仕組みです。

なお、目のいい人なら「ちょっと離れたくらいではぼやけて見えてなんかいないぞ」と感じるかもしれません。
人間の目においてパシッとピントを合わせられるのはごく狭い範囲なのですが、視野の中心のピントが合っていたら全体として「ぼやけていない」ように見えます。でも、実際のところ視野の中心以外はぼやけています。今読んでいるこの文も、今まさに読んでいる数文字以外はぼやけています。試してみてください。

「ぼけ」 の表現

人間の目（というかレンズ一般）の「ぼけ」はガウシアンフィルタの畳み込みで表現できます。そこで、目標画像と生成画像にフィルタをかけた後の両者の距離が十分に小さければ「離れたら読める」 画像であると言えると考えました。つまり、以下のような手順になります。

1. 目標画像・生成画像にガウシアンフィルタをかけ、ダウンサンプリングする³。
2. ダウンサンプル後の2つの画像のL2ノルムを距離=評価値とする。

具体的には評価関数は概ね次のように実装しました。

    def evaluate(target_image, generated_image):
        ds_rate = min(target_image.shape) // 20 # 画像の短辺が20pixになるようなダウンサンプルレート
        down_sampled_generated = scipy.ndimage.gaussian_filter(generated_image, sigma=ds_rate)[::ds_rate, ::ds_rate]
        down_sampled_target = scipy.ndimage.gaussian_filter(target_image, sigma=ds_rate)[::ds_rate, ::ds_rate]
        return np.linalg.norm(down_sampled_morph - down_sampled_target)

パラメータをどうやって最適化するか

はじめに考えたのは、最適化問題として解くことです。しかし、パラメータ->評価値の関数をちゃんと定義して適切に解くほどの数学力は僕にはありませんでした...

次に目をつけたのは「図形詰め込みアルゴリズム(参考pdf)」という最適化アルゴリズムです。与えられた空間内に与えられた図形をできるだけみっちり詰め込むというタスクで、これに落とし込むことができれば解けそうだと考えました。
一方でこのアルゴリズムは詰め込む図形同士が「接している」ということを重要な手がかりとして使っているのですが、今回配置したい円は離れていてもよく、重なっていても構いません。
その道の専門家なら難なく条件を変化させて解いてしまうのかもしれませんが、僕にはこれも難しそうだったので諦めました。

最後に、遺伝的アルゴリズムを考えました。円の配置とサイズを遺伝情報として持っておいて、その評価値さえ得られればいい⁴ので、これならできそうに思えました。遺伝的アルゴリズムには以前から興味を持っていましたが実装したことはなかったので、0から実装してみることにしました。

Pythonで0からクラス設計を行ったのも初めてなので、ご指導ご鞭撻のほどよろしくお願いいたします。

遺伝的アルゴリズムの実装

遺伝的アルゴリズム略解

遺伝的アルゴリズムについてはこのスライドが非常にわかりやすいです。
ここでは必要最低限の解説をします。（間違いの指摘をお願いいたします）

まず、遺伝的アルゴリズムには以下の構成要素があります。

• 獲得したいパラメータを格納するゲノム
• ゲノムをもつ個体
• 個体の集合である世代

これらに対して以下の操作を行います。

1. 得体パラメータをゲノムとしてエンコードし、個体に割り当てる。
2. 一定数の個体の集合を現世代とする。
3. 現世代の各個体の評価値を計算する。
4. 世代内の2個体をある確率で選び、交差を行う。
   
   交差: 2個体の持つゲノムの一部を入れ替えることによって新たに2つのゲノム(=子孫)を生成する
5. 交差によって生成した子孫にある確率で突然変異を起こす。
   
   突然変異: ゲノムの一部に変化を加える。
6. 交差を繰り返すことによって子孫を増やす。
7. 子孫と現世代の中から、次世代に残すゲノムを選択する。
8. 現世代を次世代で置換して(世代交代)以上の操作を繰り返す。

交差の親個体の選び方や選択の仕方は色々あって、それらの確率を評価値に依存させることでより良いゲノムを残していきます。今回、はじめに採用した戦略は以下の通りです。（その後変更した点もあります。結果で述べます）

• エンコーディング: 実数値エンコーディング
  • ゲノムは1d-arrayとして表現するのが普通なのかもしれませんが、今回は[[円のx座標, 円のy座標, 円の半径] * 円の数]という実数値2d-arrayをゲノムとしました。
• 評価値: 前述
• 交差方法: 二点交差
• 突然変異: 置換（ランダムに選ばれた遺伝子を、一様乱数によって書き換える）のみ
• 世代交代: 世代間最小ギャップモデル(後述)

世代間最小ギャップモデル

今回実装したのは「世代間最小ギャップモデル」というモデルです。次の手順で世代交代を行います。

1. ある世代からランダムに2個体を選んで交差し、子個体2体を得る（合計4個体）
2. 合計4個体のうち、最も評価値の良い1個体を生存個体_1として選ぶ。
3. 残った3個体のうちから、ルーレット選択によって確率的に1個体を選び、生存個体_2とする。
4. 2体の生存個体によって親個体を置換したもの（ただし、親自身も生存個体となりうる）を次世代とする。

すなわち、一回の世代交代で入れ替わる個体は最大2個体です。親個体の選択がランダムである（評価値によらない） ため、世代内の多様性の維持に優れており、計算コストも低く抑えられます。

Class構成

詳細はgithubにあります。
以下の4つのクラスによってプログラムを構成しました。Phenotype(表現型)クラスにゲノムを渡すと形質が発現する、という設計にするとかっこいいなと思って作ったのですが、あまり使い勝手のいい出来にはなりませんでした...。

• Genom(genom_length:int, genom_limits:array-like)
  • [円のx座標, 円のy座標, 円の半径] * genom_length に相当するゲノム(2d-array)を格納することを主な目的としたクラス。

• Phenotype(g:Genom, shape:array-like)
  • genomの発現を管理するクラス。コンストラクタはGenomと生成画像のshapeを引数にとる。genom内の円の位置・サイズ情報を元に、画像にエンコードする。
  • evaluate(target)メソッドによってgenomの評価値を返す。targetはtarget文字画像。

• Family(g1:Genom, g2:Genom)
  • 交差を司るクラス。コンストラクタは親となる2個体のGenomを引数にとる。
  • mgg_change(target, pm)メソッドによって世代間最小ギャップモデルに基づく世代交代における生存個体（2個体）を返す。

• Generation(args)
  • genomの集合である世代を管理するクラス。コンストラクタに世代サイズや突然変異確率などの設定を渡す。
  • mgg_change(target)メソッドで世代間最小ギャップモデルに基づく世代交代を行う。
    • 基本的にこのメソッドを繰り返すことで遺伝的アルゴリズムが実行される。

main関数はシンプルに書けて、以下のような感じです。ただしtargetは目標画像(2d-numpy.ndarray)です。

g = Generation(generation_size=GENERATION_SIZE, genom_length=circle_num, genom_limits=[height, width, min(height, width) // 2], pm=0.05)
for i in range(ITER_NUM):
    g.evaluate(255-target) # 本当は必要ない。summaryを得るには必要なので、実際のコード中では1000回に一回くらい評価してsummaryを出力している。
    min_score, max_score, ave_score = g.summary()
    if mi < THRESHOLD:
        print("score achieved, break..")
        break
    g.mgg_change(255-target)

結果と改良

さて、ここまでの内容を意気揚々と書き上げ、あとはちょっとチューニングしていい結果が出たら投稿だ〜、とたかをくくっていたのですが。
ここからが沼でした...

最初の結果

上記の設定で遺伝的アルゴリズムを回した結果です。(50000回交配)

感想としては、「おっ、ええやん」という感じです。円がだいたいいい感じに配置されていますし、ちょっと改良したら完成かな？　と思いました。
改良すべき点は以下です。

• 円の数が多くてごちゃごちゃしている。
• 円が重なりすぎている。例えば"g"の楕円部分に二つ以上の円が重なっていて、文字の再現度は高いが「円っぽさ」が失われている。

評価関数の変更

とにかく、「円の数」と「円の重なり」が問題だと思ったので二つの改良を加えました。

欠失・挿入の導入

円の数は固定で与えていましたが、これを可変にすることにしました。円の数はgenomの長さに対応します。

待てよ？　genomの長さは遺伝情報じゃないじゃん...

円の数を学習(という用語がGAにおいて正しければ)したいならgenomのどこかに円の数を書いておくべきだったと思いました。
例えばgenomの先頭に、発現させる遺伝子の位置を書いておくとか。

しかし既にgenomは2d-arrayとして実装してしまっているので、構造を抜本的に帰るのはめんどくさいコストが高すぎると考えて別の手法を考えました。

突然変異の一つである欠失・挿入の採用です。欠失はゲノム⁵の一部を削除すること（ゲノムの長さが減る）、挿入はゲノムの一部に遺伝子を追加することを言います。
これによってgenomの長さ=円の数を可変にできます。

円の重なりを評価

円の数がが少ない個体を選択したいという意味でも、出来るだけ円同士を離れさせたいという意味でも、円の重なりは小さいほうがいいと考えました。
そこで、評価関数を次のように書き換えます。依然として、小さい方がいい評価値とみなします。

evaluation = L2|Gauss\_filter(downsample(目標画像 - 生成画像))| + \mu * \max(threshold, overlap)

overlapが円の重なりを表しています。
μは円の重なりの寄与の大きさを与えるパラメータです。thresholdは、それ以下の重なりは評価しないということです。全ての円が完全に離れている必要はなく、重要なのは依然として目標画像と生成画像の距離だからです。
画像の距離の評価値(前項)とoverlap(後項)のオーダーを概算すると後者の方が一桁くらい大きかったので、前項を重視したいという気持ちでμは0.01程度に設定しました。

overlapは次のように計算しています。

class Phenotype(...)
    ...
    def overlap(self):
        total_circle_area = np.sum(self.genom[:, 2] ** 2 * np.pi)
        overlap = total_circle_area - np.sum(self.morph == 0)
        return overlap

self.morph⁶というのが生成画像で、その中の円が占める面積と、genomから得られる円の延べ面積の差を取っています。

この方法にはちょっと問題点があります。

• そもそも描画した円の面積は量子化されている(デジタルである)ので、半径くらいのオーダーで誤差が乗る
• 画像の外に円の一部がはみ出している場合、その面積がoverlapに加算されてしまう

しかし、計算が簡単なのでこれで押し切ります。上記の問題点はthresholdで吸収することを期待します。

以上の改良の結果がこれです。(ちょっと画質がいい)

...悪化してへんか????
確かにある程度円同士が離れた感じはありますが、円の数は相変わらず多いままで、総じてごちゃごちゃ感が増しています。
実行中のログを見ると、挿入・欠失によって円の数が変わることはあるものの、デフォルトの30からせいぜい28くらいまでにしか変わっていないことがわかりました。

うむ......

円の数をランダムに初期化

次にこう考えました。

.oO( 遺伝的アルゴリズムがうまくいくには多様性が必要や。円の数は初めから固定されてるねんから突然変異だけでうまく進化できるわけがない... )

そこで、円の数=genomの長さをランダムに初期化することにしてみました。具体的には30を平均とする正規分布に基づいて初期化します。

class Genom:
    def __init__(self, genom_length, ...):
        # genom_length * 0.5 ~ genom_length * 1.5の間に存在する確率がだいたい95%
        self.genom_length = np.ceil(np.random.normal(genom_length, scale=(0.25 * genom_length)))

これでどうだ。(世代の個体数を30から50に増やしました)

今度はちょっと円が足りない。特に"a"の右下あたり。というかどんどん悪化している気がする。あと"g"の下の方に小さい円が並んでいてかわいい。

初めの1000世代交代くらいで一気に円の数が少ない方に収束してしまいました。
よく考えると当たり前で、円がランダムに配置されている状態であれば円の数が少ない方がずっとoverlapが少ないし、目標画像との距離も（白黒の画素数比的に）小さくなりがちだと考えられるので、易きに流れる形でいわゆる初期収束に陥ったと考えられます。

大変異の導入

初期収束に陥ると多様性がなくなってしまって進化が起きなくなります。
多様性の消失に対して、大変異という手法があります。
これは自然界での急激な環境の変化などをモデル化したもので、変異確率を一定の周期（長い周期）でぐんと引き上げます。

ただ、今回採用している世代間最小ギャップモデルに対して大変異を施すことが妥当かどうかは自信がありません。(というか、妥当ではないと確信している)

大変異の目的は、生存個体のトレンドを大きく変えることです。例えば、世代交代ごとに世代の数10%が子孫に置換されるようなモデルにおいては、大変異によって、その子孫のトレンドは大きく変わることが期待されます。
一方、世代間最小ギャップモデルでは、世代交代で置換される個体は多くても2個体ですから、大変異が全体のトレンドに大きな影響を与えるとは思えません。

とはいえ、藁にもすがる思いでやってみました。なんにせよ、一定の周期で変異確率が上がれば多様性に寄与しなくはなさそうなので。

それから、突然変異に摂動を追加しました。
これは置換と似ていますが、置換がゲノムの一部の遺伝子を別のものに書き換えるのに対し、摂動は遺伝子の値を少しだけ変化させます。

実装としては、平均を0とする正規分布からサンプルした値を加算しています。

class Genom(...):
    ...
    def perturbate(self, pm):
        # 変異させる遺伝子を選ぶ。各遺伝子(genomの各行)は確率pmで変異する
        whichwillbechanged = np.random.choice([True, False], self.chromosome.shape, p=[pm, 1-pm])
        self.chromosome[whichwillbechanged[:,0], 0] += np.floor(np.random.normal(0, scale=10, size=np.sum(whichwillbechanged[:,0]))).astype(np.int)
        self.chromosome[whichwillbechanged[:,1], 1] += np.floor(np.random.normal(0, scale=10, size=np.sum(whichwillbechanged[:,1]))).astype(np.int)
        self.chromosome[whichwillbechanged[:,2], 2] += np.floor(np.random.normal(0, scale=10, size=np.sum(whichwillbechanged[:,2]))).astype(np.int)
        self.chromosome[self.chromosome < 0] = 0 # 負にならないように

なお、コードブロック中のコメントに書いていますが、僕は変異率について、各遺伝子が確率Pmで変異すると解釈して上のコードを書いています。(この場合、各ゲノムが変異する可能性はPmよりはるかに大きくなります)
これは誤りかもしれなくて、各ゲノムが確率Pmで変異し、変異することが決まったゲノムのうちで変異する遺伝子は別に決めるのかもしれないです。
以下のような感じで。

class Genom(...):
    ...
    def perturbate(self, pm):
        # 確率pmで変異する。
        willchange = np.random.choice([True, False], size=None, p=[pm, 1-pm])
        if willchange:
            change_ind = np.random.randint(self.length, size=3)
            self.chromosome[change_ind, :] += np.floor(np.random.normal(0, scale=10, size=3)).astype(np.int)
            self.chromosome[self.chromosome < 0] = 0 # 負にならないように

ちゃんと実装したい/知りたい方は専門家の実装とかに当たってください。

円の重なりの評価方法の変更

瑣末ですが、前述の円の評価方法の問題点がバカにならないと言う事実に直面しました。
特に、半径が大きいほど重なりが大きく評価されて、結果として円のサイズが全体で小さくなってしまっています。（さっきの結果を見てみてください）

そこで、面積ではなく距離ベースで評価し、さらに半径で標準化すると言う評価方法に切り替えました。
距離の総当たりの計算はPythonでxy座標上の2点間の距離をforループを使わずに計算する方法を参考にしましたが、すっきり書けて嬉しいです。

class Phenotype(...)
    ...
    def overlap2(self):
        # 面積ベースではなく、距離ベースで重なりを評価する。
        locations = self.genom[:, :2]
        rs = self.genom[:, 2]
        # 円同士の距離
        distances = np.sqrt(np.sum((np.expand_dims(locations, axis=0) - np.expand_dims(locations, axis=1)) ** 2, axis=2))
        # 円同士の半径
        arms = np.expand_dims(rs, axis=0) + np.expand_dims(rs, axis=1) - 2 * np.diag(rs) # 対各成分を引いている
        interfer = arms - distances
        # 各円の半径に対するinterferの値。標準化的なこと。円が大きいとペナルティが増大するのを避けたい。
        interfer = interfer / (rs + 0.1) # 半径0だとまずい
        # 十分に離れている場合は0
        interfer[interfer < 0] = 0
        return np.sum(interfer)

以上の結果（+細々したパラメータのチューニング）の結果がこれです。これで許してください！！！

展望

• 「離れると文字に見える」 と言うところで、目標画像とのL2距離を使ったが、文字のトポロジーなどを考慮して評価できると良さそう
• 円の数に関して、初めにランダマイズするよりも収束してきたときに減らしたり増やしたりした方が多様性に寄与しそう
• 上記と関連して、途中で評価関数を変えるなど、いろんなことが考えられる
• 世代間最小ギャップモデルは分散処理と相性がいいようなので、高速化してみたい
• 色々書いたけど多分やらない

まとめ

遺伝的アルゴリズム、偶然と選択のみによってまあまあいい結果が得られて気持ち悪い！！！

あと、チューニングすべきところがいっぱいあって大変（手法の選択、各確率etc...）

学ぶことはたくさんあって（ここには書いてないですが、処理があまりに遅かったのでcProfileを使ってボトルネックの処理を改善したり、numbaを導入してあんまり早くならなかったりして楽しかったです）よかった。

駄文に長々と付き合ってくださってありがとうございました！

追加情報

Pythonには遺伝的アルゴリズム用のライブラリがちゃんとあります。
Deapというのが有名なようです。
解説記事もあります。

参考文献

遺伝的アルゴリズム（Genetic Algorithm）を始めよう！

• https://www.slideshare.net/kzokm/genetic-algorithm-41617242
• このスライドは本当にわかりやすくて素晴らしかったです。
  Pythonで遺伝的アルゴリズム
• https://qiita.com/KeisukeToyota/items/0f527a72270430017d8d
  python - numpy配列の円形セクターをマスクする
• https://ja.coder.work/so/python/503264
  --GAノート--
• http://structure.cande.iwate-u.ac.jp/yamamoto/ga_note/ga_note.html#6

1. 大本のアイデアは研究室の先生から降ってきました。 ↩

2. この図は僕がマイクロソフト・パワーポイントを使って手動生成しました。 ↩

3. 畳み込みは線形操作なので、畳み込みの差=差の畳み込みです。実際のコードでは畳み込みの処理回数を減らすために、差を取ってからたたみ込んでいます。
   ガウシアンフィルタはscipyにscipy.ndimage.gaussian_filter()が備わっていたので利用しました。画像の距離がL2ノルムというのは乱暴かもしれませんが、簡単のために採用しました。
   mnistなどで事前学習したCNNによる特徴量のcos距離...とかにするとさらに本格的だとは思います。
   ↩

4. 正確な言い方は、「遺伝情報の探索空間が位相を持ち、評価値に全順序を定めることができる」だそうです。 ↩

5. 紛らわしいですが、一応genom=僕が実装した円の位置に関する遺伝情報、ゲノム=遺伝的アルゴリズム全般について述べるときの遺伝情報、と使い分けています。 ↩

6. 発現形態のアナロジーです。 ↩