结合自由能理论计算

我第一次做 结合自由能理论计算是评一个系列 8 个抑制剂，用 MM/PBSA 排序后挑了 top2 去湿实验，结果命中的是排第 5 的那个。复盘才发现我只算了单帧结构的平均，没做构象系综取样，柔性配体的熵贡献全丢了。那次让我把 结合自由能理论计算 从”算一个数”改成了”取一个分布”。

一、为什么 结合自由能理论计算先做取样

结合自由能 ΔG_bind = G_complex − G_receptor − G_ligand。终点法把每一项拆成 MM 能量 + 溶剂化（PB/GB）+ 熵（振动/准谐）。最大误差源不是公式，而是”用哪一帧算”。一个柔性配体在口袋里有十几种构象，单帧算出来的 ΔG 能差 3 kcal/mol。我做结合自由能理论计算时，默认先从一段平衡 MD 里抽 100~200 帧做系综平均，而不是挑一个”最好看”的构象。

二、核心原理：MM/PBSA 与 FEP 的取舍

MM/PBSA 用隐式溶剂近似连续介质，便宜但 PB 求解慢、GB 有参数依赖；FEP 走炼金路径，精度更高但贵。我一般先用 MM/PBSA 做排序筛，再用 FEP 校 top 候选。关键认知：MM/PBSA 的熵项（准谐分析）最容易翻车，柔性大的体系省掉它误差直接上 2~3 kcal/mol，但加上又极耗算力。我默认对刚性配体省熵、对柔性配体保留，按分子大小分级处理。

三、关键技术要点：溶剂与参数

结合自由能理论计算对隐式溶剂模型极敏感。GB 我用 GBneck2（GB-OBC2 偏早版本），PB 用非线性泊松-玻尔兹曼但要注意介电界面定义。力场上蛋白用 Amber ff19SB、配体用 GAFF2，二者要在同一 MD 轨迹里一致。曾经我用两套不同水模型跑出来的轨迹拼接算 PBSA，界面能跳变 1 kcal/mol。我现在的习惯是整套流程固定力场+水模型，绝不混搭。

四、实操流程：从轨迹到 ΔG

我的标准动作：平衡 MD（蛋白-配体 50~100 ns）→ 抽帧 → 每帧做 MM/PBSA 或跑 FEP → 系综平均 + 误差棒 → 排序并和活性数据比对。最容易被砍的是系综取样，但 结合自由能理论计算 的结论稳健性全靠它，我默认保留。

补一个实算：某激酶抑制剂系列，单帧 MM/PBSA 排序 top1 的结合能 −28 kcal/mol，系综平均后掉到 −24，而排第 4 的从 −22 升到 −25，反超。湿实验命中的恰恰是系综平均后的 top。如果我只交单帧结果，客户会浪费一轮合成在错误分子上。这之后我做结合自由能，取样帧数先按配体柔性定，不再图快取一帧。

五、常见踩坑：构象与质子化

结合自由能理论计算常见错是配体在生理 pH 下的质子化态错了，整个结合模式变。另一个坑是忽略结合诱导的蛋白构象变化，刚性受体假设让熵估错。还有人用 MM/PBSA 时不加熵直接报数，柔性配体误差被放大。金属中心配位若没特殊力场，能量也会乱飞。

六、复盘：结合自由能是排序工具

回过头看，结合自由能理论计算最值钱的是在湿实验前把候选浓缩到可验证的小集合，而不是替实验定 Kd。我接项目默认先问要排序还是要绝对定量，再定用 MM/PBSA 还是 FEP、取多少帧。把取样、溶剂、力场三件事做对，ΔG 排序才真能指导合成。被证明有用的，是愿意为多抽一百帧多等半天的清醒。下次有人只要”结合能排个序”，我会先问——取的是单帧还是系综。