""" ============================================================================= 第 3 讲 · 代码 1：蒙特卡洛策略评估（Monte Carlo Policy Evaluation） ============================================================================= 对应讲义：lectures/lecture03_model_free_evaluation.md 第 1、2 节 ----------------------------------------------------------------------------- 【这个文件要解决什么问题？】 ----------------------------------------------------------------------------- 第 2 讲的策略评估需要知道模型 P 和 R：V^π = R^π + γ P^π V^π。 现在我们**不知道模型**，只能让智能体按策略 π 和环境交互，收集一条条轨迹（回合）。 蒙特卡洛（MC）的想法极其朴素： 价值的定义是 V^π(s) = E[G_t | s_t = s]（回报的期望） 期望算不出来？那就把实际观察到的回报**取平均**。 ----------------------------------------------------------------------------- 【本文件使用的环境：随机游走（Sutton & Barto 例 6.2）】 ----------------------------------------------------------------------------- 终点(0) ← A ⇄ B ⇄ C ⇄ D ⇄ E → 终点(+1) ▲ 起点 - 5 个非终止状态 A~E（代码中编号 0~4），从 C 出发 - 每一步以 0.5 的概率向左或向右（这就是被评估的"策略 + 环境"的整体效果） - 从 E 向右走出去获得奖励 +1，其他所有转移奖励为 0；到达两端即结束 - γ = 1（回合一定会结束，不打折也没问题） 为什么用它？因为**真实价值可以精确算出**：V(A..E) = 1/6, 2/6, 3/6, 4/6, 5/6 （V(s) 恰好等于"从 s 出发最终从右边出去的概率"）。有了标准答案，才能精确衡量估计误差。 ----------------------------------------------------------------------------- 【数据格式】 ----------------------------------------------------------------------------- 一个回合（episode）是一个列表，每个元素是一次转移 (s, r, s_next)： s : 当前状态（整数） r : 这一步获得的奖励 r_t s_next : 下一个状态；如果回合在这一步结束，则为 None 例：[(2, 0, 3), (3, 0, 4), (4, 1, None)] 表示 C→D→E→从右边出去，获得 +1。 ----------------------------------------------------------------------------- 【运行方式】 ----------------------------------------------------------------------------- python 01_monte_carlo_evaluation.py """ import numpy as np np.set_printoptions(precision=3, suppress=True) STATE_NAMES = ["A", "B", "C", "D", "E"] # ============================================================================= # 第 1 部分：环境——生成随机游走的回合 # ============================================================================= def generate_random_walk_episode(rng, n_states=5, start=None, max_steps=10_000): """ 生成一条随机游走的完整回合。 【参数】 rng : numpy 随机数生成器。 n_states : int，非终止状态个数（默认 5，即 A~E）。 start : int 或 None，起点；None 表示中间状态 n_states // 2。 【Python 小知识】// 是整除：5 // 2 = 2。 max_steps : int，最多走多少步。随机游走理论上一定会结束，但加一个上限是好习惯： 万一你以后改了环境让它可能永远走不完，程序不会卡死。 【返回】 list，每个元素是 (s, r, s_next)，s_next 为 None 表示回合结束。 【边界】 - 走满 max_steps 还没结束 → RuntimeError（MC 必须要有完整的回合才能计算回报）。 """ s = n_states // 2 if start is None else start episode = [] for _ in range(max_steps): step = 1 if rng.random() < 0.5 else -1 s_next = s + step if s_next < 0: # 从左边走出去 episode.append((s, 0.0, None)) return episode if s_next >= n_states: # 从右边走出去，奖励 +1 episode.append((s, 1.0, None)) return episode episode.append((s, 0.0, s_next)) s = s_next raise RuntimeError(f"回合在 {max_steps} 步内没有结束") def true_random_walk_values(n_states=5): """ 随机游走的真实价值：V(第 i 个状态) = (i+1) / (n_states+1)。 【为什么？】γ=1、只有右出口有 +1，所以 V(s) = "从 s 出发最终从右边出去的概率"。 这是经典的"赌徒破产"问题，答案是离左端的距离除以总长度。 """ return np.arange(1, n_states + 1) / (n_states + 1) # ============================================================================= # 第 2 部分：计算回报 # ============================================================================= def returns_of_episode(episode, gamma): """ 计算一个回合中每个时刻的回报 G_t（从后往前扫，G_t = r_t + γ G_{t+1}，见第 1 讲 8.3 节）。 【参数】 episode : list of (s, r, s_next)。 gamma : float，[0, 1]。 【返回】 numpy 数组，长度 = 回合步数，第 t 个元素是 G_t。 【边界】空回合 → 返回空数组。 """ G = 0.0 returns = np.zeros(len(episode)) for t in range(len(episode) - 1, -1, -1): _, r, _ = episode[t] G = r + gamma * G returns[t] = G return returns # ============================================================================= # 第 3 部分：首次访问 MC 与每次访问 MC # ============================================================================= def first_visit_mc(episodes, n_states, gamma): """ 首次访问蒙特卡洛（First-visit MC）策略评估。 【思路】 对每个回合，找出每个状态**第一次**出现的时刻 t，记下那个时刻的回报 G_t。 最后，每个状态的价值估计 = 它所有"首次访问回报"的平均值。 为什么只用第一次？因为同一回合内对同一状态的多次访问，回报是高度相关的 （后一次的回报是前一次回报的"尾巴"）。只用第一次，每个回合对每个状态恰好贡献 一个样本，不同回合之间独立同分布，于是平均值是真实价值的**无偏**估计。 【参数】 episodes : list，每个元素是一个回合。 n_states : int，状态总数。 gamma : float，折扣因子。 【返回】 (V, counts)： - V : 长度 n_states 的价值估计；**从未被访问的状态为 NaN**（"不知道"） - counts : 每个状态被计入的样本数 【边界】 - 从未访问过的状态没有任何数据，返回 NaN 而不是 0： 0 是一个具体的数值猜测，会让人误以为"估计出来价值是 0"。 NaN 明确表示"没有信息"。（np.nanmean 等函数可以跳过 NaN。） """ total = np.zeros(n_states) counts = np.zeros(n_states, dtype=int) for episode in episodes: G = returns_of_episode(episode, gamma) seen = set() # 本回合已经访问过的状态 for t, (s, _, _) in enumerate(episode): if s in seen: continue # 不是首次访问，跳过 seen.add(s) total[s] += G[t] counts[s] += 1 # 【Python 小知识】np.where(条件, a, b)：条件为真的位置取 a，否则取 b。 # np.maximum(counts, 1) 防止除以 0（没访问过的位置反正会被替换成 NaN）。 V = np.where(counts > 0, total / np.maximum(counts, 1), np.nan) return V, counts def every_visit_mc(episodes, n_states, gamma): """ 每次访问蒙特卡洛（Every-visit MC）：状态在一个回合中每出现一次，就记一个回报样本。 【与首次访问的区别】 - 数据利用更充分（样本更多）； - 但同一回合内的样本相关，估计在样本有限时是**有偏**的（演示 3 会看到）， 不过随着回合数增加，偏差趋于 0（**一致**的）。 【返回】同 first_visit_mc。 """ total = np.zeros(n_states) counts = np.zeros(n_states, dtype=int) for episode in episodes: G = returns_of_episode(episode, gamma) for t, (s, _, _) in enumerate(episode): total[s] += G[t] counts[s] += 1 V = np.where(counts > 0, total / np.maximum(counts, 1), np.nan) return V, counts # ============================================================================= # 第 4 部分：增量式 MC # ============================================================================= def incremental_mc(episodes, n_states, gamma, alpha=None, V_init=0.0): """ 增量式（首次访问）MC：每看完一个回合就更新一次，不需要存储所有历史回报。 【更新公式】 V(s) ← V(s) + α · (G_t - V(s)) - alpha=None：α = 1/N(s)，N(s) 是 s 目前的样本数。此时结果**恰好等于**样本平均 （第 1 讲 04 文件的 incremental_mean_update 推导过）。 - alpha=常数（如 0.05）：越新的回报权重越大，旧回报的影响按 (1-α) 指数衰减。 好处：如果环境/策略在变化（非平稳），它能跟上变化； 代价：估计会一直在真实值附近抖动，不会完全收敛（因为每次都被最新样本"拉"一下）。 【括号里的 (G_t - V(s)) 是什么？】 "实际拿到的回报"减去"我原来的预期"——一种"惊讶程度"。 比预期好就调高估计，比预期差就调低。后面 TD、Q-learning 全是这个形式。 【参数】 alpha : None 或 (0, 1] 之间的常数。 V_init : float，初始估计。对 α=1/N 没有影响（第一个样本会完全覆盖它）； 对常数 α 有影响（需要一段时间才能"忘掉"初始值）。 【返回】 (V, history)：最终估计，以及每个回合结束后的估计（形状 (回合数, n_states)）。 【边界】 - alpha 不在 (0, 1] → ValueError。α > 1 会"矫枉过正"，估计会发散振荡。 """ if alpha is not None and not (0.0 < alpha <= 1.0): raise ValueError(f"alpha 必须在 (0, 1] 内或为 None，但收到了 {alpha}") V = np.full(n_states, float(V_init)) N = np.zeros(n_states, dtype=int) history = [] for episode in episodes: G = returns_of_episode(episode, gamma) seen = set() for t, (s, _, _) in enumerate(episode): if s in seen: continue seen.add(s) N[s] += 1 step = 1.0 / N[s] if alpha is None else alpha V[s] += step * (G[t] - V[s]) history.append(V.copy()) return V, np.array(history) # ============================================================================= # 第 5 部分：演示 # ============================================================================= def main(): # --------------------------------------------------------------------- print("=" * 78) print("演示 1：手算例子（复现讲义 2.4 节：火星车，一直向左，γ = 1）") print("=" * 78) # 火星车 7 个状态 s1..s7（代码 0..6），R = [1,0,0,0,0,0,10]，在 s1、s7 执行动作后回合结束。 # 一条轨迹：s3 → s2 → s2 → s1 → 结束 mars_episode = [(2, 0.0, 1), (1, 0.0, 1), (1, 0.0, 0), (0, 1.0, None)] print(" 轨迹：s3 --0--> s2 --0--> s2 --0--> s1 --1--> 结束") print(f" 各时刻的回报 G_t = {returns_of_episode(mars_episode, 1.0)}") V_fv, _ = first_visit_mc([mars_episode], 7, gamma=1.0) V_ev, cnt_ev = every_visit_mc([mars_episode], 7, gamma=1.0) print(f" 首次访问 MC：V = {V_fv} （nan = 没访问过，没有信息）") print(f" 每次访问 MC：V = {V_ev} （s2 被访问了 {cnt_ev[1]} 次，两个回报都是 1）") print(" → 一条轨迹就把终点的 +1 \"一路传回\"了 s3：MC 用的是完整的实际回报。") # --------------------------------------------------------------------- print("\n" + "=" * 78) print("演示 2：随机游走上，MC 估计随回合数增加而逼近真实值") print("=" * 78) rng = np.random.default_rng(0) V_true = true_random_walk_values() print(f" 真实价值 A~E = {V_true}") all_episodes = [generate_random_walk_episode(rng) for _ in range(10_000)] print(f" 一条回合的样子：{[(STATE_NAMES[s], r, None if sn is None else STATE_NAMES[sn]) for s, r, sn in all_episodes[0]]}") for n_ep in [10, 100, 1000, 10_000]: V, counts = first_visit_mc(all_episodes[:n_ep], 5, gamma=1.0) rms = np.sqrt(np.nanmean((V - V_true) ** 2)) print(f" {n_ep:6d} 个回合：V = {V} RMS 误差 = {rms:.4f}") print(" → 误差大体上按 1/sqrt(回合数) 缩小（回合数 ×100，误差约 ÷10）；") print(" 单次实验有随机波动，所以不会严格符合这个比例。") # --------------------------------------------------------------------- print("\n" + "=" * 78) print("演示 3：首次访问无偏，每次访问有偏（但一致）") print("=" * 78) # 一个只有 1 个状态的 MRP：每步奖励 +1，以概率 p 留下、1-p 结束，γ = 1。真实 V = 1/(1-p)。 p = 0.5 print(f" 单状态 MRP：每步 +1，以 {p} 的概率继续、{1 - p} 的概率结束。真实 V = 1/(1-p) = {1 / (1 - p)}") def self_loop_episode(rng): episode = [] while True: if rng.random() < p: episode.append((0, 1.0, 0)) else: episode.append((0, 1.0, None)) return episode rng = np.random.default_rng(1) for n_ep in [1, 10, 100]: fv_estimates, ev_estimates = [], [] for _ in range(20_000): # 重复 20000 次独立实验，看估计的"平均值" eps = [self_loop_episode(rng) for _ in range(n_ep)] fv_estimates.append(first_visit_mc(eps, 1, 1.0)[0][0]) ev_estimates.append(every_visit_mc(eps, 1, 1.0)[0][0]) print(f" 每次实验用 {n_ep:3d} 个回合：首次访问估计的均值 = {np.mean(fv_estimates):.3f}，" f"每次访问估计的均值 = {np.mean(ev_estimates):.3f}") print(" → 首次访问的平均值始终 ≈ 2（无偏）；每次访问在回合少时明显偏低（有偏），") print(" 回合多了才逐渐接近 2（一致）。讲义 2.5 节解释了原因：只有 1 个回合时，") print(" 长度为 L 的回合给出回报 L, L-1, ..., 1，每次访问取它们的平均 (L+1)/2，期望为 1.5。") # --------------------------------------------------------------------- print("\n" + "=" * 78) print("演示 4：增量式 MC——α = 1/N（精确平均） vs 常数 α") print("=" * 78) episodes = all_episodes[:2000] V_avg, hist_avg = incremental_mc(episodes, 5, 1.0, alpha=None, V_init=0.5) V_batch, _ = first_visit_mc(episodes, 5, 1.0) print(f" α = 1/N 的结果与\"存下所有回报再求平均\"是否一致？{np.allclose(V_avg, V_batch)}") for alpha in [None, 0.01, 0.05, 0.2]: _, hist = incremental_mc(episodes, 5, 1.0, alpha=alpha, V_init=0.5) # 最后 500 个回合中 V(C) 的波动（标准差） tail = hist[-500:, 2] err = np.sqrt(np.mean((hist[-1] - V_true) ** 2)) name = "1/N " if alpha is None else f"{alpha:<4}" print(f" α = {name}：最终 RMS 误差 {err:.4f}；最后 500 回合中 V(C) 的波动（标准差）{tail.std():.4f}") print(" → 常数 α 越大，估计抖动越厉害（永远被最新的样本拉来拉去）；α = 1/N 越来越稳。") print(" 但常数 α 也有好处：如果环境在变化，它能\"忘掉\"过时的数据。") # --------------------------------------------------------------------- print("\n" + "=" * 78) print("演示 5：边界情况") print("=" * 78) # 找一个没有走遍所有状态的短回合 # 【Python 小知识】next(生成器表达式) 返回第一个满足条件的元素。 short_ep = next(ep for ep in all_episodes if len({s for s, _, _ in ep}) < 5) V_short, cnt_short = first_visit_mc([short_ep], 5, 1.0) print(f" [只有 1 个短回合] 回合 = {[STATE_NAMES[s] for s, _, _ in short_ep]}") print(f" V = {V_short}，访问次数 = {cnt_short}（没到过的状态是 nan）") try: generate_random_walk_episode(np.random.default_rng(0), max_steps=2) except RuntimeError as e: print(f" [回合太长] RuntimeError: {e}（MC 必须等回合结束才能算回报）") try: incremental_mc(all_episodes[:10], 5, 1.0, alpha=1.5) except ValueError as e: print(f" [α 越界] ValueError: {e}") if __name__ == "__main__": main()