MDP策略迭代在MATLAB中的实现:从贝尔曼方程到最优策略

发布时间:2026/10/3 18:06:21
MDP策略迭代在MATLAB中的实现:从贝尔曼方程到最优策略 简介这份MDP.zip资源是一套基于MATLAB的马尔科夫决策过程MDP策略迭代算法实践包面向机器学习、强化学习入门者及相关课程学习者用来求解离散状态与动作空间下的最优策略并帮助完成从理论到代码的跨越。包体约2.32MB共93个文件以33个m函数脚本为核心覆盖策略迭代、值迭代、有限期MDP、线性规划求解等常用算法另有20个html函数文档、18个dat示例数据和8个txt说明文件便于对照帮助文档开展实验与二次开发。目前已有312人学习/下载适合作为学习策略评估与策略改进流程的动手素材。借助这套代码读者可以理解策略迭代两阶段交替收敛的过程掌握使用MDPtoolbox建模、校验与求解的方法并通过自带race示例和readme加深对MDP与强化学习的整体认识。1. MDP策略迭代在MATLAB里解决的不是一道题而是一类决策问题MDP马尔可夫决策过程配上策略迭代policy iteration是强化学习里最经典的一对组合。你手里这份“MDP.zip”如果按“mdp_policy iteration”组织那它要解决的事情很明确给你一个能用“状态-动作-转移-奖励”描述的决策环境通过反复的“评估-改进”循环输出一张最优策略表——每个状态下该选哪个动作。不管是西电机器学习期末、吴恩达机器学习作业里的迷宫题还是机器人路径规划、库存控制底层都是这套东西。这篇文章写给两类人要交策略迭代MATLAB作业的学生以及想在仿真里先验证算法再上真机的工程师。读完你能自己写、能调参、能避开那些让结果悄悄变错的坑。2. 建模先行把MDP五元组和策略迭代方程写进代码之前2.1 把决策问题翻译成(S,A,P,R,γ)动手写代码前先做一件让很多新手翻车的事建模。MDP建模的标准产物是五元组(S,A,P,R,γ)。S是状态集合有限MDP下直接编号成1到NA是动作集合编号成1到MP是状态转移函数通常存成N×N×M的三维矩阵P(i,j,a)表示在状态i选动作a后跳到状态j的概率R是奖励可以退化成只依赖(s,a)的二维表γ是折扣因子取值0到1之间越小越“短视”。这个五元组不是抽象的数学概念你后面所有函数都拿它喂数据。我写代码前一定会先把状态、动作、转移概率、奖励用表格列清楚再转成矩阵。表格列不清楚写出来的P矩阵一定错位而且这种错位极难靠调试发现——程序不报错结果却是错的。常见做法是把P存成三维矩阵第一个维度是当前状态第二个维度是下一状态第三个维度是动作。R可以按同样的方式存成N×M二维表。γ单独一个标量。这个结构MATLAB天然友好索引直观后面画热力图也方便。2.2 策略迭代核心评估一遍、改进一步、反复到不动点策略迭代的思想比名字朴素得多先随便给一个策略π0然后重复两件事。第一步“策略评估”解出当前策略下的状态价值函数V它满足贝尔曼方程V^π(s) R(s,π(s)) γ Σ_{s} P(s|s,π(s)) V^π(s)。第二步“策略改进”对每个状态重新算所有动作的Q值Q(s,a) R(s,a) γ Σ_{s} P(s|s,a) V(s)把策略更新为每个状态下argmax_a Q(s,a)。两步交替直到策略不再变化。收敛性不是玄学有限MDP、折扣因子γ1时策略改进保证每一轮的价值函数单调不减且策略总数有限所以一定在有限步到达最优策略π*。这也是你愿意把实验挂上去跑的理论底气——不会出现“越迭代越差”。但注意这个保证依赖γ1或存在吸收态后面避坑章会专门讲γ取1的事。2.3 一份策略迭代MATLAB包的常规分工拿到MDP.zip这类包我一般先看文件是不是拆成这样几块policy_evaluation.m负责解贝尔曼方程policy_improvement.m负责贪心更新策略mdp_policy_iteration.m或main.m负责把两者串起来跑主循环另外通常有一个脚本负责打印或画策略矩阵。如果缺了哪块你自己补就行——代码量不大核心三个函数加起来也就几十行。如果你手上这份包把一切都堆在一个脚本里建议你动手拆开。拆文件不是为了好看是为了调试评估阶段不收敛时你能单独跑policy_evaluation看V的变化曲线改进阶段出问题时能单独检查argmax有没有越界。一把梭的脚本最难查错。2.4 一个5状态链式问题把模型实例化下面拿一个5状态链式MDP做例子方便你拿笔和MATLAB同步验算。状态1和5是终止状态吸收态进入后奖励0并停止状态2、3、4可选动作“左”或“右”。选“右”时以0.8概率走到右边状态、0.2概率滑到左边状态选“左”时对称。到达状态5给1奖励撞回状态1给0中间转移奖励0γ取0.9。状态编号直接1到5动作编号1左、2右。这个例子规模小最优策略一眼能看出来状态2、3、4都选“右”就能最快冲向状态5。但滑移概率让策略评估很有内容足够检验代码正确性。后面所有代码都围绕这个模型展开。3. 用MATLAB手写策略迭代三个核心函数和一套可直接跑的示例3.1 策略评估迭代展开贝尔曼方程策略评估有两条路一是解线性方程组精确评估二是用迭代法逼近近似评估。对于刚接触的人来说迭代法更直观也更容易暴露问题。基本公式是反复执行 V_new(s) R(s,π(s)) γ Σ P(s|s,π(s)) V_old(s)直到两次迭代的差值小于阈值theta。这里我给出迭代式评估函数function [V, iter] policy_evaluation(P, R, gamma, policy, theta, max_iter) % 策略评估迭代求解贝尔曼方程 % P: N x N x M 状态转移矩阵 % R: N x M 奖励矩阵 % gamma: 折扣因子 % policy: N x 1 当前策略每个状态选的动作编号 % theta: 收敛阈值 % max_iter: 最大迭代次数 N size(P, 1); V zeros(N, 1); for iter 1:max_iter V_new zeros(N, 1); for s 1:N a policy(s); % 贝尔曼期望方程立即奖励 折扣 * 下一状态价值期望 V_new(s) R(s, a) gamma * sum(squeeze(P(s, :, a)) .* V); end if max(abs(V_new - V)) theta V V_new; break; end V V_new; end end逻辑说明每次迭代先用上一轮的V计算每个状态的新价值。squeeze(P(s,:,a))取出第s行第a页的N个转移概率然后和V做逐点乘法再求和就是ΣP(s|s,a)V(s)。终止判断用的是“最大绝对变化量”max(abs(delta))比平均变化更严格避免个别状态还没收敛就停下来。参数说明theta取1e-4到1e-6之间比较合适太小会让迭代次数爆炸max_iter设个500到1000的兜底防止gamma接近1时死循环。注意这里的policy是列向量动作编号必须落在1到M之间否则squeeze取出来的维度不对。3.2 策略改进贪心选择最优动作评估完当前策略的V之后要据此改进策略。改进的原则是对每个状态计算所有动作的Q值选Q值最大的那个动作作为新策略。这一步是策略迭代的“贪心”环节也是策略单调改进的引擎。function new_policy policy_improvement(P, R, gamma, V) % 策略改进按Q值贪心更新 % P: N x N x M % R: N x M % gamma: 折扣因子 % V: N x 1 当前价值函数 N size(P, 1); M size(P, 3); Q zeros(N, M); for s 1:N for a 1:M % Q(s,a) R(s,a) gamma * sum(P(s,s,a)*V(s)) Q(s, a) R(s, a) gamma * sum(squeeze(P(s, :, a)) .* V); end end [~, new_policy] max(Q, [], 2); end逻辑说明先构造一个N×M的Q表逐状态逐动作计算Q值。squeeze(P(s,:,a))取出转移概率向量后与V做点乘结果就是期望下一状态价值。max(Q,[],2)沿第二维取最大值返回每个状态的最优动作编号。这里用[~,idx]取索引而非值是MATLAB里处理argmax的标准写法。参数说明这个函数本身没有额外参数可调它的正确性完全取决于V是否可靠。如果你发现策略来回跳、不收敛先别怀疑改进函数回去看评估函数是否真的收敛了——这是策略迭代里最经典的“二八分”问题八成bug在评估阶段。3.3 主循环收敛判断与迭代记录有了评估和改进两个零件主循环就是把它们装起来。收敛判断不是看V的变化而是看策略本身有没有变化。策略是离散的所以“策略不再变化”是比“V不再变化”更硬性的停止条件。同时要记录每一轮的策略方便事后画图和分析。function [V, policy, history] mdp_policy_iteration(P, R, gamma, policy0, theta, max_iter) % MDP策略迭代主循环 % 返回收敛后的V、最优策略、以及每轮策略快照 N size(P, 1); policy policy0; history zeros(max_iter, N); for k 1:max_iter history(k, :) policy; % 先评估 [V, ~] policy_evaluation(P, R, gamma, policy, theta, 500); % 再改进 new_policy policy_improvement(P, R, gamma, V); if isequal(new_policy, policy) policy new_policy; history history(1:k, :); return; end policy new_policy; end end逻辑说明每轮先记录当前策略快照再评估、改进。isequal比较新旧策略是否完全一致一致说明到达最优策略直接返回。这里的history用来事后观察策略的演化路径——如果策略在几个版本间反复横跳说明评估精度不够或gamma设置有问题。参数说明主循环的max_iter其实很难触发因为策略迭代收敛通常在十几轮以内除非你的状态空间大到离谱或评估严重不精确。theta和max_iter传给评估函数的值要注意评估的max_iter设500是给复杂状态空间留的余量5状态问题实际上几步就收敛。3.4 一个可直接复现的最简示例把上面的函数放到同一个文件夹后写一个主脚本把5状态链式问题跑起来。这个脚本你要保存成main_example.m直接运行就能看到输出。% main_example.m % 5状态链式MDP1和5为吸收态2/3/4可选左或右 N 5; M 2; P zeros(N, N, M); % 动作1: 左。状态2/3/480%去左边20%滑到右边 for s 2:4 P(s, s-1, 1) 0.8; P(s, s1, 1) 0.2; end % 动作2: 右。状态2/3/480%去右边20%滑到左边 for s 2:4 P(s, s1, 2) 0.8; P(s, s-1, 2) 0.2; end % 吸收态任何动作都停在原地 P(1, 1, :) 1; P(5, 5, :) 1; R zeros(N, M); % 到状态5给1其余为0 R(4, 2) 0.8 * 1; % 状态4选右80%概率到达5 gamma 0.9; policy0 ones(N, 1); % 初始策略全选左 theta 1e-5; [V, policy, history] mdp_policy_iteration(P, R, gamma, policy0, theta, 50); disp(最优策略 (1左, 2右):); disp(policy); disp(价值函数:); disp(V);运行后你会看到收敛后的策略是[1,2,2,2,1]即状态2、3、4全部选右V值呈现“越靠近5越高”的单调形态。我建议你手动验算一次把V(5)0当作已知反推V(4)再反推V(3)看和程序输出是否吻合。这一步能帮你确认自己对贝尔曼方程的理解是否到位。参数说明这里R(4,2)0.8是有意设计的——在状态4选右有80%概率直接到状态5拿1所以期望奖励就是0.8。如果你把奖励设计为“到达5之后额外加1”那R的写法会变。奖励的定义直接决定最优策略是建模时最容易和队友吵起来的地方。4. 策略迭代避坑五个翻车现场与修复方案4.1 状态编号从1开始你的转移矩阵全部错位现象程序能跑、不报错但输出策略明显不合理比如某个状态选了动作后总是走到错误方向。检查逻辑看起来也没错。原因MATLAB没有0索引所有数组从1开始。如果你在纸上把状态标成0到4或者借鉴C代码改写成0起步编号那么P(s, :, a)取出来的永远是错的那一行结果就是“看起来在算实际上全错位”。解决统一所有状态编号从1到N转移矩阵、奖励矩阵、策略向量全部用同一套编号。如果一定要对应外部数据的0编号在建模阶段统一做一次“1”映射s_matlab s_external 1。别在程序里东一个1西一个-1迟早漏一处。我自己的习惯是建模表格里就写1到N彻底绕开这个问题。4.2 gamma取1且没有吸收态策略评估永远不收敛现象策略评估函数跑到max_iter上限被强制退出V值和预期差异巨大。策略改进后策略乱跳。原因γ1意味着未来奖励不打折扣在无吸收态的无限持续MDP里价值可能无穷大贝尔曼方程没有唯一解迭代自然不收敛。很多人学的时候直接把γ写成1感觉“少一个参数省事”结果翻车。解决要么显式建模吸收态状态1和5那样进入后停在原地要么把γ设到0.9到0.99之间。教材里的作业题通常会给出γ别自作主张改。如果是自己的仿真项目γ0.95是个常用起点任务越长期γ越要接近1但不要等于1。4.3 评估阈值theta设太小迭代次数爆炸现象5状态问题评估阶段跑了1000步还没收敛主循环几十秒不结束。调小theta反而更慢。原因theta设成1e-12这类“高精度”数值在双精度浮点下不但无益还会让迭代逼近一个永远达不到的精度。线性收敛的迭代法后期收敛速度就是慢阈值越严步数翻倍增长。解决theta取1e-4到1e-6就足够支持策略迭代的正确决策。因为策略改进只关心Q值的相对大小微小误差不会改变argmax结果。如果你想确认精度影响可以同一问题分别跑theta1e-3和1e-6对比策略是否一致——通常一致。4.4 转移概率没归一化Q值越迭代越大现象V值不断膨胀后期直接超出合理范围策略反复横跳。检查gamma0.9没错吸收态也建了。原因很多新手写滑移概率时漏掉一种情况比如状态2选“左”0.8概率到状态10.2概率滑到状态3但这两个概率加起来正好是1吗如果不是或者你额外写了留在原地的概率却忘了把它加进归一化那么ΣP(s|s,a)不为1贝尔曼方程变成带放大器的不动点V就会被不断放大。解决建模后立刻用assert检查归一化——对每个状态和动作sum(P(s,:,a))必须严格等于1浮点数容差内。for s 1:N for a 1:M assert(abs(sum(P(s, :, a)) - 1) 1e-10, ... [转移概率未归一化: s, num2str(s), a, num2str(a)]); end end这段代码放在主脚本开头能把你从“结果悄悄错”的泥潭里拉出来。血泪经验这类bug在结果图里很难一眼看出来但策略落地到真机上一定会出事。4.5 中文注释乱码与脚本路径带空格现象MATLAB 2023及更新版本里从网盘或GitHub拉下来的脚本打开后中文注释一堆乱码linxu/mac上运行还可能出现路径找不到的错误。原因旧版MATLAB默认GBK编码新版切到UTF-8代码文件编码不匹配就会乱码。脚本放在带空格的路径下某些版本对文件路径解析不完整导致脚本调脚本时找不到文件。解决统一用UTF-8保存代码文件MATLAB里设置Preferences General Encoding选择UTF-8。路径规范上文件夹名避免空格和中文用连字符或下划线替代。这两个问题看着小但真能卡你一下午——尤其是期末ddl前从网盘拖下来的包乱码注释让你完全看不懂原作者的意图。5. 验证与进阶从网格世界到路径规划的正确性检查5.1 用贝尔曼最优性方程校验收敛结果策略迭代跑完后别急着信结果。独立算一遍最优贝尔曼方程残差对每个状态验证 V(s) 是否等于 max_a [ R(s,a) γ Σ P(s|s,a) V(s) ]。残差小于1e-6说明代码正确如果残差大先怀疑评估精度再怀疑P矩阵有没有归一化。这个校验脚本应该和你主程序分开每次跑完都执行一次形成习惯。5.2 把策略和价值画到网格上策略表是一堆数字不直观。把V值用imagesc画成热力图把策略画成箭头一眼就能看出“最优策略是否符合直觉”。比如5状态链式问题热力图应该显示V从状态1到状态5单调上升箭头全部指向右。如果图像里某个箭头方向反了就回去查那个状态的转移概率和奖励定义。这一步在课程作业里也是加分项。5.3 把策略导出成.mat供下游使用最优策略最终要落地。我会把policy和V保存成policy_result.mat用类似save(policy_result.mat, policy, V)的方式。后续在Simulink或机器人路径规划模块里直接load读取按状态查表决策。如果你的状态空间很大V没必要保存只存policy即可一个N×1的整型数组内存开销极其有限。在Linux服务器上跑大规模实验的话注意脚本编码和路径问题和4.5节说的一致。我自己的收尾习惯是每次跑完策略迭代先看一眼策略是否有突然跳变的孤立状态——那种和邻居策略方向相反的状态十有八九是P矩阵某一行没归一化。这个习惯替我省下的排查时间远比当初写策略迭代代码的时间多。希望这篇笔记能帮你把这套MDP策略迭代在MATLAB里跑通、跑对也让你在作业和仿真里少走几个我走过的弯路。本文还有配套的精品资源点击获取