用Python实战推导贝尔曼方程强化学习核心原理的代码化理解许多初学者在接触强化学习时都会被贝尔曼方程这个看似复杂的数学概念吓退。但今天我要告诉你一个秘密贝尔曼方程的本质其实就是一个自我更新的价值计算过程。让我们暂时放下那些令人头疼的数学符号用Python代码和可视化来直观感受这个强化学习的核心思想。1. 从网格世界开始建立强化学习直觉想象一个4×4的网格世界我们的小机器人可以在其中移动。某些格子有奖励比如5某些有惩罚比如-5大部分格子则是中性的0。机器人的目标是通过尝试不同的移动策略上、下、左、右学习如何最大化长期累积奖励。import numpy as np import matplotlib.pyplot as plt # 定义网格世界 grid_size 4 rewards np.zeros((grid_size, grid_size)) rewards[0, 3] 10 # 右上角是目标 rewards[1, 3] -10 # 陷阱 gamma 0.9 # 折扣因子 # 可视化初始网格 def plot_grid(values): fig, ax plt.subplots() im ax.imshow(values, cmapRdYlGn) for i in range(grid_size): for j in range(grid_size): text ax.text(j, i, f{values[i, j]:.1f}, hacenter, vacenter, colorblack) ax.set_title(网格世界状态值分布) plt.colorbar(im) plt.show() plot_grid(rewards)这个简单的网格世界已经包含了强化学习的几个关键要素状态(State): 每个网格位置动作(Action): 移动方向奖励(Reward): 进入特定网格的即时反馈策略(Policy): 决定在每个状态下采取什么动作的规则2. 状态价值的计算从直觉到公式状态价值(State Value)是贝尔曼方程的核心概念表示从某个状态出发按照当前策略能获得的预期累积奖励。让我们先用最直观的方式计算# 初始化状态值为0 values np.zeros((grid_size, grid_size)) # 简单计算相邻状态的影响 for i in range(grid_size): for j in range(grid_size): # 考虑四个方向的移动上、下、左、右 neighbors [] if i 0: neighbors.append(values[i-1, j]) # 上 if i grid_size-1: neighbors.append(values[i1, j]) # 下 if j 0: neighbors.append(values[i, j-1]) # 左 if j grid_size-1: neighbors.append(values[i, j1]) # 右 if neighbors: values[i, j] rewards[i, j] gamma * np.mean(neighbors) plot_grid(values)这个简单的计算已经体现了贝尔曼方程的核心思想当前状态的价值 即时奖励 折扣后的未来状态价值。但你会发现这种单次计算并不准确因为状态价值之间相互依赖。3. 迭代求解实现贝尔曼方程的完整形式贝尔曼方程的完整形式考虑了所有可能的动作和状态转移概率。让我们实现一个更精确的版本# 定义策略随机选择动作上、下、左、右概率均等 policy np.ones((grid_size, grid_size, 4)) / 4 # 每个状态有4个动作 # 迭代求解贝尔曼方程 def bellman_update(values): new_values np.zeros_like(values) for i in range(grid_size): for j in range(grid_size): if (i, j) (0, 3) or (i, j) (1, 3): # 终止状态 new_values[i, j] rewards[i, j] continue total 0 # 四个动作上(0)、下(1)、左(2)、右(3) for action in range(4): # 计算执行该动作后的新位置 new_i, new_j i, j if action 0 and i 0: new_i i - 1 elif action 1 and i grid_size - 1: new_i i 1 elif action 2 and j 0: new_j j - 1 elif action 3 and j grid_size - 1: new_j j 1 # 贝尔曼方程计算 total policy[i, j, action] * (rewards[i, j] gamma * values[new_i, new_j]) new_values[i, j] total return new_values # 迭代计算直到收敛 for _ in range(100): values bellman_update(values) plot_grid(values)这个实现展示了贝尔曼方程的几个关键点Bootstrapping: 使用当前估计的值函数来更新自身期望计算: 对所有可能的动作和状态转移取期望折扣因子γ: 平衡即时奖励和未来奖励的重要性4. 价值迭代从状态价值到最优策略有了准确的状态价值我们就能改进策略。价值迭代算法就是不断交替进行策略评估和策略改进# 价值迭代实现 def value_iteration(max_iter100, theta1e-4): values np.zeros((grid_size, grid_size)) for _ in range(max_iter): delta 0 new_values np.zeros_like(values) for i in range(grid_size): for j in range(grid_size): if (i, j) (0, 3) or (i, j) (1, 3): # 终止状态 new_values[i, j] rewards[i, j] continue # 找出能带来最大价值的动作 max_value -np.inf for action in range(4): new_i, new_j i, j if action 0 and i 0: new_i i - 1 elif action 1 and i grid_size - 1: new_i i 1 elif action 2 and j 0: new_j j - 1 elif action 3 and j grid_size - 1: new_j j 1 current_value rewards[i, j] gamma * values[new_i, new_j] if current_value max_value: max_value current_value new_values[i, j] max_value delta max(delta, abs(new_values[i, j] - values[i, j])) values new_values if delta theta: break # 从最终的价值函数推导最优策略 optimal_policy np.zeros((grid_size, grid_size), dtypeint) for i in range(grid_size): for j in range(grid_size): if (i, j) (0, 3) or (i, j) (1, 3): continue best_action 0 max_value -np.inf for action in range(4): new_i, new_j i, j if action 0 and i 0: new_i i - 1 elif action 1 and i grid_size - 1: new_i i 1 elif action 2 and j 0: new_j j - 1 elif action 3 and j grid_size - 1: new_j j 1 current_value values[new_i, new_j] if current_value max_value: max_value current_value best_action action optimal_policy[i, j] best_action return values, optimal_policy optimal_values, optimal_policy value_iteration() plot_grid(optimal_values) # 可视化最优策略 action_symbols [↑, ↓, ←, →] policy_grid np.array([[action_symbols[a] for a in row] for row in optimal_policy]) policy_grid[0, 3] G # 目标 policy_grid[1, 3] T # 陷阱 print(最优策略) print(policy_grid)这个实现展示了如何通过最大化状态价值来改进策略从收敛的状态价值函数中提取最优策略可视化学习到的最优移动方向5. 贝尔曼方程的数学本质与代码对应现在让我们将数学公式与代码实现对应起来加深理解贝尔曼方程的标准形式v(s) E[R γv(s)|s]我们的代码实现对应关系数学概念代码实现状态s(i, j)网格坐标动作aaction(0-3对应上下左右)奖励Rrewards[i, j]折扣因子γgamma 0.9状态转移概率由policy数组表示期望E[...]对动作和状态的循环与加权平均贝尔曼最优方程v*(s) max_a E[R γv*(s)|s,a]价值迭代的核心代码正是实现了这个方程max_value -np.inf for action in range(4): # 计算状态转移 current_value rewards[i, j] gamma * values[new_i, new_j] if current_value max_value: max_value current_value new_values[i, j] max_value6. 高级话题扩展与优化理解了基础实现后我们可以考虑几个优化方向动态规划的高效实现# 使用向量化操作加速计算 def vectorized_bellman_update(values): new_values np.zeros_like(values) # 上移 up np.roll(values, -1, axis0) up[-1, :] values[-1, :] # 下移 down np.roll(values, 1, axis0) down[0, :] values[0, :] # 左移 left np.roll(values, -1, axis1) left[:, -1] values[:, -1] # 右移 right np.roll(values, 1, axis1) right[:, 0] values[:, 0] # 计算四个方向的最大值 new_values rewards gamma * np.maximum.reduce([up, down, left, right]) # 固定终止状态的值 new_values[0, 3] rewards[0, 3] new_values[1, 3] rewards[1, 3] return new_values异步动态规划 不同于同步更新所有状态异步方法每次只更新部分状态可以更快收敛def async_value_iteration(): values np.zeros((grid_size, grid_size)) state_order [(i, j) for i in range(grid_size) for j in range(grid_size)] for _ in range(100): np.random.shuffle(state_order) for i, j in state_order: if (i, j) in [(0, 3), (1, 3)]: continue max_value -np.inf for action in range(4): new_i, new_j i, j if action 0 and i 0: new_i i - 1 elif action 1 and i grid_size - 1: new_i i 1 elif action 2 and j 0: new_j j - 1 elif action 3 and j grid_size - 1: new_j j 1 current_value rewards[i, j] gamma * values[new_i, new_j] if current_value max_value: max_value current_value values[i, j] max_value return values处理随机环境 现实中的状态转移往往不是确定的。我们可以修改代码考虑转移概率# 定义转移概率有80%概率按预期移动20%概率随机移动 def stochastic_bellman_update(values): new_values np.zeros_like(values) for i in range(grid_size): for j in range(grid_size): if (i, j) in [(0, 3), (1, 3)]: new_values[i, j] rewards[i, j] continue total 0 for intended_action in range(4): # 80%概率执行预期动作 actual_actions [intended_action] * 8 [a for a in range(4) if a ! intended_action] * 2 for action in actual_actions: new_i, new_j i, j if action 0 and i 0: new_i i - 1 elif action 1 and i grid_size - 1: new_i i 1 elif action 2 and j 0: new_j j - 1 elif action 3 and j grid_size - 1: new_j j 1 prob 0.8 if action intended_action else 0.2/3 total prob * (rewards[i, j] gamma * values[new_i, new_j]) new_values[i, j] total / len(actual_actions) return new_values在实现这些扩展时我发现随机环境的策略收敛速度明显慢于确定环境这解释了为什么在实际强化学习问题中我们常常需要更复杂的算法如Q-learning或策略梯度方法。