はじめに
粒子群最適化（Particle Swarm Optimization, PSO）は、鳥や魚の群れの集合行動から着想を得たメタヒューリスティクス最適化手法です。Kennedy と Eberhart によって 1995 年に提案され、勾配情報を必要としないブラックボックス最適化アルゴリズムとして、連続値の多峰性問題で広く用いられてきました。
PSO の特徴は、各**粒子（particle）**が「位置」と「速度」を持ち、自分自身が見つけたこれまでの最良解（個人最良 pbest）と、群れ全体が共有する最良解（全体最良 gbest）の双方に引き寄せられながら探索空間を移動する点にあります。これにより、ランダム探索だけでは到達しにくい大域解の周辺へ、群れ全体の知能で素早く収束していきます。
本稿では PSO の更新式を慣性重み \(w\) 、認知係数 \(c_1\) 、社会係数 \(c_2\) の役割から丁寧に整理し、Sphere 関数と Rastrigin 関数の最小化を題材に Python 実装と収束の様子を可視化します。
アルゴリズム
粒子の状態
集団サイズを \(N\) 、探索空間の次元を \(D\) とします。各粒子 \(i \in \{1, \dots, N\}\) は、
- 位置 \(\mathbf{x}_i^t \in \mathbb{R}^D\)
- 速度 \(\mathbf{v}_i^t \in \mathbb{R}^D\)
- 個人最良位置 \(\mathbf{p}_i^t\) （粒子 \(i\) がこれまでに到達した最良位置）
を保持します。さらに群れ全体で共有する全体最良位置 \(\mathbf{g}^t = \arg\min_{i} f(\mathbf{p}_i^t)\) を持ちます。
速度更新
速度は次の式で更新されます。
$$ \mathbf{v}_i^{t+1} = w , \mathbf{v}_i^{t}
- c_1 \mathbf{r}_1 \odot (\mathbf{p}_i^{t} - \mathbf{x}_i^{t})
- c_2 \mathbf{r}_2 \odot (\mathbf{g}^{t} - \mathbf{x}_i^{t}) \tag{1} $$
ここで \(\mathbf{r}_1, \mathbf{r}_2 \sim U[0,1]^D\) は次元ごとに独立な一様乱数ベクトル、\(\odot\) は要素積です。右辺の各項は次の意味を持ちます。
- 第1項 \(w \, \mathbf{v}_i^t\) : 慣性項。直前の速度をどれだけ引き継ぐか
- 第2項 \(c_1 \mathbf{r}_1 \odot (\mathbf{p}_i^t - \mathbf{x}_i^t)\) : 認知項。自分自身の経験（pbest）への引力
- 第3項 \(c_2 \mathbf{r}_2 \odot (\mathbf{g}^t - \mathbf{x}_i^t)\) : 社会項。群れ全体の経験（gbest）への引力
位置更新
位置は速度を用いて単純に進めます。
\[ \mathbf{x}_i^{t+1} = \mathbf{x}_i^{t} + \mathbf{v}_i^{t+1} \tag{2} \]pbest と gbest の更新
各イテレーションで、粒子の現在位置の評価値が個人最良を更新したら \(\mathbf{p}_i\) を置き換えます。
\[ \mathbf{p}_i^{t+1} = \begin{cases} \mathbf{x}_i^{t+1} & \text{if } f(\mathbf{x}_i^{t+1}) < f(\mathbf{p}_i^{t}) \\ \mathbf{p}_i^{t} & \text{otherwise} \end{cases} \tag{3} \]そして全体最良も同様に更新します。
\[ \mathbf{g}^{t+1} = \arg\min_{i} f(\mathbf{p}_i^{t+1}) \tag{4} \]ハイパーパラメータの役割
慣性重み \(w\)
\(w\) は探索の 多様性（exploration） と 集中（exploitation） のバランスを制御します。
- \(w\) が大きい（例: \(w \approx 0.9\) ）: 粒子が直前の速度を強く保つため、広域探索が促進される
- \(w\) が小さい（例: \(w \approx 0.4\) ）: 速度が急速に減衰し、局所的な精緻化が進む
実用上は、序盤に広く探索し終盤に収束させるため、線形に減衰させる線形減衰慣性重み（linearly decreasing inertia weight, LDIW） が定石です。
\[ w_t = w_{\max} - (w_{\max} - w_{\min}) \cdot \frac{t}{T_{\max}} \tag{5} \]典型値は \(w_{\max} = 0.9, w_{\min} = 0.4\) です。
認知係数 \(c_1\) と社会係数 \(c_2\)
\(c_1\) と \(c_2\) は、それぞれ自分自身の経験と群れの経験への引力の強さを決めます。
- \(c_1 \gg c_2\) : 各粒子が自分の経験を重視し、独立的な探索が進む
- \(c_1 \ll c_2\) : 群れの最良解への収束が支配的となり、局所解に集中しやすい
Clerc & Kennedy (2002) は、収束性の解析から \(c_1 = c_2 = 2.05\) 、慣性重み \(w = 0.7298\) という**収束係数（constriction coefficient）**を提案しています。実装時のデフォルトとしては \(w = 0.729, c_1 = c_2 = 1.4944\) がよく使われます。
速度のクリッピング
粒子が探索領域外へ発散することを防ぐため、各次元の速度を \(|v_{i,d}| \leq v_{\max}\) にクリップする運用が一般的です。\(v_{\max}\) は探索範囲幅の 10〜20% 程度が目安です。
ベンチマーク関数
Sphere 関数（単峰性）
\[ f_{\text{sphere}}(\mathbf{x}) = \sum_{i=1}^{D} x_i^2 \tag{6} \]大域最小値は \(f(\mathbf{0}) = 0\) 。単峰性なので PSO は容易に収束します。
Rastrigin 関数（多峰性）
\[ f_{\text{rastrigin}}(\mathbf{x}) = 10 D + \sum_{i=1}^{D} \left[ x_i^2 - 10 \cos(2\pi x_i) \right] \tag{7} \]大域最小値は \(f(\mathbf{0}) = 0\) ですが、無数の局所解を持つため、慣性重みや係数の設定によっては粒子が早期に局所解に集中して停滞することがあります。
Python 実装
NumPy のみを用いた基本的な PSO 実装です。線形減衰慣性重みと速度クリッピングを組み込んでいます。
import numpy as np
import matplotlib.pyplot as plt
def sphere(x):
"""Sphere 関数（単峰性）"""
return np.sum(x**2, axis=-1)
def rastrigin(x):
"""Rastrigin 関数（多峰性）"""
return 10 * x.shape[-1] + np.sum(x**2 - 10 * np.cos(2 * np.pi * x), axis=-1)
class PSO:
def __init__(self, dim, n_particles=40, bounds=(-5.12, 5.12),
w_max=0.9, w_min=0.4, c1=1.4944, c2=1.4944,
v_max_ratio=0.2):
self.dim = dim
self.n_particles = n_particles
self.bounds = bounds
self.w_max = w_max
self.w_min = w_min
self.c1 = c1
self.c2 = c2
# 速度上限は探索範囲幅の v_max_ratio 倍
self.v_max = v_max_ratio * (bounds[1] - bounds[0])
def initialize(self):
low, high = self.bounds
x = np.random.uniform(low, high, (self.n_particles, self.dim))
v = np.random.uniform(-self.v_max, self.v_max,
(self.n_particles, self.dim))
return x, v
def run(self, objective, max_iter=200):
x, v = self.initialize()
f = objective(x)
# pbest / gbest 初期化
pbest = x.copy()
pbest_f = f.copy()
g_idx = int(np.argmin(pbest_f))
gbest = pbest[g_idx].copy()
gbest_f = float(pbest_f[g_idx])
history = [gbest_f]
for t in range(max_iter):
# 線形減衰慣性重み
w = self.w_max - (self.w_max - self.w_min) * (t / max_iter)
# 速度更新（式 1）
r1 = np.random.random((self.n_particles, self.dim))
r2 = np.random.random((self.n_particles, self.dim))
v = (w * v
+ self.c1 * r1 * (pbest - x)
+ self.c2 * r2 * (gbest - x))
v = np.clip(v, -self.v_max, self.v_max)
# 位置更新（式 2）
x = x + v
x = np.clip(x, *self.bounds)
# 評価と pbest / gbest 更新
f = objective(x)
improved = f < pbest_f
pbest[improved] = x[improved]
pbest_f[improved] = f[improved]
g_idx = int(np.argmin(pbest_f))
if pbest_f[g_idx] < gbest_f:
gbest = pbest[g_idx].copy()
gbest_f = float(pbest_f[g_idx])
history.append(gbest_f)
return gbest, gbest_f, history
# --- 実行 ---
np.random.seed(42)
pso_sphere = PSO(dim=10)
_, best_sphere, hist_sphere = pso_sphere.run(sphere, max_iter=200)
np.random.seed(42)
pso_rastrigin = PSO(dim=10)
_, best_rastrigin, hist_rastrigin = pso_rastrigin.run(rastrigin, max_iter=200)
print(f"Sphere 最良適応度: {best_sphere:.6e}")
print(f"Rastrigin 最良適応度: {best_rastrigin:.6e}")
# --- 収束曲線のプロット ---
plt.figure(figsize=(10, 5))
plt.plot(hist_sphere, label='Sphere (10D)')
plt.plot(hist_rastrigin, label='Rastrigin (10D)')
plt.xlabel('Iteration')
plt.ylabel('Best Objective Value')
plt.title('PSO Convergence on Sphere and Rastrigin')
plt.yscale('log')
plt.grid(True, alpha=0.3)
plt.legend()
plt.tight_layout()
plt.show()
Sphere 関数では数十イテレーションで \(10^{-10}\) 程度まで収束する一方、Rastrigin 関数では多峰性により早期に停滞しがちで、慣性重みの初期値・粒子数・係数の調整が結果を大きく左右します。
慣性重みが収束に与える影響
固定の \(w\) を変化させた場合の比較も簡単に再現できます。
import numpy as np
import matplotlib.pyplot as plt
results = {}
for w in [0.3, 0.5, 0.7, 0.9]:
np.random.seed(0)
pso = PSO(dim=10, w_max=w, w_min=w) # 固定慣性重み
_, _, history = pso.run(rastrigin, max_iter=200)
results[w] = history
plt.figure(figsize=(10, 5))
for w, history in results.items():
plt.plot(history, label=f'w = {w}')
plt.xlabel('Iteration')
plt.ylabel('Best Objective Value')
plt.title('Effect of Inertia Weight on Rastrigin (10D)')
plt.yscale('log')
plt.grid(True, alpha=0.3)
plt.legend()
plt.tight_layout()
plt.show()
一般に \(w = 0.3\) では序盤から大域探索が弱く局所解に捕まりやすく、\(w = 0.9\) では収束が遅くなる傾向があります。LDIW のように動的に調整する戦略が頑健です。
他の最適化手法との比較
| 手法 | 集団ベース | 主なメカニズム | 多峰性への耐性 | 評価コストが高い問題 |
|---|---|---|---|---|
| PSO | Yes | 速度・位置の更新（pbest/gbest） | 中 | 不向き |
| GA | Yes | 選択・交叉・突然変異 | 高 | 不向き |
| CEM | Yes | 確率分布の反復更新 | 中 | 不向き |
| SA | No | 温度付き確率的受容 | 高 | 不向き |
| ベイズ最適化 | No | ガウス過程による代理モデル | 中 | 得意 |
PSO は実装が極めて簡潔でハイパーパラメータも少なく、連続値の中規模次元問題に対して高い実用性を持ちます。一方で、評価が極端に高コストな問題（実機実験など）では、評価回数を抑えられるベイズ最適化の方が有利です。
まとめ
- PSO は粒子の位置・速度更新（式 1, 式 2）と pbest/gbest の更新（式 3, 式 4）から構成される非常に単純なアルゴリズム
- 慣性重み \(w\) 、認知係数 \(c_1\) 、社会係数 \(c_2\) がそれぞれ多様性・自己経験・群れ経験の重みを決める
- 線形減衰慣性重み（式 5）や速度クリッピングは、Rastrigin のような多峰性関数で実用上ほぼ必須
- Sphere では非常に高速に収束し、Rastrigin では設定次第で挙動が大きく変わる
最適化シリーズの他の記事と読み比べると、それぞれの手法が得意とする問題クラスがより明確になります。
関連記事
- 遺伝的アルゴリズム（GA）の基礎とPython実装 - 進化計算ベースのメタヒューリスティクス。PSO は速度更新で連続的に動くのに対し、GA は交叉・突然変異で離散的に世代交代します。
- 焼きなまし法（Simulated Annealing）の仕組みとPython実装 - 単一解を保持する SA との比較で、集団ベース手法である PSO の特徴が見えてきます。
- クロスエントロピー法：モンテカルロ最適化の実践的手法 - 確率分布の反復更新で探索する CEM。PSO の「群れ全体が gbest に引き寄せられる」過程と通底する点があります。
- ベイズ最適化の基礎とPython実装 - 評価コストが高い問題で威力を発揮する手法。PSO とは想定する問題設定が大きく異なるため使い分けが重要です。
参考文献
- Kennedy, J., & Eberhart, R. (1995). “Particle Swarm Optimization”. Proceedings of IEEE International Conference on Neural Networks, 4, 1942-1948.
- Shi, Y., & Eberhart, R. (1998). “A modified particle swarm optimizer”. Proceedings of IEEE International Conference on Evolutionary Computation, 69-73.
- Clerc, M., & Kennedy, J. (2002). “The particle swarm – explosion, stability, and convergence in a multidimensional complex space”. IEEE Transactions on Evolutionary Computation, 6(1), 58-73.