この記事について
多目的最適化入門シリーズ④で「これで一通り説明できたかと思います」と締めたのですが、後から読み返すと思いっきり説明を飛ばしていた箇所がありました。
②のNSGA-IIの実装のここの箇所や
# 子集団の生成
offspring = generate_offspring(population, objectives)
③のMOEA/Dの実装のここの箇所です。
# 近傍からランダムに2つ選んで交叉
j, k = np.random.choice(neighborhoods[i], 2, replace=False)
y = crossover(population[j], population[k])
y = mutation(y)
generate_offspring、crossover、mutation……この辺が全部ブラックボックスのままでした。。。
NSGA-IIもMOEA/Dも「集団をどう評価して、どれを残すか」が主役のアルゴリズムなので、記事ではそちらを厚く書いたのですが、いざ動かすとなるとこの中身が結果を大きく左右します。実際、後半で実験しますが同じNSGA-IIでも交叉を変えるだけでIGDが275倍悪化しました（0.0050 → 1.3761）。実際に論文レベルでもこの辺を詳しく調査しているも野もあるので、さすがにこれを説明しないのはまずいなと。
というわけでシリーズを2記事だけ延長して、 交叉（Crossover） と 突然変異（Mutation） の中身を説明します。今回は交叉編ということで、実数値GAの定番である SBX（Simulated Binary Crossover） を扱います。
多目的最適化入門シリーズの他の記事はこちら。
| 記事 | 内容 |
|---|---|
| ① | 多目的最適化とは？基礎概念 |
| ② | NSGA-IIの説明 |
| ③ | MOEA/Dの説明 |
| ④ | ベンチマーク問題の説明 |
| ⑤（本記事） | 交叉（SBX交叉）の説明 |
| ⑥ | 突然変異と選択の説明 |
交叉・突然変異はアルゴリズムのどこにいるのか
改めて、②で示した進化型アルゴリズムの基本フローを再掲します。
① 初期集団の生成
② 各個体の評価（目的関数値の計算）
③ 選択（良い個体を親として選ぶ） ← ⑥で説明
④ 交叉・突然変異（新しい個体を生成） ← 本記事と⑥で説明
⑤ 次世代集団の形成 ← ②③で説明済み（非劣解ソート・スカラー化）
⑥ ②に戻って繰り返す（収束するまで）
②のNSGA-IIも③のMOEA/Dも、説明していたのはステップ⑤（次世代集団の形成）だけで、ステップ③の選択とステップ④の交叉・突然変異はまるごと飛ばしていたわけです。
3つのオペレータの役割分担は以下のようになっています。
| オペレータ | やること | 探索と活用でいうと |
|---|---|---|
| 選択（Selection） | 良い個体を親に選ぶ | 探索を良い領域に集中させる |
| 交叉（Crossover） | 2つの親を組み合わせて子を作る | 活用（exploitation）：今ある良い情報の組換え |
| 突然変異（Mutation） | 一部をランダムに変える | 探索（exploration）：情報がないところへ飛ぶ |
交叉は「すでに集団が持っている情報を組み替える」操作なので、集団の中にない情報は基本的に生み出せません。これが後ほど効いてきます。（ちなみに厳密には交叉でも集団にない情報を生み出すような交叉を用いることが多いですが、概念・考え方的には集団の中にない情報は基本的に生み出さないと考えていただいて遜色ないです。）
解の表現が変われば演算子も変わる
まず押さえておきたいのが、交叉は「解をどう表現するか」とセットで決まるということです。
| 解の表現 | 例 | 代表的な交叉 | 代表的な突然変異 |
|---|---|---|---|
| バイナリ列 | 10110100 |
1点交叉・一様交叉 | ビット反転 |
| 実数値ベクトル | $(0.31, 0.85, \ldots)$ | SBX・ブレンド交叉・DE | 多項式突然変異・ガウス変異 |
| 順列 | [0, 3, 1, 4, 2] |
OX・PMX | スワップ・2-opt |
多目的最適化のベンチマーク問題（④で説明したZDT・DTLZ・WFG）はすべて実数値ベクトルの問題です。決定変数が $\mathbf{x} \in [0,1]^{30}$ みたいな形をしているので、実数値向けの交叉が必要になります。
そういうわけで、この記事では実数値ベクトル向けの交叉を扱います。
バイナリの1点交叉を実数値にそのまま持ってくると？
なぜわざわざ「SBX」という専用の交叉が必要なのか。素朴なやり方でどう失敗するかを見るのが良さそうなので、それを見ていきましょう。
そもそも1点交叉とは
GAの原点はバイナリ列です。切断点を1つ決めて、そこから後ろを入れ替えます。
切断点：3ビット目と4ビット目の間
親1: 1 0 1 | 1 0 1 0 → 10進で 90
親2: 0 1 1 | 0 1 1 1 → 10進で 55
子1: 1 0 1 | 0 1 1 1 → 10進で 87 ← 親1の前半 + 親2の後半
子2: 0 1 1 | 1 0 1 0 → 10進で 58 ← 親2の前半 + 親1の後半
ここで、この1点交叉が持っている性質を2つ確認しておきます。あとでSBXがまねをする性質です。
性質1：子の平均は親の平均に一致する
$$\frac{87 + 58}{2} = 72.5 = \frac{90 + 55}{2}$$
同じビット位置を入れ替えているだけなので、$c_1 + c_2 = p_1 + p_2$ が常に成り立ちます。
性質2：子の間隔は親の間隔とだいたい同じ
「親どうしの距離」に対する「子どうしの距離」の比を 広がり係数（spread factor） $\beta$ と呼びます。
\beta = \frac{\vert c_1 - c_2 \vert}{\vert p_1 - p_2 \vert}
= \frac{\vert 87 - 58 \vert}{\vert 90 - 55 \vert} = \frac{29}{35} \approx 0.83
$\beta < 1$ なら子は親より縮んで（contracting）生まれ、$\beta > 1$ なら親より広がって（expanding）生まれます。1点交叉ではこの $\beta$ が 1の近くに集中し、たまに大きく外れるという分布になります。
実数値にそのまま移植するとどうなるか
さて、これを実数値ベクトルに持ってくるとどうなるでしょうか。素朴な発想は2つあります。
① 変数の値をそのまま入れ替える（バイナリの一様交叉と同じ発想）
親1: (2.0, 5.0)
親2: (7.0, 2.0)
子の候補は (2.0, 5.0) (2.0, 2.0) (7.0, 5.0) (7.0, 2.0) の4通りだけ
変数の中身が変わらないので、集団に存在する値の組み合わせしか作れません。$x_1 = 3.5$ みたいな新しい値は永遠に生まれないわけです。連続最適化としては致命的です。
② 親の平均を取る（算術交叉）
$c = \frac{p_1 + p_2}{2}$ とすれば新しい値は作れます。ただしこれは必ず親の内側、しかも1点に潰れます。世代を重ねるほど集団は中心に集まっていき、探索範囲が縮む一方になります。
図にすると違いが一目瞭然です。
- 赤い星：2つの親
- 青い点：生成された子（交叉3000回分。③は子2体ともプロットしています）
- ①：4隅のどれかにしかならない（新しい値が作れない）
- ②：必ず線分の中点1点に潰れる（多様性が消える）
- ③：SBX交叉。親のまわりに雲状に広がり、しかも親の外側にも出られる
欲しいのは③です。「親の近くを重点的に、でも時々は外側にも飛べる」交叉。これを設計したのがSBXです。
上図を生成したコードは記事の後半にまとめて載せます。
SBX交叉（Simulated Binary Crossover）
SBXは1995年にDebとAgrawalが提案した実数値GA向けの交叉です。
K. Deb and R. B. Agrawal, "Simulated Binary Crossover for Continuous Search Space," Complex Systems, vol. 9, pp. 115-148, 1995.
名前の通り「バイナリ（Simulated Binary）の1点交叉をまねる」交叉です。何をまねるかというと、さきほど確認した性質1・性質2、つまり広がり係数 $\beta$ のふるまいです。
設計の方針
SBXは以下の3つを満たすように設計されています。
| 条件 | 意味 |
|---|---|
| $c_1 + c_2 = p_1 + p_2$ | 子の平均＝親の平均（性質1） |
| $P(\beta < 1) = P(\beta > 1) = 0.5$ | 縮む交叉と広がる交叉が半々（性質2） |
| $\beta = 1$ 付近の確率が高い | 子は親の近くに生まれやすい |
βの確率密度
この条件を満たすように、$\beta$ の確率密度は次のように定められています。
P(\beta) = \begin{cases}
\dfrac{1}{2}(\eta_c + 1)\, \beta^{\eta_c} & (\beta \le 1) \quad \text{縮む側} \\[10pt]
\dfrac{1}{2}(\eta_c + 1)\, \dfrac{1}{\beta^{\eta_c + 2}} & (\beta > 1) \quad \text{広がる側}
\end{cases}
$\eta_c$ は 分布指数（distribution index） と呼ばれるパラメータです。この分布を図にすると以下のようになります。
- 左：$\eta_c$ を変えたときの $P(\beta)$。$\eta_c$ が大きいほど $\beta = 1$（＝親と同じ間隔）に鋭く集中します
- 右：後述するサンプリング式で20万個の $\beta$ を生成して理論式と重ねたもの。ぴったり一致しています
念のため、$\beta \le 1$ 側と $\beta > 1$ 側の積分がそれぞれ $0.5$ になっていることを確認しておきます。
\int_0^1 \frac{1}{2}(\eta_c+1)\beta^{\eta_c}\, d\beta
= \frac{1}{2}\Big[\beta^{\eta_c+1}\Big]_0^1 = \frac{1}{2}
\int_1^{\infty} \frac{1}{2}(\eta_c+1)\beta^{-(\eta_c+2)}\, d\beta
= \frac{1}{2}\Big[-\beta^{-(\eta_c+1)}\Big]_1^{\infty} = \frac{1}{2}
ちゃんと「縮む半分・広がる半分」になっていますね。
逆関数法でβを作る
分布が決まっても、そこから乱数を作れないと実装できません。ここで使うのが逆関数法です。
一様乱数 $u \sim U[0,1]$ を用意して、累積分布 $F(\beta) = u$ を $\beta$ について解きます。
$u \le 0.5$ のとき（縮む側）：
u = \int_0^{\beta} \frac{1}{2}(\eta_c+1)t^{\eta_c}\, dt = \frac{1}{2}\beta^{\eta_c+1}
\quad \Longrightarrow \quad
\beta = (2u)^{\frac{1}{\eta_c+1}}
$u > 0.5$ のとき（広がる側）：
u = \frac{1}{2} + \int_1^{\beta} \frac{1}{2}(\eta_c+1)t^{-(\eta_c+2)}\, dt
= 1 - \frac{1}{2}\beta^{-(\eta_c+1)}
\quad \Longrightarrow \quad
\beta = \left(\frac{1}{2(1-u)}\right)^{\frac{1}{\eta_c+1}}
まとめると以下の式になります。
\beta = \begin{cases}
(2u)^{\frac{1}{\eta_c+1}} & (u \le 0.5) \\[8pt]
\left(\dfrac{1}{2(1-u)}\right)^{\frac{1}{\eta_c+1}} & (u > 0.5)
\end{cases}
子の生成式
$\beta$ が決まれば、条件「$c_1 + c_2 = p_1 + p_2$」と「$c_1 - c_2 = -\beta(p_1 - p_2)$」の連立を解いて子が決まります。
\begin{aligned}
c_1 &= \frac{1}{2}\big[(1+\beta)\,p_1 + (1-\beta)\,p_2\big] \\
c_2 &= \frac{1}{2}\big[(1-\beta)\,p_1 + (1+\beta)\,p_2\big]
\end{aligned}
これは次のように書き直すとイメージしやすいです。
c_1 = \underbrace{\frac{p_1+p_2}{2}}_{\text{親の中点}}
- \underbrace{\beta \cdot \frac{p_2-p_1}{2}}_{\text{中点からのずれ}}
「親の中点から、$\beta$ 倍だけ左右に振り分ける」 という操作です。$\beta = 1$ なら子は親そのもの、$\beta < 1$ なら内側、$\beta > 1$ なら外側に生まれます。
実際に生成してみる
理屈だけだとピンと来ないので、親を $p_1 = 2.0$、$p_2 = 5.0$ に固定して30万組の子を作ってみます。
import numpy as np
def sample_beta(eta_c, size, rng):
"""逆関数法による広がり係数betaのサンプリング"""
u = rng.random(size)
return np.where(u <= 0.5,
(2 * u) ** (1 / (eta_c + 1)),
(1 / (2 * (1 - u))) ** (1 / (eta_c + 1)))
def sbx_simple(p1, p2, eta_c, rng):
"""境界処理なしの素のSBX"""
b = sample_beta(eta_c, p1.shape, rng)
c1 = 0.5 * ((1 + b) * p1 + (1 - b) * p2)
c2 = 0.5 * ((1 - b) * p1 + (1 + b) * p2)
return c1, c2
rng = np.random.default_rng(0)
P1, P2 = 2.0, 5.0
for eta_c in (2, 15, 50):
a, b_ = np.full(300_000, P1), np.full(300_000, P2)
c1, c2 = sbx_simple(a, b_, eta_c, rng)
both = np.concatenate([c1, c2])
outside = np.mean((both < P1) | (both > P2)) # 親の外側に出た割合
dev = np.median(np.abs(c1 - P1)) # 親からのずれ（中央値）
print(f'eta_c={eta_c:2d}: (c1+c2)/2 の平均={np.mean((c1 + c2) / 2):.4f}'
f'（親の平均=3.5）, 親の外側に出た子={outside*100:.1f}%, '
f'親からのずれの中央値={dev:.3f}')
eta_c= 2: (c1+c2)/2 の平均=3.5000（親の平均=3.5）, 親の外側に出た子=49.9%, 親からのずれの中央値=0.342
eta_c=15: (c1+c2)/2 の平均=3.5000（親の平均=3.5）, 親の外側に出た子=50.1%, 親からのずれの中央値=0.065
eta_c=50: (c1+c2)/2 の平均=3.5000（親の平均=3.5）, 親の外側に出た子=49.9%, 親からのずれの中央値=0.020
読み取れることが3つあります。
- 子の平均は3.5000。$\eta_c$ が何であっても親の平均と完全に一致しています（性質1が保たれている）
- 親の外側に出る子はほぼちょうど50%。これも $\eta_c$ によらず一定です（設計どおり）
- $\eta_c$ が変えるのは「どれくらい外側に出るか」だけ。ずれの中央値が $0.342 \to 0.065 \to 0.020$ と縮んでいきます
「外側に出る確率」と「外側に出る距離」がきれいに分離されているのが、この設計の上手いところだと思います。分布を図にすると以下のようになります。
$\eta_c = 2$ では子が広い範囲にばらまかれ、$\eta_c = 50$ ではほぼ親そのものになっているのがわかります。
ηc の選び方
まとめると $\eta_c$ の役割はこうなります。
| $\eta_c$ | 子の生まれ方 | 性格 |
|---|---|---|
| 小さい（2〜5） | 親から遠くまで飛ぶ | 探索寄り。収束は速いが解が粗くなりがち |
| 標準（15〜20） | 親のそこそこ近く | NSGA-II論文以来の定番 |
| 大きい（50〜） | ほぼ親のコピー | 活用寄り。終盤の微調整向き |
NSGA-IIの原論文をはじめ、多目的最適化の文献ではだいたい $\eta_c = 15$ か $20$ が使われています。迷ったら15でいいかと思います。
実装の細かいところ
論文の式をそのまま実装しただけだと、実は困ることが2つあります。pymoo や jMetal などのライブラリ実装ではきちんと手当てされている部分です。
① 変数ごとに確率0.5で適用する
30変数の問題で全変数にSBXを適用すると、子は親から30次元すべての方向に同時にずれます。これはずれが大きすぎて、親の良さが壊れてしまいます。
そこでNSGA-IIのオリジナル実装では、各変数に確率0.5でSBXを適用し、残りの変数は親の値をそのまま引き継ぎます。バイナリの一様交叉と同じ発想ですね。
② 変数の定義域からはみ出す問題
$\beta > 1$ で外側に出られるということは、定義域の外にも出てしまうということです。$x \in [0, 1]$ の問題で $c_1 = -0.3$ が生まれたら困ります。
素朴な対処は np.clip で境界に丸めることですが、これをやると境界に解が張り付くという副作用があります（この副作用は⑥の突然変異の方で実際に測ってみます）。
Debの実装では、そもそも定義域からはみ出さないように $\beta$ の分布を切り詰めるという方法を取っています。親 $x_1 < x_2$ と定義域 $[\text{lb}, \text{ub}]$ に対して、
\beta_{\text{lower}} = 1 + \frac{2(x_1 - \text{lb})}{x_2 - x_1}, \qquad
\alpha = 2 - \beta_{\text{lower}}^{-(\eta_c+1)}
としてから、
\beta_q = \begin{cases}
(u\alpha)^{\frac{1}{\eta_c+1}} & \left(u \le \dfrac{1}{\alpha}\right) \\[8pt]
\left(\dfrac{1}{2 - u\alpha}\right)^{\frac{1}{\eta_c+1}} & (\text{それ以外})
\end{cases}
とします。$\beta_{\text{lower}}$ は「親 $x_1$ が下限 $\text{lb}$ にちょうど届くときの $\beta$」で、それより広がる確率を $\alpha$ で再配分している、という感じです。上限側も同様に計算します。
実装コード
上の①②を入れたSBXの実装です。以降の実験ではこれを使います。
def sbx(p1, p2, eta_c, rng, lb=0.0, ub=1.0):
"""境界処理つきSBX交叉（変数ごとに確率0.5で適用）"""
n_var = len(p1)
c1, c2 = p1.copy(), p2.copy()
# ① 変数ごとに確率0.5で適用。値が同じ変数は交叉しても意味がないので除く
do = (rng.random(n_var) <= 0.5) & (np.abs(p1 - p2) > 1e-14)
x1 = np.minimum(p1, p2)[do]
x2 = np.maximum(p1, p2)[do]
if len(x1) == 0:
return c1, c2
u = rng.random(len(x1))
p = 1.0 / (eta_c + 1)
def betaq(beta): # ② 境界処理つきのbeta
alpha = 2.0 - beta ** -(eta_c + 1)
return np.where(u <= 1.0 / alpha,
(u * alpha) ** p,
(1.0 / (2.0 - u * alpha)) ** p)
bq1 = betaq(1.0 + 2.0 * (x1 - lb) / (x2 - x1)) # 下限側
bq2 = betaq(1.0 + 2.0 * (ub - x2) / (x2 - x1)) # 上限側
y1 = 0.5 * ((1 + bq1) * x1 + (1 - bq1) * x2)
y2 = 0.5 * ((1 - bq2) * x1 + (1 + bq2) * x2)
swap = rng.random(len(x1)) < 0.5 # どちらの子にどちらを渡すかを入れ替える
c1[do] = np.where(swap, y2, y1)
c2[do] = np.where(swap, y1, y2)
return np.clip(c1, lb, ub), np.clip(c2, lb, ub)
最後の swap は、子1が必ず小さい方の親側に生まれてしまう偏りをなくすための処理です。地味ですが入れないと集団に妙な偏りが出ます。
もう一つの選択肢：DE（Differential Evolution）
多目的最適化、特にMOEA/Dの文脈では DE（Differential Evolution：差分進化） の演算子もよく使われます。
H. Li and Q. Zhang, "Multiobjective optimization problems with complicated Pareto sets, MOEA/D and NSGA-II," IEEE Trans. on Evolutionary Computation, vol. 13, no. 2, pp. 284-302, 2009.
3つの解 $\mathbf{x}{r_1}, \mathbf{x}{r_2}, \mathbf{x}_{r_3}$ を選んで、差分ベクトルを足すというやり方をします。
\mathbf{v} = \mathbf{x}_{r_1} + F \cdot (\mathbf{x}_{r_2} - \mathbf{x}_{r_3})
そのあと、対象個体 $\mathbf{x}_i$ と変数ごとに確率 $CR$ で混ぜます（二項交叉）。
y_j = \begin{cases}
v_j & (\text{確率 } CR) \\
x_{i,j} & (\text{それ以外})
\end{cases}
SBXとの違いは、集団の広がりそのものがステップ幅になることです。$\mathbf{x}{r_2} - \mathbf{x}{r_3}$ は集団の散らばりを反映しているので、集団が広がっているうちは大きく動き、収束してくると自動的に小さく動きます。ステップ幅の自動調整が入っているようなものです。
| SBX | DE | |
|---|---|---|
| 親の数 | 2 | 3（+ 対象個体） |
| ステップ幅の決まり方 | $\eta_c$ で人が指定 | 集団の広がりから自動的に決まる |
| 変数間の相関 | 扱えない（変数ごとに独立） | 差分ベクトルの向きに沿って動ける |
| 主なパラメータ | $\eta_c$, $p_c$ | $F$, $CR$ |
| よく使われる場面 | NSGA-II系の標準 | MOEA/D-DE、変数間に相関がある問題 |
「変数間の相関を扱えるか」がけっこう本質的な違いです。SBXは変数ごとに独立に処理するので、$x_1 + x_2 = 1$ のような斜めの構造に沿って動くのが苦手だったりします。
実験：交叉を変えるとパレートフロントはどう変わるか
ここまでの話を実際に確かめます。②のNSGA-IIに交叉だけ差し替えて、④で紹介したベンチマーク問題のひとつ ZDT1（30変数）を解きます。
ZDT1はパレートフロントが $f_2 = 1 - \sqrt{f_1}$ と解析的にわかっているので、④で説明した IGD（小さいほど良い）で評価できます。
f_1(\mathbf{x}) = x_1, \qquad
g(\mathbf{x}) = 1 + \frac{9}{n-1}\sum_{i=2}^{n} x_i, \qquad
f_2(\mathbf{x}) = g\left(1 - \sqrt{f_1/g}\right)
真のパレートフロントは $g = 1$、つまり $x_2 = \cdots = x_{30} = 0$ のときです。
集団100・250世代・5シードの平均を取りました。結果はこうなりました。
| 交叉 | IGD（50世代時点） | IGD（250世代時点） | $x_2 \ldots x_{30}$ の平均（最適は0） | $f_1$ の範囲（理想は1.0） |
|---|---|---|---|---|
| SBX（$\eta_c=15$） | 0.1169 | 0.0050 | 0.0003 | 1.000 |
| 算術交叉（親の平均） | 2.0680 | 1.3761 | 0.2029 | 0.435 |
| DE（$F=0.5, CR=0.5$） | 0.0373 | 0.0046 | 0.0002 | 1.000 |
| DE（$F=0.5, CR=1.0$） | 1.2592 | 0.2837 | 0.0703 | 1.000 |
※突然変異はすべて多項式突然変異（$\eta_m = 20$, $p_m = 1/30$）で統一
※「$f_1$ の範囲」は得られた解集合の $f_1$ の最大値−最小値。ZDT1の真のパレートフロントは $f_1 \in [0, 1]$ なので、1.0に近いほどフロント全体をカバーできている＝多様性が保たれている、という指標です
得られた解とパレートフロントを重ねると以下のようになります。
- 灰色の破線：真のパレートフロント
- 青い点：250世代後に得られた解（シード0）
- 真ん中の算術交叉だけ縦軸の位置がまるで違います。$f_2$ が2以上の場所で止まってしまっているわけです
収束の推移を見るとさらにはっきりします。
算術交叉はなぜ止まるのか
IGDが1.3761、つまりほぼ何も最適化できていません。理由は $x_2 \ldots x_{30}$ の平均が0.2029のままであることから読み取れます。
ZDT1で $g$ を小さくするには $x_2 \ldots x_{30}$ を0に近づける必要がありますが、算術交叉は親の平均しか作れないので、集団の最小値より小さい値を絶対に作れません。子は常に親たちが張る範囲の内側にしか生まれないからです。あとは突然変異の細かいランダム性だけが頼りになるので、まったく進まないというわけです。
ここでさっきの「SBXの子は50%が親の外側に生まれる」が効いてきます。外側に出られるからこそ、集団は今いる範囲の外へ進めるわけです。この一点だけでもSBXを使う理由になるかと思います。
$f_1$ の範囲（多様性の指標）が0.435しかないのも、集団が中心に潰れていっていることの表れです。
DEはCRの設定で結果がひっくり返った
個人的に一番おもしろかったのがここです。同じDEでも、
- $CR = 0.5$：IGD 0.0046 でSBXより良い。しかも50世代時点で0.0373と収束が圧倒的に速い
- $CR = 1.0$：IGD 0.2837 で全然ダメ
$CR = 1.0$ は「変数を全部、差分ベクトルで置き換える」という設定です。これがなぜダメかというと、$\mathbf{x}{r_1}, \mathbf{x}{r_2}, \mathbf{x}_{r_3}$ を集団全体からランダムに選んでいるからです。無関係な3つの解から作った差分ベクトルで30変数すべてを一気に書き換えたら、そりゃ壊れますよね。。。
ここで③のMOEA/Dの話に戻ります。MOEA/D-DEは $CR = 1.0$ を使うのですが、親を近傍 $B(i)$ の中からしか選びません。近傍の解は互いに似ているので差分ベクトルが小さく、$CR = 1.0$ でも壊れないわけです。
③で「近傍の情報を共有して効率よく探索」と書きましたが、近傍制限は演算子を壊さないための仕掛けでもあった、ということです。MOEA/DとDE演算子の相性の良さは、こういうところから来ているのだと思います。
順列表現の場合（配送最適化入門⑥との対応）
ちなみに解が順列の場合、SBXは使えません。[0,3,1,4,2] と [2,0,4,1,3] の平均を取っても順列になりませんし、変数を入れ替えると都市が重複します。
順列向けにはOX（Order Crossover）やPMX（Partially Mapped Crossover）という専用の交叉があります。こちらは配送最適化入門シリーズ⑥「遺伝的アルゴリズム（GA）でルートを進化」で説明しているので、興味があれば覗いてみてください。
| 実数値（本記事） | 順列（配送⑥） | |
|---|---|---|
| 交叉 | SBX・DE | OX・PMX |
| 突然変異 | 多項式突然変異 | スワップ・2-opt |
| 守るべき制約 | 変数の定義域 $[\text{lb}, \text{ub}]$ | 「全都市がちょうど1回ずつ」 |
「選択 → 交叉 → 突然変異」という骨格は共通で、表現に合わせて演算子だけ差し替えるというのが、GAという枠組みの強さかなと思います。
図を生成したコード
記事中の図を生成したコードは以下のとおりです。
import numpy as np
import matplotlib
matplotlib.use('Agg')
import matplotlib.pyplot as plt
plt.rcParams['font.family'] = 'Yu Gothic'
def beta_pdf(b, eta_c):
"""betaの理論確率密度"""
return np.where(b <= 1,
0.5 * (eta_c + 1) * b ** eta_c,
0.5 * (eta_c + 1) / b ** (eta_c + 2))
rng = np.random.default_rng(0)
# --- βの確率密度とサンプリングの検証 ---
fig, axes = plt.subplots(1, 2, figsize=(11, 4))
b_grid = np.linspace(0.001, 3.0, 600)
for eta_c, c in {2: '#3498db', 5: '#27ae60',
15: '#e67e22', 30: '#c0392b'}.items():
axes[0].plot(b_grid, beta_pdf(b_grid, eta_c), color=c, lw=2,
label=f'ηc = {eta_c}')
axes[0].axvline(1.0, color='gray', ls='--', lw=1)
axes[0].set_xlabel('広がり係数 β'); axes[0].set_ylabel('確率密度 P(β)')
axes[0].set_ylim(0, 4.2); axes[0].legend(); axes[0].grid(True, alpha=0.3)
b_sample = sample_beta(15, 200_000, rng)
axes[1].hist(b_sample[b_sample < 3.0], bins=120, density=True,
color='#95a5a6', alpha=0.75, label='サンプリング結果（20万個）')
axes[1].plot(b_grid, beta_pdf(b_grid, 15), color='#c0392b', lw=2,
label='理論式 P(β)')
axes[1].set_xlabel('広がり係数 β'); axes[1].set_ylabel('確率密度')
axes[1].legend(); axes[1].grid(True, alpha=0.3)
plt.tight_layout()
plt.savefig('20260812_sbx_beta_pdf.png', dpi=150)
# --- 2変数での子の広がり方（3方式の比較）---
P = np.array([2.0, 5.0]); Q = np.array([7.0, 2.0]); n_s = 3000
mask = rng.random((n_s, 2)) < 0.5
swap_child = np.where(mask, P, Q) # ① 変数の入れ替えだけ
mean_child = np.tile((P + Q) / 2, (n_s, 1)) # ② 算術交叉
b = sample_beta(15, (n_s, 2), rng) # ③ SBX（子2体とも描く）
sbx_child = np.vstack([0.5 * ((1 + b) * P + (1 - b) * Q),
0.5 * ((1 - b) * P + (1 + b) * Q)])
fig, axes = plt.subplots(1, 3, figsize=(13, 4.2), sharex=True, sharey=True)
data = [(swap_child, '① 変数の入れ替えだけ', 90, 0.5),
(mean_child, '② 算術交叉：親の平均', 90, 0.5),
(sbx_child, '③ SBX交叉（ηc=15）', 8, 0.2)]
for ax, (child, title, size, al) in zip(axes, data):
ax.scatter(child[:, 0], child[:, 1], s=size, color='#3498db',
alpha=al, zorder=3)
ax.plot([P[0], Q[0]], [P[1], Q[1]], color='#7f8c8d', ls='--', lw=1.2)
ax.scatter(*P, s=160, marker='*', color='#c0392b', zorder=5)
ax.scatter(*Q, s=160, marker='*', color='#c0392b', zorder=5)
ax.set_title(title, fontsize=11); ax.set_xlabel('x₁')
ax.grid(True, alpha=0.3)
axes[0].set_ylabel('x₂'); axes[0].set_xlim(-0.5, 9.5); axes[0].set_ylim(-0.5, 8.0)
plt.tight_layout()
plt.savefig('20260812_sbx_2d.png', dpi=150)
NSGA-IIに各交叉を組み込んだ実験のコードは、次回⑥で全体をまとめて載せます。
まとめ
今回は交叉、特にSBX交叉を説明しました。
- 交叉は解の表現とセット：実数値ベクトルにバイナリGAの交叉をそのまま持ってきても機能しない
- SBXはバイナリの1点交叉の性質をまねた交叉：「子の平均＝親の平均」「縮む／広がるが半々」を満たすように $\beta$ の分布を設計している
- $\beta$ の式は逆関数法から導出できる：実装でよく見る式には出どころがある
- $\eta_c$ は「どれだけ親の近くに子を作るか」だけを制御する。外側に出る確率は常に50%で一定
- 親の外側に出られることが本質：算術交叉のように内側しか作れない交叉ではIGDが0.0050→1.3761と壊滅した
- DEも有力な選択肢。ただし $CR=1.0$ は近傍制限（MOEA/D-DE）とセットで使うもので、集団全体から親を選ぶ設定だと壊れる
次回⑥は突然変異（多項式突然変異）と選択（トーナメント選択・混雑比較演算子）を説明して、オペレータ一式を揃えます。「突然変異なんて微調整でしょ」と思われがちですが、外すとNSGA-IIが途中で完全に止まったので、その辺も含めて実験結果を出します。
参考にした記事や本、論文等
- K. Deb and R. B. Agrawal, "Simulated Binary Crossover for Continuous Search Space," Complex Systems, vol. 9, pp. 115-148, 1995.
- K. Deb, A. Pratap, S. Agarwal, and T. Meyarivan, "A fast and elitist multiobjective genetic algorithm: NSGA-II," IEEE Trans. on Evolutionary Computation, vol. 6, no. 2, pp. 182-197, 2002.
- H. Li and Q. Zhang, "Multiobjective optimization problems with complicated Pareto sets, MOEA/D and NSGA-II," IEEE Trans. on Evolutionary Computation, vol. 13, no. 2, pp. 284-302, 2009.
- E. Zitzler, K. Deb, and L. Thiele, "Comparison of multiobjective evolutionary algorithms: Empirical results," Evolutionary Computation, vol. 8, no. 2, pp. 173-195, 2000.（ZDT問題群）