ベルマン方程式は強化学習の理論的中核をなす方程式です。価値関数の再帰的な関係を記述し、「現在の価値 = 即時報酬 + 将来の価値の割引」という直感的な構造を持ちます。本記事では、ベルマン期待方程式と最適方程式を導出し、動的計画法（価値反復法・方策反復法）を実装します。
本記事の内容
- ベルマン期待方程式の導出
- ベルマン最適方程式の導出
- 最適方策と最適価値関数
- 価値反復法（Value Iteration）
- 方策反復法（Policy Iteration）
- Pythonでの実装
前提知識
この記事を読む前に、以下の記事を読んでおくと理解が深まります。
ベルマン期待方程式
状態価値関数のベルマン期待方程式
方策 $\pi$ のもとでの状態価値関数 $V^\pi(s)$ について、収益 $G_t$ の再帰的定義 $G_t = R_{t+1} + \gamma G_{t+1}$ を利用します。
$$ \begin{align} V^\pi(s) &= E_\pi[G_t | S_t = s] \\ &= E_\pi[R_{t+1} + \gamma G_{t+1} | S_t = s] \\ &= E_\pi[R_{t+1} | S_t = s] + \gamma E_\pi[G_{t+1} | S_t = s] \end{align} $$
行動 $a$ を方策 $\pi$ に従って選び、遷移確率 $P(s’|s,a)$ に従って次状態 $s’$ に遷移することを展開します。
$$ \begin{align} V^\pi(s) &= \sum_{a} \pi(a|s) \sum_{s’} P(s’|s, a) \left[R(s, a, s’) + \gamma E_\pi[G_{t+1} | S_{t+1} = s’]\right] \\ &= \sum_{a} \pi(a|s) \sum_{s’} P(s’|s, a) \left[R(s, a, s’) + \gamma V^\pi(s’)\right] \end{align} $$
$$ \boxed{V^\pi(s) = \sum_{a} \pi(a|s) \sum_{s’} P(s’|s, a) \left[R(s, a, s’) + \gamma V^\pi(s’)\right]} $$
行動価値関数のベルマン期待方程式
同様に $Q^\pi(s, a)$ について:
$$ \boxed{Q^\pi(s, a) = \sum_{s’} P(s’|s, a) \left[R(s, a, s’) + \gamma \sum_{a’} \pi(a’|s’) Q^\pi(s’, a’)\right]} $$
行列形式
全状態をベクトルとして書くと:
$$ \bm{V}^\pi = \bm{R}^\pi + \gamma \bm{P}^\pi \bm{V}^\pi $$
ここで $\bm{P}^\pi$ は方策 $\pi$ のもとでの状態遷移行列、$\bm{R}^\pi$ は期待即時報酬ベクトルです。
これを解くと:
$$ \bm{V}^\pi = (\bm{I} – \gamma \bm{P}^\pi)^{-1} \bm{R}^\pi $$
小さなMDPでは直接解けますが、状態空間が大きい場合は反復法を使います。
最適方策と最適価値関数
最適状態価値関数 $V^*(s)$ は全ての方策の中で最大の価値をとる価値関数です。
$$ V^*(s) = \max_\pi V^\pi(s), \quad \forall s \in S $$
最適行動価値関数 $Q^*(s, a)$ は:
$$ Q^*(s, a) = \max_\pi Q^\pi(s, a), \quad \forall s \in S, a \in A $$
最適方策 $\pi^*$ は最適行動価値関数から得られます。
$$ \pi^*(s) = \arg\max_a Q^*(s, a) $$
ベルマン最適方程式
状態価値関数のベルマン最適方程式
最適方策は各状態で最良の行動を貪欲に選ぶので:
$$ \begin{align} V^*(s) &= \max_a Q^*(s, a) \\ &= \max_a \sum_{s’} P(s’|s, a) \left[R(s, a, s’) + \gamma V^*(s’)\right] \end{align} $$
$$ \boxed{V^*(s) = \max_a \sum_{s’} P(s’|s, a) \left[R(s, a, s’) + \gamma V^*(s’)\right]} $$
行動価値関数のベルマン最適方程式
$$ \boxed{Q^*(s, a) = \sum_{s’} P(s’|s, a) \left[R(s, a, s’) + \gamma \max_{a’} Q^*(s’, a’)\right]} $$
ベルマン期待方程式との違いは、方策 $\pi$ による重み付き平均の代わりに $\max$ を使う点です。
価値反復法（Value Iteration）
価値反復法は、ベルマン最適方程式を反復的に適用して $V^*$ を求めるアルゴリズムです。
アルゴリズム
- $V(s) = 0$ で初期化
- 全状態 $s$ について更新:
$$ V(s) \leftarrow \max_a \sum_{s’} P(s’|s, a) \left[R(s, a, s’) + \gamma V(s’)\right] $$
- $\max_s |V_{\text{new}}(s) – V_{\text{old}}(s)| < \varepsilon$ なら収束
収束後、最適方策を抽出:
$$ \pi^*(s) = \arg\max_a \sum_{s’} P(s’|s, a) \left[R(s, a, s’) + \gamma V^*(s’)\right] $$
方策反復法（Policy Iteration）
方策反復法は、方策評価と方策改善を交互に繰り返して最適方策を求めます。
アルゴリズム
- 適当な方策 $\pi$ で初期化
- 方策評価: $\pi$ のもとで $V^\pi$ を計算（ベルマン期待方程式を反復解法）
- 方策改善: 新しい方策を $\pi'(s) = \arg\max_a Q^\pi(s, a)$ で更新
- $\pi’ = \pi$ なら収束。さもなくば $\pi \leftarrow \pi’$ としてステップ2へ
方策改善定理
$\pi'(s) = \arg\max_a Q^\pi(s, a)$ とすると、$V^{\pi’}(s) \geq V^\pi(s)$ が全ての $s$ で成り立ちます。
Pythonでの実装
import numpy as np
import matplotlib.pyplot as plt
class GridWorldMDP:
"""4x4 グリッドワールド MDP"""
def __init__(self, size=4, gamma=0.9):
self.size = size
self.n_states = size * size
self.n_actions = 4
self.gamma = gamma
self.actions = [(-1, 0), (1, 0), (0, -1), (0, 1)]
self.action_names = ['↑', '↓', '←', '→']
self.terminal_states = [0, self.n_states - 1]
self.P = np.zeros((self.n_states, self.n_actions, self.n_states))
self.R = np.zeros((self.n_states, self.n_actions, self.n_states))
self._build_transitions()
def _state_to_pos(self, s):
return s // self.size, s % self.size
def _pos_to_state(self, r, c):
return r * self.size + c
def _build_transitions(self):
for s in range(self.n_states):
if s in self.terminal_states:
for a in range(self.n_actions):
self.P[s, a, s] = 1.0
continue
r, c = self._state_to_pos(s)
for a, (dr, dc) in enumerate(self.actions):
nr, nc = r + dr, c + dc
if 0 <= nr < self.size and 0 <= nc < self.size:
s_next = self._pos_to_state(nr, nc)
else:
s_next = s
self.P[s, a, s_next] = 1.0
self.R[s, a, s_next] = -1.0
def value_iteration(mdp, tol=1e-8, max_iter=1000):
"""価値反復法"""
V = np.zeros(mdp.n_states)
history = [V.copy()]
for iteration in range(max_iter):
V_new = np.zeros(mdp.n_states)
for s in range(mdp.n_states):
q_values = np.zeros(mdp.n_actions)
for a in range(mdp.n_actions):
for s_next in range(mdp.n_states):
q_values[a] += mdp.P[s, a, s_next] * (
mdp.R[s, a, s_next] + mdp.gamma * V[s_next]
)
V_new[s] = np.max(q_values)
delta = np.max(np.abs(V_new - V))
V = V_new
history.append(V.copy())
if delta < tol:
print(f"価値反復法: {iteration + 1} 回で収束")
break
# 最適方策の抽出
pi = np.zeros(mdp.n_states, dtype=int)
for s in range(mdp.n_states):
q_values = np.zeros(mdp.n_actions)
for a in range(mdp.n_actions):
for s_next in range(mdp.n_states):
q_values[a] += mdp.P[s, a, s_next] * (
mdp.R[s, a, s_next] + mdp.gamma * V[s_next]
)
pi[s] = np.argmax(q_values)
return V, pi, history
def policy_iteration(mdp, tol=1e-8, max_iter=100):
"""方策反復法"""
# ランダム方策で初期化
pi = np.zeros(mdp.n_states, dtype=int)
V_history = []
for iteration in range(max_iter):
# 方策評価
pi_prob = np.zeros((mdp.n_states, mdp.n_actions))
for s in range(mdp.n_states):
pi_prob[s, pi[s]] = 1.0
V = np.zeros(mdp.n_states)
for _ in range(1000):
V_new = np.zeros(mdp.n_states)
for s in range(mdp.n_states):
for a in range(mdp.n_actions):
for s_next in range(mdp.n_states):
V_new[s] += pi_prob[s, a] * mdp.P[s, a, s_next] * (
mdp.R[s, a, s_next] + mdp.gamma * V[s_next]
)
if np.max(np.abs(V_new - V)) < tol:
break
V = V_new
V_history.append(V.copy())
# 方策改善
pi_new = np.zeros(mdp.n_states, dtype=int)
policy_stable = True
for s in range(mdp.n_states):
q_values = np.zeros(mdp.n_actions)
for a in range(mdp.n_actions):
for s_next in range(mdp.n_states):
q_values[a] += mdp.P[s, a, s_next] * (
mdp.R[s, a, s_next] + mdp.gamma * V[s_next]
)
pi_new[s] = np.argmax(q_values)
if pi_new[s] != pi[s]:
policy_stable = False
pi = pi_new
if policy_stable:
print(f"方策反復法: {iteration + 1} 回で収束")
break
return V, pi, V_history
# --- 実行と比較 ---
mdp = GridWorldMDP(size=4, gamma=0.9)
V_vi, pi_vi, history_vi = value_iteration(mdp)
V_pi, pi_pi, history_pi = policy_iteration(mdp)
# --- 結果の可視化 ---
fig, axes = plt.subplots(1, 3, figsize=(18, 5))
# 価値反復法の結果
V_grid = V_vi.reshape(mdp.size, mdp.size)
im1 = axes[0].imshow(V_grid, cmap='RdYlGn', interpolation='nearest')
for i in range(mdp.size):
for j in range(mdp.size):
axes[0].text(j, i, f'{V_grid[i, j]:.2f}',
ha='center', va='center', fontsize=11, fontweight='bold')
axes[0].set_title('Value Iteration: V*(s)', fontsize=13)
plt.colorbar(im1, ax=axes[0])
# 最適方策
policy_grid = np.array([mdp.action_names[a] for a in pi_vi]).reshape(mdp.size, mdp.size)
axes[1].imshow(np.zeros((mdp.size, mdp.size)), cmap='Greys', vmin=0, vmax=1, alpha=0.1)
for i in range(mdp.size):
for j in range(mdp.size):
axes[1].text(j, i, policy_grid[i, j],
ha='center', va='center', fontsize=20, fontweight='bold')
axes[1].set_title('Optimal Policy π*(s)', fontsize=13)
axes[1].set_xticks(range(mdp.size))
axes[1].set_yticks(range(mdp.size))
# 収束過程
for s in [5, 6, 9, 10]: # 中央付近の状態
values = [h[s] for h in history_vi]
axes[2].plot(values, linewidth=2, label=f'V(s={s})')
axes[2].set_xlabel('Iteration', fontsize=12)
axes[2].set_ylabel('V(s)', fontsize=12)
axes[2].set_title('Value Iteration Convergence', fontsize=13)
axes[2].legend(fontsize=10)
axes[2].grid(True, alpha=0.3)
plt.tight_layout()
plt.show()
まとめ
本記事では、ベルマン方程式の導出と動的計画法について解説しました。
- ベルマン期待方程式: $V^\pi(s) = \sum_a \pi(a|s) \sum_{s’} P(s’|s,a)[R + \gamma V^\pi(s’)]$
- ベルマン最適方程式: $V^*(s) = \max_a \sum_{s’} P(s’|s,a)[R + \gamma V^*(s’)]$
- 価値反復法: ベルマン最適方程式を反復適用して $V^*$ を求める
- 方策反復法: 方策評価と方策改善を交互に繰り返す
- どちらも最適方策に収束することが保証されている
次のステップとして、以下の記事も参考にしてください。