分子模拟/伞形采样

伞形采样是一种在计算物理和化学中使用的采样方法。这种采样可以采样正常分子动力学采样忽略的稀有状态。因此，当系统发生系统性变化时，伞形采样可以改善自由能计算。

有偏分子动力学模拟

普通 MD 模拟在平衡状态下对系统进行采样。在正丁烷 (aq) 的 C-C-C-C 二面角时间序列的 MD 模拟中，只对 gauche 状态和 trans 状态进行采样。因为这种模拟只进行了 2 纳秒，所以自由能较高的状态（例如 cis 状态）不太可能发生。这些配置被忽略了，无法从这种模拟中计算这些状态的自由能。在这种情况下，需要一个人工偏差势来帮助分子越过能量势垒。使用偏差势，可以有效地采样稀有状态。

来自 2 纳秒分子动力学模拟的正丁烷 (aq) 的 C-C-C-C 二面角时间序列。

在这种情况下，需要一个谐波偏差势 $w(\phi )=-{\mathcal {V}}(\phi )$ 来抵消二面角势垒。

$w(\phi )=-k_{\phi }(1+\cos(n\phi -\delta ))$

来自使用偏差势来稳定 gauche 构象的 2 纳秒分子动力学模拟的正丁烷 (aq) 的 C-C-C-C 二面角时间序列。

有偏模拟捕获了高自由能状态。为了计算这些状态的自由能曲线，必须将有偏概率分布转换为无偏概率分布。

正丁烷在液态水中进行等温等压分子动力学模拟。对 C-C-C-C 二面角施加偏差势，以促进正丁烷旋转势垒的采样。

正丁烷在液态水中的等温等压分子动力学模拟

从有偏模拟中获取自由能曲线

势能 $V'(r)$ 包括反应坐标 $s$ 上的偏差势 $w(s)$ 是 ${\mathcal {V}}'(r)={\mathcal {V}}(r)+w(s)$

此势的概率分布为

$P'(s)={\frac {\int {\exp \left({\frac {-{\mathcal {V}}'(r)}{k_{B}T}}\right)\delta (f_{\alpha }(r_{1},...,r_{N})-s)d^{N}r}}{\int {\exp \left({\frac {-{\mathcal {V}}'(r)}{k_{B}T}}\right)}d^{N}r}}$

无偏势的概率分布为

$P(s)={\frac {\int {\exp \left({\frac {-{\mathcal {V}}(r)}{k_{B}T}}\right)\delta (f_{\alpha }(r_{1},...,r_{N})-s)d^{N}r}}{\int {\exp \left({\frac {-{\mathcal {V}}(r)}{k_{B}T}}\right)}d^{N}r}}$

根据这个公式，我们可以推导出，

$P'(s)=\exp \left({\frac {w(s)}{k_{B}T}}\right)P(s)\left\langle \exp \left({\frac {-w(s)}{k_{B}T}}\right)\right\rangle ^{-1}$

自由能曲线可以通过概率分布计算得出，公式为：

$A(s)=k_{B}T\ln P(s)$

使用此关系，偏置模拟的PMF可以转换为无偏PMF，公式为：

$A(s)=A'(s)-w(s)-k_{B}T\ln \left\langle \exp \left({\frac {-w(s)}{k_{B}T}}\right)\right\rangle$

$-k_{B}T\ln \left\langle \exp \left({\frac {-w(s)}{k_{B}T}}\right)\right\rangle$ 项记为 $F_{i}$。它通常是一个常数，在某些情况下不会影响相对能量，因此无需计算。可以通过 计算得出。

$\exp \left({\frac {-F_{i}}{k_{B}T}}\right)=\int {P(s)\exp \left({\frac {-w(s)}{k_{B}T}}\right)}ds=\int {\exp \left({\frac {-w(s)-A(s)}{k_{B}T}}\right)}ds$

文件：正丁烷伞形采样PMF.png

通过分子动力学模拟，并在模拟中引入偏置势，计算出的正丁烷C-C-C-C二面角的平均力势。偏置概率分布得出的PMF用红色表示。偏置势用蓝色表示。

这种方法可以得到所有可能状态的自由能曲线。在正丁烷（aq）的伞形采样中，所选偏置势涵盖了所有的反应坐标。一般情况下更加复杂，这会导致偏置势的确定更加复杂。

偏置势的选择

上一节讨论了正丁烷(aq)的偏置分子动力学模拟。反应坐标是一维且周期性的，偏置势能被选择为正丁烷二面角势能的负值^[2]。最佳偏置势能是自由能 $A(s)$ 的反向^[1]。然而，对于大多数情况，$A(s)$ 是未知的。对于一般情况，偏置势能需要沿反应坐标进行调整。因此，在参考点 $s_{0,i}$ 上施加的谐波偏置势能相对于反应坐标上的窗口 $i$ 被引入^[2]。

$w_{i}(s)={\frac {1}{2}}k_{i}(s-s_{0,i})^{2}$

因此，可以通过一系列在反应坐标上不同参考点进行的偏置MD模拟来获得完整的伞形抽样。

从伞形抽样数据计算平均力势

加权直方图分析法(WHAM)^[3]将一系列偏置分布函数转换为单个分布函数。$F_{i}$ 的值需要进行估计以给出 $A(s)$ 的正确值：$A(s)=A'(s)-w(s)+F_{i}$

真实分布P(s) 是每一步的加权平均值^[1]

$P(s)=\sum _{i=1}^{window}p_{i}(s)P_{i}(s)$

并且 $p_{i}=N_{i}\exp \left({\frac {-w(s)+F_{i}}{k_{B}T}}\right)$，其中 $N_{i}$ 是窗口 $i$ 采样的总步数^[3]。

结合 $\exp \left({\frac {-F_{i}}{k_{B}T}}\right)=\int {P(s)\exp \left({\frac {-w(s)}{k_{B}T}}\right)}ds=\int {\exp \left({\frac {-w(s)-A(s)}{k_{B}T}}\right)}ds$，$P(s)$ 和 $F_{i}$ 都可以得到。

分析伞形抽样的另一种方法是伞形积分，参见。

另请参阅

维基百科：伞形抽样

有关伞形抽样的更多信息，请参见