分布推定アルゴリズムを学ぶ (1) ～概要～
分布推定アルゴリズムを学ぶ (2) ～PBIL～
分布推定アルゴリズムを学ぶ (3) ～UMDA～
分布推定アルゴリズムを学ぶ (4) ～CGA～の続き

Extended Compact Genetic Algorithm(ECGA)とは

目的関数の変数をクラスタリングすることで，変数間に依存関係を持つ問題に対処した離散ブラックボックス関数最適化手法．
これまでに紹介したPBILやCGAでは問題対象がビット列なのに対し，ECGAは離散空間全体が対象(もちろんビット列も対象)．

概要

ECGAは，確率モデルとして周辺分布モデル(Marginal Product Model; MPM)を用い，変数間の依存関係を考慮しながら最適化を行う．
ECGAでは，周辺分布モデルを構築するために，評価指標として最小記述長(Minimum Description Length; MDL)を用いた貪欲探索を行う．

以下では，ECGAを構成する周辺分布モデルと最小記述長について説明し，貪欲探索を含むECGAのアルゴリズムについてまとめる．

周辺分布モデル(Marginal Product Model; MPM)

以下に，4次元ビット列の場合の周辺分布モデルの例を示す．

この例では，2，3次元目の変数は独立であり，1，4次元目の変数が1つの変数としてまとめられている．
このように，任意個の変数を1つの変数として扱うモデルを周辺分布モデルと言う．

上記の[1, 4]の確率分布の場合，00または11の値しかサンプリングされない．
この確率分布は集団(左図)を用いて推定しているが，01と10の確率が0であることから，集団内の個体は全て(0,x,x,0)または(1,x,x,1)の形になっているということが表現されている(xは0または1を取る)．
逆に言えば，集団内には(0,x,x,1)や(1,x,x,0)といった個体は存在しないということである．

仮に，各変数を独立に扱い，確率分布を推定したとすると上記のような表現ができない．
以下に，各変数を独立に扱い，確率を推定した結果を示す．

各変数が独立の場合，[1]と[4]はそれぞれ0または1の値を取りうるため，(0,x,x,1)や(1,x,x,0)といった個体が生成される可能性があり，集団内の個体は(0,x,x,0)または(1,x,x,1)だけであるという特徴を表現できていない．
1，4次元目の値は00または11で良いということが分かっているのに，01や10の個体も生成してしまうということは，それだけ非効率的な探索になってしまうということである．

以上より，各変数を独立に扱うより，複数の変数を適切に1つの変数として扱った方が，与えられた集団にフィットする確率モデルを構築でき，効率的な探索に繋がることが分かる．

次の問題は，どの変数を1つの変数として扱うべきかという問題である．
上記の例では，[1,4]という情報を事前に与えていたが，ブラックボックス関数ではどの変数間に依存関係があるかは不明である．
そのため，変数間の依存関係の度合いを定量化できる指標が必要となる．
これが，次に説明する最小記述長である．

これ以降，複数の変数を1つの変数としてまとめることをクラスタリング，生成された変数をクラスタと呼ぶ．

最小記述長(Minimum Description Length; MDL)

周辺分布モデルの説明では，クラスタリングを行うことで集団によりフィットした確率モデルが構築でき，効率的な探索が実現できると述べた．
それならば全ての変数を1つのクラスタとしてまとめることが，集団にフィットさせるという観点では最良だと考えられる．
しかし，これは(現実的な最適化の文脈では)誤りであり，全ての変数を1つのクラスタとして扱うとモデルが非常に複雑になり，与えられた集団のみに適応した状態，つまり過学習を引き起こす．
最適化途中で与えられる集団は常に問題構造を表現した正しいものであるとは限らないため，その集団のみに適応してしまうことは危険である．
そのため，集団にフィットしつつも，ある程度の汎化性能を持った周辺分布モデルになるようなクラスタリングが求められる．

これを実現するため，最小記述長$C$は，集団へのフィット具合を表す項$C_p$とモデルの複雑さに対するペナルティー項$C_m$から構成される．

C = C_p + C_m, \\
C_p = \lambda\sum_{i} E(M_i), \\
C_m = \log\lambda \sum_{i}2^{S[i]},

ここで，$\lambda$は集団サイズ，$M_i$は$i$番目のクラスタの確率分布，$E(M_i)$は$M_i$に対するエントロピー，$S[i]$は$i$番目のクラスタのサイズ(上記の周辺分布モデルの例では，[1,4]はサイズ2，[2]はサイズ1)を表す．
1つのクラスタのサイズが大きくなると，特定の同時確率の値が大きくなるため，$C_p$は小さくなる(=集団にフィット)．
一方で，$C_m$はクラスタのサイズに指数的に比例して大きくなる(=モデルの複雑さに対するペナルティー)．
そのため，$C$を最小化するように変数のクラスタリングを行うことで，汎化性能を保ったまま集団にフィットした周辺分布モデルを構築することができる．

アルゴリズム

Algorithm 4にECGAの全体の流れを示す．
流れは，EDAsの概要で示したアルゴリズムと同じである．

5行目の確率モデルの構築時，Algorithm 5に示す貪欲探索によって周辺分布モデルを構築する．

実装

プログラム全体はgithubに置いてある．

環境

Ubuntu18.04
Python 3.6.5
anaconda Command line client (version 1.7.2)
conda 4.8.3
numpy 1.18.1

githubの方にSingularityの環境ファイルを置いてあるので，それをビルドしても動く(動作確認済み 2019/12/19時点)．

ECGAの実装を以下に示す．
基底クラスのEDABaseクラスに関しては分布推定アルゴリズムを学ぶ (2) ～PBIL～を参照

import numpy as np

from eda.optimizer.eda_base import EDABase
from eda.optimizer.util import SubSet, Cache


class ECGA(EDABase):
    def __init__(self, categories, replacement,
                 selection=None, lam=500, theta_init=None):
        super(ECGA, self).__init__(categories, lam=lam, theta_init=theta_init)
        self.replacement = replacement
        self.selection = selection

        self.population = None
        self.fitness = None
        self.cluster = None

    def update(self, x, evals, range_restriction=False):
        x, evals = self._preprocess(x, evals)
        if self.selection is not None:
            x, evals = self.selection(x, evals)
        x = np.argmax(x, axis=2)
        if self.population is None:
            self.population = x
            self.fitness = evals
            self.lam = int(self.lam * self.replacement.replace_rate)
        else:
            self.population, self.fitness = self.replacement(self.population,
                                                             self.fitness,
                                                             x,
                                                             evals)
        self.cluster = self.greedy_mpm_search(self.population)

    def greedy_mpm_search(self, population):
        # initialize subset of cluster
        cluster = [SubSet(i, population[:, i], self.Cmax) for i in range(self.d)]
        # initialize cache
        cache = self.initialize_mpm(cluster)
        # clustering according to CCO
        while True:
            pos_i, pos_j = cache.argmax_cc()
            if cache.cc_list[pos_i, pos_j] <= 0:
                break
            cluster[pos_i] = cache.subsets[pos_i, pos_j]
            cluster.pop(pos_j)
            cache.remove(pos_j)
            for k in range(len(cluster)):
                if k == pos_i:
                    continue
                i, k = (k, pos_i) if k < pos_i else (pos_i, k)
                merge = cluster[i].merge(cluster[k])
                cache.add(i, k, cluster[i].cc + cluster[k].cc - merge.cc, merge)
        return cluster

    def initialize_mpm(self, cluster):
        cache = Cache(self.d)
        for i in range(len(cluster)-1):
            for j in range(i+1, len(cluster)):
                subset1 = cluster[i]
                subset2 = cluster[j]
                merge = subset1.merge(subset2)
                cache.add(i, j, subset1.cc + subset2.cc - merge.cc, merge)
        return cache

    def sampling(self):
        # random sampling, only first generation
        if self.cluster is None:
            rand = np.random.rand(self.d, 1)
            cum_theta = self.theta.cumsum(axis=1)

            c = (cum_theta - self.theta <= rand) & (rand < cum_theta)
            return c
        # sample by using each probability of cluster that clustered by CCO
        else:
            c = np.zeros((self.d, self.Cmax), dtype=bool)
            for cl in self.cluster:
                rand = np.random.rand()
                cum_theta = cl.theta.cumsum().reshape(cl.theta.shape)
                _c = (cum_theta - cl.theta <= rand) & (rand < cum_theta)
                if len(cl) > 1:
                    _c = np.unravel_index(np.argmax(_c), _c.shape)
                c[cl.idx_set, _c] = True
            return c

    def convergence(self):
        if self.cluster is None:
            return 0.5
        return np.mean([np.max(c.theta) for c in self.cluster])

実験：ECGAの性能評価

設定

目的関数はOne-Max問題，3-Deceptive問題を用いる．
One-Max問題は$D=100$，3-Deceptive問題は$D=30$とする．
また，3-Deceptive問題のユーザーパラメータは$d=0.1$とする．

その他のユーザーパラメータは以下の通り．

パラメータ名	値
集団サイズ$\lambda$	{100, 300, 500}
選択手法	トーナメント選択
選択割合	集団サイズの50%
置換方法	制限付きトーナメント置換
置換割合	集団サイズの50%

制限付きトーナメント置換(Restricted Tournament Replacement; RTR)は，論文[2]においてECGAに適用されている集団内の多様性を保つ工夫が導入された置換方法である．

以下のいずれかを満たした場合，最適化を終了する．

大域的最適解をサンプルできたとき(成功)
確率分布が収束したとき(失敗)
上限世代数($2\times 10^3$)に達したとき(失敗)

以上の設定で，乱数のシードを変更して10試行行う．

観察指標

以下の指標を記録する．

成功試行数(大域的最適解を見つけられた回数)
平均評価回数($\pm$標準偏差)
最適化終了時の最良評価値の平均($\pm$標準偏差)

ここで平均評価回数は，成功試行時のみのものとする．

結果

One-Max問題($D=100$)をECGAで最適化したときの結果を以下の表に示す．

集団サイズ	成功試行数	平均評価回数($\pm$標準偏差)	平均最良評価値($\pm$標準偏差)
100	8/10	2600$\pm$206.16	99.8$\pm$0.4
300	10/10	7560$\pm$640.62	100$\pm$0.0
500	10/10	12350$\pm$550.00	100$\pm$0.0

3-Deceptive問題($D=30$)をECGAで最適化したときの結果を以下の表に示す．

集団サイズ	成功試行数	平均評価回数($\pm$標準偏差)	平均最良評価値($\pm$標準偏差)
100	0/10	-	9.56$\pm$0.17
300	7/10	4242.86$\pm$493.39	9.95$\pm$0.08
500	10/10	6200.00$\pm$1288.41	10.0$\pm$0.0

考察

One-Max問題では，集団サイズ300と500で全試行で最適解を発見できている．
評価回数の観点では集団サイズが小さい方が優れている．
これは，One-Max問題では各変数が独立であり，ECGAのクラスタリングに意味が無いからである．
一般に正確なクラスタリングのためには大きな集団サイズが必要だが，クラスタリング自体が無意味なので小さな集団サイズでも解けてしまう．

3-Deceptive問題では，集団サイズが500のときに全試行で最適解を発見できている．
集団サイズが小さいと最適解を見つけられない原因としては，変数間のクラスタリングがうまく行えていないことが挙げられる．
クラスタリングの際，MDLを用いて変数間に依存関係が存在するかを判定するが，集団サイズが小さいと依存関係によって生まれるバイアス(ある変数とある変数が同じ値になりやすいなど)が現れにくくなり，誤った依存関係を検知してしまう．
その結果，問題構造をうまく表現できず局所解に収束してしまう．

以下に，3-Deceptive問題を集団サイズ500で最適化したときのクラスタリング結果の一例を示す．

世代	クラスタリング結果
1世代目	[0,9], [1,2], [3,5], [4,13], [6,15], [7,18], [8,20], [10,24], [11,21], [12,27], [14,28], [16,17], [19,23], [22,25], [26,29]
5世代目	[0,2,1], [3,5], [4,13], [6,14], [7,18], [8,12], [9,11,10], [15,17,16], [19,20], [21,23,22], [24,25,26], [27,28,29]
10世代目	[0,1,2], [3,5,4], [6,8,7], [9,10,11], [12,14,13], [15,17,16], [18,19,20], [21,23,22], [24,25,26], [27,28,29]
13世代目	[0,1,2], [3,5,4], [6,8,7], [9,11,10], [12,14,13], [15,17,16], [18,19,20], [21,22,23], [24,25,26], [27,28,29]

クラスタリング結果より，世代を追うごとに3ビットずつの変数間依存性を正しくクラスタリングできている様子が観察できる．
この試行では，13世代目に最適解を発見しているが，おおよそ10世代目で問題構造に対する最適化が終了し，13世代目まででパラメータの最適化が行われたと考えられる．

まとめ

ECGAは問題内の変数をクラスタリングすることで，変数間に依存関係を持った問題に対処
変数間依存性のある問題では，少ない反復数で最適解を見つけられる(ただし，適切な依存関係を検知するために集団サイズを大きく設定する必要がある)
各変数が独立な問題に対しての性能は微妙(これまでに紹介した手法と比べて)

参考文献

Linkage Learning via Probabilistic Modeling in the ECGA
Investigating Restricted Tournament Replacement in ECGA for Non-Stationary Environments
遺伝的アルゴリズム -その理論と先端的手法-