面向高通量的绝对结合自由能的计算

tink

随着计算机硬件和算法发展，绝对结合自由能计算不再昂贵。今天给大家介绍一种可以应用于高通量筛选的绝对结合自由能计算。

常见的绝对结合自由能计算一般是由Benoit Roux提出的DDM，该方法可以分为两大步
i）配体在溶液中的消失/生长
ii）配体在结合位点的消失/生长
其中，第二步中，配体在结合位点的消失生成，一般在分为配体的生长以及约束的消失。通常情况下，配体的生长这一步的炼金术计算收敛较为困难，需要极长的模拟时间。
为了解决这一困难，我们提出了一种新的热力学循环：LDDM,一种 “0”受力的热力学循环路线。

如图A实线所示，这总方法将分子的生长和约束的消失在同一步实现，极大的缩小了模拟的所需要的步骤。该方法的核心在于，
在如图B的关键自由度上，获取每一个cv上的力常数，随着alchemy的变化，力常数逐渐的减小，我们对分子的限制势逐渐增加。
这样就保证了，每一个自由度上配体所受到的力不变，这样就最大限度上的减少了在模拟中体系的扰动。
同时，由于模拟的步数缩减，加上使用HMR和DWS算法，模拟的时间也得到的极大的缩减。

模拟所需要的时间，仅仅是其他算法的一半的时间。

此外，在真实的计算中，只需要三天就可以计算得到一个蛋白配体的绝对结合自由能。（如果体系较小，可以更快）

在收敛速度的测试上，也是全面优于传统的DMM方法。
那么在准确度上，我们所有的测试体系都表现的很好，


在众多的数据测试中，R达到了0.96，AUE达到了0.7。
在精度的统计中，

绝对误差：计算值与实验值的差别，可能受很多因素（力场、算法系统误差，模拟收敛误差等）影响，对应Accuracy
粘滞误差：FEP正向和反向模拟的差别，反映正向和反向模拟是否收敛到了相同的构象空间，对应Accuracy（算法部分）
统计误差：BAR estimator的统计误差，反应已经得到采样的构象空间的采样充分程度，对应Precision
最后，也挑战了一下，fep+中一些较为困难的算法，

我们在模拟时间仅有一半的情况下，每一个体系都优于fep+。

最后介绍如何实现，该算法已经半自动化实现了，首先用户需要下载BFEE3（测试版）（github fhh2626）

输入文件准备
在使用 BFEE3 实现的 LDDM 方法时，用户干预的主要部分是提供分子系统的结构和拓扑文件，并向工具指明蛋白质和配体。具体来说，用户需要准备用于平衡的拓扑、坐标和力场文件。这些文件可以使用 CHARMM-GUI 或 AmberTools 自动生成。如果使用 Amber 力场，仅需提供拓扑和坐标文件（即 parm7 和 pdb 文件）（见 Fig. A）。对于 CHARMM 力场，还必须上传力场文件（即 prm 文件）（见 Fig. B）。随后，用户需要使用 MDAnalysis 语法指明蛋白质和配体（见 Fig. A/B 中的“Select protein”和“Select ligand”选项）。
完成上述步骤后，可以通过点击“Quick settings->Protein:Ligand / LDDM”菜单，然后点击“Generate Inputs”按钮生成所有自由能计算所需的输入文件（见 Fig. C）。
对于高级用户，可以通过提供修改过氢原子质量的 parm7/psf 文件，并在点击“Quick settings->Protein:Ligand / LDDM”后、点击“Generate Inputs”前，在“Advanced Settings”中将时间步长设置为 4 fs 来应用氢质量重新分配（hydrogen mass repartitioning）。如果需要，还可以调整窗口数量（number of windows）。

模拟与后处理
用户需要按以下步骤执行模拟或依次运行脚本：

复合物的平衡 (Equilibration of the complex):
1.1. 提交由 000_eq/000.1_eq.conf 和 000_eq/000.2_eq_ligandOnly.conf 配置的模拟。
1.2. 执行 000_eq/000.3_updateCenters.py。
1.3. 提交由 000_eq/000.1_eq2_re-eq.conf 配置的模拟。

分子在结合状态下 (Molecule in Bound State):
2.1. 执行 001_MoleculeBound/000_updateForceConstant.py。
2.2. 提交由 001_MoleculeBound/001_fep_doubleWide.conf 配置的模拟。

分子在非结合状态下 (Molecule in Unbound State):
3.1. 提交由 003_MoleculeUnbound/003_fep_doubleWide.conf 配置的模拟。

后处理 (Post-treatments)
BFEE3 会自动执行所有后处理，通过读取 001_MoleculeBound 中的 colvars.in.tmp，以及 001_MoleculeBound/output 中的 fep_doubleWide.colvars.traj 和 fep_doubleWide.fepout，以及 003_MoleculeUnbound/output 中的 fep_doubleWide.fepout（见 Fig.D）。一般推荐使用 BAR 估算器（BAR estimator）。

jxt

有搞头哇，竟然超过了fep+

zhouoh

想问下能贴一下GitHub仓库的链接吗？我这里搜索没有找到

Peng14

zhouoh 发表于 2025-4-25 05:13
想问下能贴一下GitHub仓库的链接吗？我这里搜索没有找到

可能是这个？https://github.com/fhh2626/BFEE2/tree/3.0-dev
也可以直接：git clone -b 3.0-dev https://github.com/fhh2626/BFEE2.git

tink

zhouoh 发表于 2025-4-25 05:13
想问下能贴一下GitHub仓库的链接吗？我这里搜索没有找到

google 直接搜关键词 应该可以

tink

Peng14 发表于 2025-4-25 08:48
可能是这个？https://github.com/fhh2626/BFEE2/tree/3.0-dev
也可以直接：git clone -b 3.0-dev https: ...

感谢帮帖

lonewolf659

哇，免费的ABFE工具，👍👍👍

吹梦到西洲

请问一下是只支持AMBER系列和CHARMM系列的力场吗？体系除了蛋白配体以外还有一个用UFF描述的基底，可以用同样的方法做吗？

ErikZhang

能贴一下paper么？想学习学习

tink

ErikZhang 发表于 2025-5-30 11:45
能贴一下paper么？想学习学习

目前还没有published，还要再等一会

hdjjys

请问只有NAMD才支持LDDM吗？如果是Gromacs用户，似乎只有Geometrical route，和BFEE2一样。

tink

hdjjys 发表于 2025-6-2 15:58
请问只有NAMD才支持LDDM吗？如果是Gromacs用户，似乎只有Geometrical route，和BFEE2一样。

我们前期的测试都是在namd上进行的，目前是完成了算法部分，接下来的bfee3，我们在努力实在gmx的实现。

hdjjys

tink 发表于 2025-6-3 09:40
我们前期的测试都是在namd上进行的，目前是完成了算法部分，接下来的bfee3，我们在努力实在gmx的实现。

谢谢，期待

tink

ErikZhang 发表于 2025-5-30 11:45
能贴一下paper么？想学习学习

doi.org/10.1038/s43588-025-00821-w

ErikZhang

tink 发表于 2025-6-18 21:12
doi.org/10.1038/s43588-025-00821-w

感谢，我去学习学习