逆强化学习IRL最小实现：专家示范反推奖励与策略-育师

简介：围绕逆强化学习（IRL）搭建的示例代码包，面向具备Java基础、希望从演示轨迹反推奖励函数并理解IRL框架实现的学习者与开发者，可用于课程设计、论文复现和算法原理解析。包体以484个Java源文件为主体，配合XML与pom构建配置、若干依赖JAR以及许可证和说明文档，共497个文件、约2.29MB，导入IDE后可按模块阅读与调试。内容包含BURLAP代码库快照及为创建IRL框架所做的修改，并涉及网格世界、单阶段正规形博弈、学徒学习、性能绘图等实验场景，读者可据此理解环境建模、求解器接入、算法流程与结果评估的组织方式，把逆强化学习的抽象概念落到可运行代码上。已有1048人学习下载，适合需要对照源码梳理IRL实现细节、搭建实验基线的中高级学习者。

1. 从示范轨迹反推奖励：一套不依赖大框架的 IRL 最小实现

最尴尬的一次调参，是给一个路径规划任务手写奖励函数。每改一次代价系数，智能体就换一种绕开目标的方式，直到某个版本它学会了贴着墙走——奖励写错了，强化学习会非常诚实地把错误放大。后来我把思路反过来：不写奖励，只给几段专家示范，让算法自己把奖励反推出来。这套教程示例做的就是这件事，逆强化学习（Inverse Reinforcement Learning）。代码是能跑通的最小实现：一个网格环境、一个用真实奖励产出示范轨迹的专家、两套求解器（最大熵 IRL 与学徒学习），输入是示范轨迹，输出是一组奖励权重以及由它导出的策略。适合已经懂值迭代、Q-learning 这些基础概念、想把模仿学习真正落到项目里的人，也适合带学生做课题时找一个不依赖大型框架的起点——跑一遍只要几分钟，所有中间量都能打印出来看。

2. 项目骨架拆解：环境层、专家层与求解器层的三段式结构

这套代码最值得学的地方不是算法本身，而是它把「环境」「示范数据」「求解器」切得干干净净。三者之间只通过几个数组约定接口，任何一层都可以单独替换：换一个更复杂的环境，求解器不用改；换一个求解器，环境不用改。工程代码看多了会发现，能不能被替换，往往比跑得多快更重要。

2.1 目录骨架：从 envs 到 runners 的分层

先把目录铺开看，结构本身就是一份设计说明。

irl_tutorial/ ├── envs/ │ ├── __init__.py │ └── gridworld.py # 网格环境：动力学、特征、可视化 ├── experts/ │ └── expert_policy.py # 在真实奖励上做值迭代并采样示范轨迹 ├── irl/ │ ├── base.py # IRLSolver 抽象：fit / reward / feature_expectation │ ├── maxent.py # 最大熵 IRL │ ├── apprentice.py # 学徒学习 │ └── utils.py # logsumexp、行归一化、特征期望统计 ├── runners/ │ ├── run_maxent.py # 命令行入口，读配置、训练、落盘 │ └── run_apprentice.py ├── configs/ │ └── grid_5x5.yaml # 环境参数与训练超参 ├── outputs/ # 权重文件、热力图、收敛曲线 └── tests/ └── test_transition.py # 断言转移矩阵每行和为 1

依赖不复杂，numpy承担全部数值计算，scipy用来做 log-sum-exp 和线性规划，matplotlib负责把权重画成热力图，配置读取靠pyyaml。没有深度学习框架，也没有 RL 平台，装完四个库就能跑。这一点对二次开发很友好——你不会被某个框架的版本锁死，把核心循环搬到别的语言里也不难。

runners和irl分开是个好习惯。runners只做三件事：解析参数、调用求解器、把结果写进outputs/；irl里全是纯函数式的数学，不碰文件、不碰绘图。以后要把训练塞进调度系统，只替换runners那一层就行。

2.2 环境接口：一个 5×5 网格怎么说清自己的动力学

环境层要对外提供的东西非常少：状态数、动作数、转移概率、特征函数、折扣因子。这四个量一旦确定，后面的求解器就可以完全通用。下面是环境类里最关键的部分。

# envs/gridworld.py import numpy as np class GridWorld: """确定性 4 动作网格：0=上 1=下 2=左 3=右""" def __init__(self, size=5, walls=(), goal=(4, 4), gamma=0.9): self.size = size self.n_states = size * size self.n_actions = 4 self.walls = {self.idx(*w) for w in walls} # 不可达格集合 self.goal = self.idx(*goal) self.gamma = gamma self.transition = self._build_transition() # 形状 (S, A, S) def idx(self, r, c): """二维坐标压成一维状态编号，热力图还原时再 divmod 拆回来""" return r * self.size + c def _build_transition(self): T = np.zeros((self.n_states, self.n_actions, self.n_states)) moves = [(-1, 0), (1, 0), (0, -1), (0, 1)] for s in range(self.n_states): r, c = divmod(s, self.size) for a, (dr, dc) in enumerate(moves): nr, nc = r + dr, c + dc # 越界或撞墙时停在原地，保证每行概率和为 1 if not (0 <= nr < self.size and 0 <= nc < self.size) \ or self.idx(nr, nc) in self.walls: T[s, a, s] = 1.0 else: T[s, a, self.idx(nr, nc)] = 1.0 return T def feature(self, s): """单状态 one-hot，权重向量就能直接画成一张 S 格子的热力图""" return np.eye(self.n_states)[s]

transition的维度顺序是(S, A, S)，也就是T[s, a, s']：给定状态和动作，下一个状态的概率分布。这个顺序在后面写值迭代时能省掉一堆转置，代价是记住一点——对s'求和永远是最右边那个轴。

撞墙处理成自环而不是「丢弃该动作重采样」，是这段代码里容易被忽略但很关键的一处。如果写成丢弃，转移矩阵每行的和就不再是 1，值迭代和状态访问频率统计都会悄悄偏掉，最后表现成「训练看起来在跑，但学出来的奖励永远不对劲」。tests/test_transition.py里那条np.allclose(T.sum(axis=-1), 1.0)的断言就是为了卡住这类问题。

feature用的是状态 one-hot，好处是简单直白，权重直接就是「每个状态值多少分」，画出来一目了然。真实项目里这一层要换掉，换成距离、拥堵度、能耗这类业务特征，第 6 章会细说。但要提醒一句：特征一旦不是 one-hot，权重就不能直接当热力图看了，得先算feature_matrix @ w才能解释。

2.3 求解器接口：为什么两个算法能共用同一个 fit

irl/base.py只定义了三个方法，两个算法各实现一遍。

# irl/base.py class IRLSolver: def fit(self, trajectories): """输入示范轨迹列表，每条轨迹是状态编号序列，返回奖励权重 w""" raise NotImplementedError def reward(self): """返回形状为 (S,) 的奖励向量，供画图和下游策略求解使用""" raise NotImplementedError def feature_expectation(self, trajectories): # 专家特征期望：所有状态特征的平均，轨迹长短不一时按步数加权 feats = [self.env.feature(s) for traj in trajectories for s in traj] return np.mean(feats, axis=0)

feature_expectation是整个 IRL 的核心统计量，很多论文里写作 μ_E。它的含义很直白：专家在状态空间里「停留分布的加权平均」。这里用所有轨迹所有步的特征做简单平均，等价于按步数加权，轨迹越长贡献越大。如果你的任务里每个示范的长度差异极大，就得改成先对每条轨迹求平均、再对轨迹求平均，否则一条超长轨迹会淹没其他所有示范。

fit的签名里没有超参，超参放在构造函数里。这样切算法只改一个类名：

solver = MaxEntIRL(env, lr=0.05, n_iters=60, l2=1e-3) w, history = solver.fit(expert_trajectories)

换成ApprenticeIRL(env, ...)时，runners/run_maxent.py里除了 logger 的名字，其他一行都不用动。这种设计在小项目里看着有点「过度抽象」，但一旦要做对照实验，两种方法的差异就能被压到最小——影响结果的只有算法本身，不是脚本差异。

提示：两个求解器共享env对象而不是重新构造。环境里有转移矩阵和特征矩阵，5×5 时无所谓，换成 100×100 网格，重复构造会明显拖慢实验。

3. 跑通第一个实验：从生成专家轨迹到画出奖励热力图

实验流程是三步：先造一个「真实奖励」，用它求最优策略并采样一批示范轨迹；再把示范喂给求解器，让它反推权重；最后把反推的权重和真实的比一比。这个对比不是为了作弊，而是检验实现有没有写错——如果连自己造的数据都恢复不出来，别指望它能处理真实示范。

3.1 依赖安装与最小运行命令

建虚拟环境这一步别省，numpy和scipy在不同 Python 版本下的行为差异足够让你排查半天。

# Linux / macOS python -m venv .venv source .venv/bin/activate pip install numpy scipy matplotlib pyyaml tqdm # Windows 用下面这一行激活 # .venv\Scripts\activate python -m runners.run_maxent --config configs/grid_5x5.yaml

--config是唯一必须传的参数，其余都能用命令行覆盖，比如--n-iters 100 --lr 0.02，覆盖逻辑在runners/run_maxent.py里做，优先级是「命令行 > 配置文件 > 代码默认值」。这套约定和大多数训练脚本一致，不用额外学。

跑起来之后标准输出会打印几行关键信息：专家轨迹的平均步数、初始特征期望的范数、每一轮梯度的模长、最终权重的范数。这几项是判断实验是否健康的体检指标，比盯着损失曲线更直接。

3.2 输出物怎么读：热力图、权重收敛曲线与策略回放

跑完以后outputs/下会有三样东西。

第一是weights_maxent.npy，一维数组，长度等于状态数。用np.save存是为了下次直接加载，不用重训。

第二是reward_heatmap.png，把权重按size × size重新排布画成热力图。看这张图有个诀窍：不要只看最大值在哪，要看「负值区域连不连成一片」。如果负值散成孤立的点，多半是迭代次数不够；如果负值在目标周围围了一圈，那是典型的学习不充分——模型还没意识到接近目标是有价值的。

第三是history.json，每轮记录梯度模长、权重模长、专家期望与当前期望的 L2 距离。这三个量的正常形态是：梯度模长先快速下降，然后在某个量级上下抖动，如果它单调趋近于 0，很可能是lr太小或者提前收敛到了退化解（权重全 0）。

策略回放这一步容易被跳过，但对新手最有价值。加载权重后跑一次值迭代，把最优策略的箭头画在网格上，和专家的策略并排看。如果两边在「绕墙走」这件事上选路不同但代价接近，说明奖励学到了，只是策略不唯一——这是可以接受的；如果一边直奔目标另一边在兜圈子，那就是奖励没学到，回去看第 5 章。

3.3 参数调优：轨迹条数、学习率与正则项

配置文件里参数不多，但每一个都有讲究。

参数	默认值	作用	调大之后	调小之后
`n_trajectories`	30	专家示范条数	特征期望更准，但采样成本上升	方差大，学出的奖励抖动
`horizon`	40	每条轨迹最大步数	覆盖更多状态，长尾信息进来	状态覆盖不足，权重集中在一小片
`n_iters`	60	外层次数	更充分，但可能开始过拟合噪声	权重没长开
`lr`	0.05	权重更新步长	收敛快但来回横跳	太慢，几十轮看不出变化
`l2`	1e-3	权重 L2 正则	权重整体被压小，热力图变淡	权重可能发散
`gamma`	0.9	折扣因子	更看重长期回报	退化成只看一步

这几行配置在configs/grid_5x5.yaml里，改完直接重跑。我一般的调参顺序是：先把n_trajectories和horizon加起来让专家数据足够密，再定lr，最后用l2收一下权重的尺度。反过来调很容易陷入「怎么调都不收敛」的循环。

# runners/run_maxent.py 里读取配置的核心几行 cfg = { "n_trajectories": 30, # 示范轨迹条数，太少会让特征期望方差过大 "horizon": 40, # 单条轨迹最大步数 "n_iters": 60, # 奖励权重的更新轮数 "lr": 0.05, # 步长，5x5 网格上 0.05 是比较稳的起点 "l2": 1e-3, # 正则强度，防止权重随迭代无限膨胀 "gamma": 0.9, } cfg.update({k: v for k, v in vars(args).items() if v is not None}) # 命令行覆盖

这里的cfg.update那行用了字典推导过滤None，是因为argparse里没传的参数默认是None。如果直接整字典覆盖，配置文件里的值会被全部冲掉——这是个很常见的翻车点，写训练脚本时值得留意。

lr的量级跟特征尺度强相关。one-hot 特征下每个状态的特征范数是 1，0.05 合适；换成连续特征（比如把坐标归一化后当特征），特征范数可能小到 0.01，这时候lr得放大到 1 以上才有明显更新。换句话说，换特征就得重调lr，别想着抄配置。

4. 核心算法逐行读：最大熵 IRL 的软值迭代与梯度更新

到这一层，代码量其实不大，maxent.py主体加起来一百多行。但里面每一行都对应一个数学定义，读的时候不能跳。最大熵 IRL 的思路可以概括成一句话：在所有能解释专家行为的奖励函数里，选熵最大的那个。落到实现上，就是「用当前奖励算出软最优策略，再比较这个策略的访问频率和专家的访问频率，差多少就把奖励往哪个方向推」。

4.1 前向：软值迭代把奖励摊成状态访问频率

普通值迭代取max，软值迭代取log-sum-exp。这个替换是整个算法的心脏。

# irl/maxent.py import numpy as np from scipy.special import logsumexp def soft_value_iteration(reward, transition, gamma, n_iters=200, tol=1e-6): """reward: (S,) 奖励向量；transition: (S, A, S) 返回软值函数 V: (S,)""" V = np.zeros(reward.shape[0]) for _ in range(n_iters): # Q(s,a) = R(s) + gamma * sum_s' T(s,a,s') * V(s') Q = reward[:, None] + gamma * (transition @ V) # (S, A) V_new = logsumexp(Q, axis=1) # 对动作做 soft-max if np.max(np.abs(V_new - V)) < tol: return V_new V = V_new return V def soft_policy(reward, transition, gamma, V=None): """策略 pi(a|s) = exp(Q(s,a) - V(s))，每个状态上动作概率和为 1""" V = soft_value_iteration(reward, transition, gamma) if V is None else V Q = reward[:, None] + gamma * (transition @ V) return np.exp(Q - V[:, None])

transition @ V这一步是(S, A, S) @ (S,)，结果是(S, A)，正好是每个状态动作对的「下一步期望值」。这个写法比手写einsum('sat,t->sa', T, V)快得多，因为走的是 BLAS 矩阵乘法。

logsumexp而不是np.log(np.sum(np.exp(Q), axis=1))，是防溢出的。网格环境下 Q 值通常只有个位数，看不出差别；一旦奖励尺度上去（比如把真实代价放进来），exp直接炸成inf，log出来就是nan，然后整个权重向量一夜之间变成nan数组。这种问题排查起来极其耗时，因为报错的地方和出错的地方隔了好几层。

有了策略，下一步算状态访问频率。

def state_visitation(policy, transition, start_dist, horizon): """policy: (S, A)；start_dist: (S,) 初始状态分布；返回平均每步访问频率""" # 策略诱导的转移矩阵 P_pi(s, s') = sum_a pi(s,a) * T(s,a,s') P_pi = np.einsum("sa,sat->st", policy, transition) mu = np.zeros(policy.shape[0]) d = start_dist.astype(float).copy() for _ in range(horizon): mu += d # 累加每一步的访问分布 d = d @ P_pi # 往前推一步 return mu / horizon # 归一化成“平均每步”

P_pi那行用einsum在策略和转移矩阵上做加权，循环外算一次就够，别放进 for 里。horizon决定往前看多远，不要超过专家轨迹的实际平均长度太多——超出部分d已经衰减到接近 0，纯属浪费算力。

4.2 反传：梯度就是两条访问频率曲线之间的差

整个训练循环只有八行，但每一行都值得盯住。

def fit(self, expert_trajectories): mu_expert = self.feature_expectation(expert_trajectories) # 专家的状态访问频率 w = np.zeros(self.n_features) history = [] for it in range(self.n_iters): reward = self.feature_matrix @ w # R(s) = w^T phi(s) V = soft_value_iteration(reward, self.transition, self.gamma) pi = soft_policy(reward, self.transition, self.gamma, V) mu = state_visitation(pi, self.transition, self.start_dist, self.horizon) grad = mu_expert - mu - self.l2 * w # 关键一行 w += self.lr * grad history.append(float(np.linalg.norm(grad))) return w, history

grad = mu_expert - mu - self.l2 * w：专家经常去的状态，mu_expert大，梯度为正，权重被推高；模型当前策略经常去而专家不去的状态，mu大，梯度为负，权重被压低。就这么简单。

注意这里是「加」不是「减」。最大熵 IRL 的似然函数对权重的梯度是mu_expert - mu（专家期望减模型期望），做梯度上升所以是w += lr * grad。很多复现踩坑就踩在这——写成了w -= lr * grad，结果权重符号整个反过来，热力图看着像模像样但完全反了。

self.feature_matrix是(S, d)的特征矩阵，第 s 行是状态 s 的特征。feature_matrix @ w一次算出所有状态的奖励，比循环调env.feature(s)快一个数量级。feature_expectation里那条路径是给示范数据用的，两条路径算出来的必须是同一个phi，否则梯度方向就错了——这也是第 5 章第一条要说的。

4.3 学徒学习版本：用投影和 margin 代替梯度

apprentice.py走的是另一条路。它不优化似然，而是维护一个「已经能解释的历史策略集合」，每一轮找当前权重下最违反约束的那个策略，然后把权重往缩小 violation 的方向推。

# irl/apprentice.py def apprentice_step(w, mu_expert, mu_history, lr=0.1, eps=1e-8): """mu_history: 历史策略的特征期望列表，每条形状 (d,)""" diffs = np.stack([mu_expert - m for m in mu_history]) # (K, d) margins = diffs @ w # 每个历史策略的 margin i = int(np.argmin(margins)) # 最违反约束的那一个 w_new = w + lr * diffs[i] norm = np.linalg.norm(w_new) return w_new / norm if norm > eps else w_new # 投影到单位球

margins的含义是：当前权重下，专家比第 i 个历史策略好多少。argmin找最小的那个，也就是「专家优势最不明显」的对手，往它身上推一把。最后那步归一化对应论文里的||w|| <= 1约束，不归一化的话每轮都在放大权重，margin 数字会越跑越大但方向没变。

两个方法在同一份数据上的表现差别，大致是这样：

维度	最大熵 IRL	学徒学习
输出	稠密奖励，每个状态一个值	稠密奖励，但只有方向可靠
收敛速度	慢，通常要几十轮	快，十几轮就有形状
对示范数量	敏感，少了方差大	相对宽松
权重可解释性	好，能直接看热力图	一般，尺度靠归一化定
典型用途	需要精确奖励做下游规划	只需要策略像专家

我自己的用法是：先用学徒学习快速确认「示范数据里有没有可学的东西」，有形状了再切最大熵去精修。反过来先跑最大熵，常常因为初始化附近梯度太小，看不出任何变化就放弃了。

4.4 数值安全与收敛判据

有几个数要在训练循环里盯着，不然出了问题是黑匣子。

一是np.linalg.norm(grad)。正常形态是前几轮快速下降，然后在一个小量级抖动。如果它一直在涨，lr太大；如果它从第一轮就接近 0，多半是mu_expert和mu几乎一样，也就是示范数据不够分散，策略随便走走就撞上了专家的分布。

二是np.linalg.norm(w)。one-hot 特征下权重范数会随迭代缓慢增长，这是正常的，因为奖励的绝对尺度不影响策略。但如果它指数级涨，说明正则没起作用，把l2从 1e-3 提到 1e-2 看看。

三是软值迭代的内部收敛。soft_value_iteration的tol默认 1e-6，n_iters200。如果它每次都跑满 200 轮还没收敛，说明gamma * max|V|这一项在放大，通常出现在gamma设成 0.99 以上、奖励又比较大的时候。把gamma降到 0.9 附近，或者把奖励先做个归一化，都能缓解。

5. 避坑与排查：五个让实验原地打转的常见问题

这一节的内容基本都来自我自己踩过的坑，按「现象 → 原因 → 解决」整理，遇到对得上的直接照着查。

5.1 训练相关：权重不收敛、学不出形状

现象一：梯度模长在几轮内来回横跳，热力图的红蓝区域每轮翻转一次。

原因通常是lr相对特征尺度偏大。one-hot 特征下梯度每个分量最大是 1，lr取到 0.5 就足以让单个状态的权重在两轮之间从 -1 冲到 +1，看上去像在震荡。解决是把lr降到 0.05 甚至 0.02，同时把n_iters加上去，用更多小步代替少量大步。顺带检查一下特征有没有归一化——如果你的特征某一维量级是 100，那一维上的有效学习率就是其他维的 100 倍。

现象二：权重全为 0，或者所有状态一个值。

这是最典型的「专家期望和模型期望没差别」。先打印mu_expert：如果它几乎是均匀分布，说明示范轨迹太短或太少，专家还没走出有信息量的路径就结束了。把horizon从 20 提到 40、n_trajectories从 10 提到 30，立刻能看到分布出现峰。另一种可能是start_dist设错了——如果所有示范都从同一个角落出发，而你训练时的start_dist是均匀分布，两边访问频率从一开始就对不上，梯度会把你往错误方向推。

现象三：损失（梯度模长）单调下降得很好看，但学出的策略明显不合理。

常见原因是奖励只依赖状态，动作偏好表达不出来。网格环境里这通常不是问题，因为动作的效果体现在状态转移上；但一旦环境里有「同一状态不同动作导致不同后果」（比如带惯性的移动），纯状态奖励就是不够的。改法是把特征从phi(s)换成phi(s, a)，feature_matrix的维度从(S, d)变成(S*A, d)，奖励也相应变成(S, A)。改动量不大，但要把soft_value_iteration里reward[:, None]那行换成直接传(S, A)的奖励矩阵。

5.2 数值与随机性：NaN、行和不为 1、结果不可复现

现象四：跑着跑着权重变成nan，前面几轮完全正常。

八成是exp溢出。检查两处：soft_policy里的np.exp(Q - V[:, None])——减V这一步已经是在做归一化，理论上不会溢出，但如果V本身是nan，减法救不回来；以及soft_value_iteration里有没有偷懒用np.log(np.sum(np.exp(Q), axis=1))。把这两处统一改成scipy.special.logsumexp，并且在整个循环开头加一句assert np.isfinite(reward).all()，能在第一轮就抓到问题，而不是等到二十轮后。

现象五：同样的代码，两次跑出来的热力图不一样。

三个随机源没固定：专家采样、初始状态分布、以及某些实现里值迭代的初始化。解决是在runners入口处统一设种子：

import numpy as np, random def set_seed(seed=0): random.seed(seed) np.random.seed(seed) # 单进程够用；多进程还要额外处理

设完种子后再跑两次，weights_maxent.npy应该逐位相同。如果还是不同，检查有没有把结果写进同一个文件导致互相覆盖——outputs/下带时间戳命名就没有这个问题。

还有一条不在上面五条里但很常见：转移矩阵行和不为 1。表现是训练能跑、指标也在动，但学出来的奖励总是偏，而且换随机种子结论也不变。tests/test_transition.py里那条断言就是为了这个，pytest tests/跑一遍只要一秒，比事后排查省事得多。

6. 二次开发：把奖励推断接到 Java 业务侧的两种改法

这套代码是脚本形态，做实验很舒服，但很多团队的业务系统主栈是 Java——这时候硬要用脚本跑线上推理是不现实的。我实践下来有两条路，选哪条取决于「奖励多久更新一次」。

6.1 方案 A：Python 离线训练，Java 只加载权重

绝大多数场景该走这条。奖励函数的变化频率远低于请求频率，一天甚至一周重训一次都够。Python 侧把权重导出成一个朴素的文本文件，一行一个 double，Java 侧读进来做个线性组合就行。

// RewardModel.java —— 只负责加载权重和查询奖励 import java.io.IOException; import java.nio.file.*; import java.util.List; public final class RewardModel { private final double[] weights; // 与特征维度一一对应 private RewardModel(double[] weights) { this.weights = weights; } public static RewardModel load(Path path) throws IOException { List<String> lines = Files.readAllLines(path); double[] w = new double[lines.size()]; for (int i = 0; i < lines.size(); i++) { w[i] = Double.parseDouble(lines.get(i).trim()); // 空行会抛异常，提前清掉 } return new RewardModel(w); } /** 状态 one-hot 场景下直接取对应维；换成业务特征时改成点积 */ public double reward(int stateIdx) { return weights[stateIdx]; } /** 业务特征版本：feature 由调用方按同一套规则算好再传进来 */ public double reward(double[] feature) { double sum = 0.0; for (int i = 0; i < weights.length; i++) sum += weights[i] * feature[i]; return sum; } }

加载和查询都很轻，weights是不可变数组，多线程读没有同步问题。要注意的是特征顺序必须两边一致——Python 侧feature_matrix的列顺序和 Java 侧feature的构造顺序一旦错位，奖励会算得似是而非，而且不会报错。

6.2 方案 B：Java 侧重写推理，Python 只留训练

如果延迟要求高到不能跨进程，就把soft_value_iteration和state_visitation搬到 Java。这两段没有依赖特殊的库，logsumexp自己写十行就够，转移矩阵用double[S][A][S]或者展平成一维数组都能跑。真正的成本不在算法，在于两边要对齐。

我的验证方法是：用同一批专家轨迹，Python 侧和 Java 侧各算一次mu_expert，把两个数组打出来做差，要求最大绝对值小于 1e-6。这一步过不了，后面的一切都不用谈。之后再对同一组权重跑一次软值迭代，比较V向量的差，通常能压到 1e-9 量级。做完这两步对齐，剩下的就是普通的服务封装。