通过包络分布采样实现更平滑的炼金术变换用于自由能估计

《Journal of Chemical Theory and Computation》：Smoother Alchemical Transformations via Enveloping Distribution Sampling for Free-Energy Estimation

　　

相对自由能（例如，用于溶剂化或蛋白质-配体结合）的计算精度取决于两个炼金术终态之间相空间变换的平滑性。平滑变换确保了在平衡（EQ）模拟中连接两个终态的相邻中间态之间有足够的相空间重叠，并在非平衡（NEQ）模拟中产生更少的耗散功。传统的能量插值（EI）耦合方案通过线性组合终态势能来构建中间态。研究人员表明，包络分布采样（EDS）耦合方案（EI的推广形式，其中相应的玻尔兹曼因子被线性组合）是一种更灵活的选择。通过使用负平滑参数，EDS方案沿变换轴增加了采样相空间的局部曲率，从而避免了相变并创造了更平滑的变换。研究人员在从谐波振荡器和Ising模型系统到FreeSolv数据集上的绝对水化自由能（AHFE）计算等日益复杂的设置中验证了这种行为。与传统的EI方案相比，EDS在模型系统计算中始终产生更准确和统计上更稳健的自由能估计，而在NEQ机制下，AHFE计算中观察到明显优势，其中耗散较小的转变导致更可靠的自有能估计。

自由能是热力学系统的基本平衡性质，在生物学、材料科学和药物发现中具有高度相关性。然而，实际应用中两个系统的相空间重叠往往不足，因此需要通过耦合参数λ∈[0,1]沿炼金术坐标构建非物理中间态，这一过程称为炼金术变换。传统能量插值（EI）耦合方案通过线性组合终态势能（U_EI(r,λ)=(1-λ)U_A(r)+λU_B(r)）构建中间态，但常面临相变（如有序/无序转变）和端点灾难等问题，需借助软核势能来缓解。相比之下，包络分布采样（EDS）方案通过线性组合玻尔兹曼因子（U_EDS(r,λ;s,E(λ))=-(sβ)^-1ln[(1-λ)e^-sβU_A(r)+λe^{-sβ(U_B(r)-E(λ))}]）引入平滑参数s和能量偏移E(λ)，提供更灵活的控制。研究人员首次系统研究了负s参数在炼金术变换中的使用，发现其能增加局部曲率从而避免相变，并在谐波振荡器、Ising模型以及FreeSolv数据集（642个小分子）的绝对水化自由能（AHFE）计算中验证了EDS方案的优势。该研究发表在《Journal of Chemical Theory and Computation》。研究结论表明，EDS方案在非平衡（NEQ）模拟中通过配置感知的有效耦合参数λ_EDS显著减少耗散功，提高自由能估计精度，尤其在分子尺寸较大的溶质中效果更明显，为自由能计算提供了简单而有效的改进方法。

研究人员主要采用以下关键技术方法：EDS耦合方案（线性组合玻尔兹曼因子，引入平滑参数s和能量偏移E(λ)）；负s参数（通过公式s_α=αβ[|?ΔU?₀-?ΔU?₁|+2σ⁰_ΔU+2σ¹_ΔU]^-1选择，α=-10）；EQ和NEQ模拟（分别使用Metropolis-Hastings蒙特卡洛和分子动力学模拟，其中MD使用OpenMM 8.1软件、GAFF 1.8力场和TIP3P水模型）；样本队列来源为FreeSolv数据集（版本0.52，含642个小分子）。AHFE计算采用解耦方法，使用软核势能处理Lennard-Jones和静电相互作用，并采用反应场处理长程静电。

**Analytical Demonstration Using Harmonic Oscillators as Toy Systems**
研究人员通过解析计算分析了谐波振荡器（HO1、HO2、iHO）中不同s和E(λ)选择对局部曲率和ΔU分布的影响。结果表明：负s参数使得坐标分布变窄、ΔU分布向E(λ)偏置，从而平滑了相空间变换；在力常数不同的HO2系统中，EDS方案使坐标分布沿λ更均匀地传输；在反相谐波振荡器iHO中，EDS方案通过增加曲率避免了边界处积累并允许中间位置采样。

**Comparison of EQ and NEQ Methods on Toy Systems**
研究人员在五个玩具系统（HO1、HO2、iHO、2D-Ising1、2D-Ising2）中比较了EQ和NEQ方法下EI与EDS方案的表现。EQ结果显示：在HO1中两种方案收敛相似；在HO2和iHO中，EI方案因缺乏相空间重叠而失败，而EDS方案成功收敛（RMSE<1k_BT）。NEQ结果显示：在HO2和iHO中，EDS方案比EI方案更快收敛、耗散更小，且λ_EDS表现出非单调适应行为以分配更多资源采样过渡区。在2D-Ising系统中，EDS方案避免了EI方案因相变导致的磁化状态采样不足，并提供了更稳定的自由能估计。

**Comparison on Absolute Hydration Free Energies**
研究人员在FreeSolv数据集的642个小分子上计算AHFE，比较EQ和NEQ协议下EI与EDS方案。EQ条件下，两种方案随模拟时间收敛相似，RMSE在1-2 ns/窗口时均低于0.5 kJ/mol。NEQ条件下，EDS方案在所有变换速度下均优于EI，RMSE更低且化学精度（绝对误差<4.184 kJ/mol）分子比例更高；EDS方案在98.5%的变换中产生更少的耗散功，尤其在分子量较大时效果显著。通过分析最重分子（decachlorobiphenyl）的NEQ过程，研究人员发现EDS方案通过配置感知的λ_EDS在反向切换（再耦合）中提前排斥溶剂，使空腔形成更平滑，从而减少耗散功。进一步研究表明，耗散功与分子量及2D分子描述符（如chi0v、labuteASA）高度相关（Pearson系数ρ=0.93-0.99），而EDS方案的改进主要源于空腔形成过程的调节。

在讨论部分，研究人员指出NEQ条件下EDS方案的优势源于其对迟滞根本原因（构型-耦合不匹配）的针对性解决，通过自适应有效耦合λ_EDS独立于具体的软核势能形式，并能稳健处理不同的软核参数选择。结论部分翻译如下：EDS方法最初作为耦合方案引入，用于在单个模拟中采样多个系综。本研究首次系统研究和验证了负平滑参数s在成对炼金术变换中的使用。研究人员证明，这一特定应用将EDS方案转化为自由能计算中强大而灵活的工具，提供了对两个终态之间中间态构建的优越控制。研究结果一致表明，通过采用负s参数，EDS耦合方案有效平滑了相空间变换，即防止了变换过程中相对终态能量（ΔU(r)分布）的突变。这增加了势能景观的局部曲率，从而避免了传统EI耦合方案模拟中出现的相变和端点灾难。该方法的优势在EI耦合引入相空间重叠突变的系统中最为明显。在谐波和Ising测试系统中，EDS始终提高了EQ和NEQ自由能估计的稳定性和准确性。在AHFE计算中，EQ采样下EDS相对于EI没有显著优势，两种方案收敛到参考值相似；但在NEQ机制下，EDS表现出明显优势，这归因于配置感知的λ_EDS，它根据瞬时ΔU(r)与能量偏移E(λ)之间的不匹配来调节溶质插入速率，从而减少耗散并改善自由能估计。特别是，研究人员发现AHFE估计的“难度”与溶质分子大小几乎完全相关，而溶质插入速率被λ_EDS调节。总之，EDS耦合方案为自由能计算提供了简单、基本但显著的优势。在不增加计算开销的情况下，它带来了更平滑、更稳定的炼金术变换。研究人员预计这项研究将为其他高级自由能方法创造未来机会。