""" ============================================================================= 第 5 讲 · 代码 2：四种优化器 —— 梯度上升、牛顿法、MM 算法、Adam ============================================================================= 对应讲义：lectures/lecture05_mle.md 第 4 节 ----------------------------------------------------------------------------- 【这个文件要解决什么问题？】 ----------------------------------------------------------------------------- BT 的对数似然是凹的，只有一个最高点（代码 01）。怎样最快地爬上去？ 演示 1：同一个问题（100 个选项、5000 条随机配对的比较），四种优化器各走多少步、花多少时间， 目标函数离最优值还差多少（"最优性差距"）随步数的变化。 演示 2：步长的意义 —— 梯度上升的步长取 1/L（L 是光滑度上界）一定收敛；取太大就会来回震荡甚至发散。 演示 3：条件数 —— 换成"链式"的比较设计（只和邻居比），梯度上升慢了十几倍，牛顿法几乎不受影响。 ----------------------------------------------------------------------------- 【运行方式】 ----------------------------------------------------------------------------- python3 02_optimizers.py # < 1 秒 需要 numpy。只依赖本目录的 bt_numpy.py。 """ import time import numpy as np from bt_numpy import (fit_adam, fit_gd, fit_mm, fit_newton, objective, pad, print_header, random_pairs, simulate_bt, smoothness_bound) def chain_pairs(M, n, rng): """“链式”设计：每个选项只和编号相差 1 或 2 的选项比较。""" i = rng.integers(0, M, size=n) step = rng.choice(np.array([-2, -1, 1, 2]), size=n) j = np.clip(i + step, 0, M - 1) keep = i != j return i[keep], j[keep] def run_all(i, j, y, M, tol=1e-6): """ 用四种优化器拟合同一份数据（都不加正则化），记录步数、时间和每一步的目标函数值。 【返回】dict {名字: (分数, 步数, 秒数, 目标函数值的历史)} """ out = {} for name, fn in [("梯度上升（步长 1/L）", lambda h: fit_gd(i, j, y, M, tol=tol, max_iters=200000, history=h)), ("牛顿法", lambda h: fit_newton(i, j, y, M, tol=tol, history=h)), ("MM 算法", lambda h: fit_mm(i, j, y, M, tol=tol, max_iters=200000, history=h)), ("Adam（lr = 0.05）", lambda h: fit_adam(i, j, y, M, lr=0.05, tol=tol, max_iters=20000, history=h))]: hist = [] t0 = time.perf_counter() s, steps = fn(hist) out[name] = (s, steps, time.perf_counter() - t0, hist) return out def report(results, best): """打印每种优化器的步数、时间和若干步时的最优性差距。""" checkpoints = [1, 10, 100, 1000, 10000] head = " ".join(pad(f"第{c}步", 9, "right") for c in checkpoints) print(f" {pad('优化器', 22)} {pad('收敛步数', 9, 'right')} {pad('耗时(秒)', 9, 'right')} 最优性差距：{head}") for name, (s, steps, secs, hist) in results.items(): gaps = [] for c in checkpoints: gaps.append(f"{best - hist[c - 1]:>9.1e}" if c - 1 < len(hist) else f"{'—':>9s}") steps_text = f"{steps:>9d}" if steps < 200000 and not (name.startswith("Adam") and steps >= 20000) else f"{'>' + str(steps):>9s}" print(f" {pad(name, 22)} {steps_text} {secs:>9.2f} {'':12s}{' '.join(gaps)}") # ============================================================================= # 演示 1：四种优化器 # ============================================================================= def demo_compare(rng): print_header("演示 1：同一个问题，四种优化器（100 个选项、5000 条随机配对的比较）") M = 100 true_s = rng.normal(0, 1, size=M) i, j, y = simulate_bt(true_s, 5000, rng) results = run_all(i, j, y, M) best = max(objective(r[0], i, j, y) for r in results.values()) print("收敛的标准：梯度的最大分量 < 1e-6。“最优性差距” = 最优的平均对数似然 − 当前的平均对数似然（越小越好）；") print("“—”表示在那一步之前已经收敛了。\n") report(results, best) s_ref = results["牛顿法"][0] print("\n 各优化器的解和牛顿法的解的最大差别：" + "，".join(f"{name.split('（')[0]} {np.max(np.abs(r[0] - s_ref)):.1e}" for name, r in results.items() if name != "牛顿法")) print(""" 观察： · 牛顿法只用了几步：在最优点附近，误差每一步大约平方一次（二次收敛）。代价是每一步要解一个 100 × 100 的 线性方程组；选项有几百万个时（例如给每条推文打分），这一步就太贵了。 · 梯度上升每步很便宜，但要走几十步：误差每一步只按一个固定的比例缩小（线性收敛）。 · MM 算法（Zermelo 1929 年就用它给棋手排名）不需要选步长，每一步都保证似然不下降，步数和梯度上升差不多。 · Adam 的步数最多：它的步长是固定的“绝对大小”（lr），在这种光滑的凸问题上没有优势，靠近最优点时还会小幅震荡。 它的优势在神经网络这种大规模、随机梯度、非凸的问题（第 6 讲）。 · “收敛”到梯度 < 1e-6 时，各方法的解和牛顿法仍差 1e-3 左右：梯度很小不等于参数很准，平缓的方向上 梯度小、参数却还能差一截（演示 3 会看到更极端的情形）。 · 所有方法收敛到同一个解（凹函数只有一个最高点）—— 选哪个优化器只影响速度，不影响答案。""") # ============================================================================= # 演示 2：步长 # ============================================================================= def demo_step_size(rng): print_header("演示 2：梯度上升的步长 —— 1/L 一定收敛，太大就震荡或发散") M = 100 true_s = rng.normal(0, 1, size=M) i, j, y = simulate_bt(true_s, 5000, rng) L = smoothness_bound(i, j, M, len(y), 0.0) s_star, _ = fit_newton(i, j, y, M) best = objective(s_star, i, j, y) print(f"光滑度上界 L = {L:.4f}，“安全”步长 1/L = {1 / L:.1f}。\n") print(f" {pad('步长', 16)} {pad('200 步后的最优性差距', 22, 'right')} {pad('结果', 6)}") for mult in [0.5, 1.0, 4.0, 16.0, 64.0]: lr = mult / L hist = [] s, _ = fit_gd(i, j, y, M, lr=lr, max_iters=200, tol=0.0, history=hist) last = hist[-1] if np.isfinite(hist[-1]) else float("nan") gap = best - last verdict = "收敛" if gap < 1e-3 else ("发散" if not np.isfinite(gap) or gap > 0.5 else "震荡 / 慢") print(f" {pad(f'{mult:g} × (1/L)', 16)} {gap:>22.2e} {verdict}") print(""" 观察： · 1/L 是“保证不会变差”的步长：L 是曲率的上界，按最陡的方向估计，一步不会跨过山顶。 · 1/L 是一个保守的上界，实际的曲率往往小得多，所以稍大的步长（几倍）也能收敛，甚至更快。 · 步长太大就会一步跨过山顶，来回震荡，甚至越跳越远（发散）。深度学习里“学习率太大，损失变成 NaN”就是这个现象。 · 这就是为什么第 1~4 讲的拟合函数都用 1/L 自动设定步长：不用调参，也保证收敛。""") # ============================================================================= # 演示 3：条件数 # ============================================================================= def demo_conditioning(rng): print_header("演示 3：比较设计决定了问题有多“难爬” —— 链式设计") M = 100 true_s = np.linspace(-2, 2, M) for name, pairs in [("随机配对", random_pairs(M, 5000, rng)), ("链式（只和邻居比）", chain_pairs(M, 5000, rng))]: i, j, y = simulate_bt(true_s, len(pairs[0]), rng, pairs=pairs) res = {} for opt, fn in [("梯度上升", lambda: fit_gd(i, j, y, M, tol=1e-6, max_iters=200000)), ("牛顿法", lambda: fit_newton(i, j, y, M, tol=1e-6))]: t0 = time.perf_counter() _, steps = fn() res[opt] = (steps, time.perf_counter() - t0) print(f" {pad(name, 20)} 梯度上升 {res['梯度上升'][0]:>6d} 步（{res['梯度上升'][1]:.2f} 秒）；" f"牛顿法 {res['牛顿法'][0]:>3d} 步（{res['牛顿法'][1]:.2f} 秒）") print(""" 观察： · 链式设计里，相隔很远的选项只能通过一长串中间人间接比较，海森矩阵的最小非零特征值（代数连通度，代码 01）很小， 最大特征值却不变：两者之比（“条件数”）很大。梯度上升的收敛速度由条件数决定 —— 最平缓的方向上每步只挪一点点， 要很多步才走到头。 · 牛顿法用海森矩阵的逆“矫正”了各个方向的曲率差异，几乎不受条件数影响。 · 实践意义：比较数据的**设计**不只影响统计精度（代码 01 演示 3），也影响计算。Chatbot Arena 这样的平台 让新模型多和各种水平的模型对战，两方面都受益。""") def main(): rng = np.random.default_rng(1929) demo_compare(rng) demo_step_size(rng) demo_conditioning(rng) if __name__ == "__main__": main()