WHAM目的是把有偏采样的统计结果转换为无偏采样的统计结果,在伞型采样或副本交换等模拟中广泛应用于计算potentials of mean force (PMF)。其大致原理思路如下:
假设试图计算函数 $ f(x) $ 的期望, 其中 $ x $ ~ $ p(x) $, 对 $ E[f(x)] $ 有如下估计:
\begin{equation} E[f(x)]=\int f(x) p(x) d x \approx \frac{1}{n} \sum_i f\left(x_i\right) \end{equation}
蒙特卡罗抽样法就是简单地从分布 $ p(x) $ 中抽取 $ x $ 样本, 然后取所有样本的平均值, 得到期望值的估计值。那么问题来了, 如果 $ p(x) $ 很难取样呢? 我们是否能根据已知的, 容易取样的分布来估计期望值呢?
我们做如下变换:
\begin{equation} E[f(x)]=\int f(x) p(x) d x=\int f(x) \frac{p(x)}{q(x)} q(x) d x \approx \frac{1}{n} \sum_i f\left(x_i\right) \frac{p\left(x_i\right)}{q\left(x_i\right)} \end{equation}
通过这种方法,估计期望值就能从另一个分布 $ q(x) $ 中采样, 而 $ p(x)/q(x) $ 称为采样比或采样权重, 它作为校正权重来抵消从不同分布中采样的概率。在实际应用中可以从有偏采样中推导无偏置的自由能分布。
模拟后我们首先得到 $ p_{ij}^b $ 为第 $ i $ 次模拟中第 $ j $ 个bin中(有偏差的)概率估计值。我们假设 $ p_{ij}^b $ 与第 $ j $ 个bin的无偏估计 $ p_{j}^u $ 的关系为:
\begin{equation} p_{ij}^b=f_i c_{ij} p_{j}^u \end{equation}
其中 $ c_{ij} $ 是biasing factor, 例如对于温度偏倚, $ c_{ij}=\exp \left[-\left(\beta_i-\beta_0\right) E_j\right] $, 对于结构偏倚如US, $ c_{ij}=\exp \left[-\beta V_i\left(x_j\right)\right] $。 $ f_i $ 是normalizing constant,使得 $ \sum_j p_{ij}^b=1 $ , 即为第 $ i $ 次模拟的配分函数。
使用最大似然法推导,第 $ i $ 次模拟,观察到直方图 $ n_{ij} $ ($ j = 1… M $)的概率的多项式分布如下:
\begin{equation} P\left(n_{i 1}, \ldots, n_{i M}\right)=\frac{N_{i} !}{\prod_{k=1}^M n_{i k} !} \prod_{j=1}^M\left(p_{ij}^b\right)^{n_{i j}} \end{equation}
其中 $ M $ 是bin数目, $ N_{i} $ 是第 $ i $ 次模拟产生的总样本数。所有模拟的概率分布如下:
\begin{equation} P\left(n_{i 1}, \ldots, n_{i M}, \ldots, n_{S 1}, \ldots, n_{S M}\right)=\prod_{i=1}^S\left[\frac{N_{i} !}{\prod_{k=1}^M n_{i k} !} \prod_{j=1}^M\left(p_{ij}^b\right)^{n_{ij}}\right] \end{equation}
其log-likelihood function如下:
\begin{equation} L=\sum_{i=1}^S \sum_{j=1}^M n_{i j} \ln \left(f_i c_{i j} p_j^{u}\right) + \text{const.} \end{equation}
对这个式子进行优化, 其中要满足 $ \sum_{j=1}^M f_i c_{i j} p_j^{u}=1 $, 引入Lagrange算子 $ \lambda $, 如下有:
\begin{equation} F=\sum_{i=1}^S \sum_{j=1}^M\left[n_{i j} \ln \left(f_i c_{i j} p_j^{u}\right)+\lambda_i f_i c_{i j} p_j^{u}\right] \end{equation}
对 $ f_i $ 求导得:
\begin{equation} \begin{aligned} \frac{\partial F}{\partial f_k} & =\sum_{j=1}^M\left[\frac{n_{k j}}{f_k}+\lambda_k c_{k j} p_j^{u}\right] \\ & =\frac{N_k}{f_k}+\lambda_k \sum_{j=1}^M c_{k j} p_j^{u}\end{aligned} \end{equation}
对 $ p_j^{u} $ 求导得:
\begin{equation} \begin{aligned} \frac{\partial F}{\partial p_k^{u}} & =\sum_{i=1}^S\left[\frac{n_{k j}}{p_k^{u}}+\lambda_i f_i c_{i k}\right] \\ & =\frac{\sum_{i=1}^S n_{i k}}{p_k^{u}}+\sum_{i=1}^S \lambda_i f_i c_{i k}\end{aligned} \end{equation}
将 $ \frac{\partial F}{\partial f_k} = 0 $ 有:
\begin{equation} \lambda_k=\frac{-N_k}{\sum_{j=1}^M f_k c_{k j} p_j^{u}}=-N_k \end{equation}
将上式代入 $ \frac{\partial F}{\partial p_k^{u}} = 0 $ 有:
\begin{equation} p_j^{u}=\frac{\sum_{i=1}^S n_{i j}}{\sum_{i=1}^S N_i f_i c_{i j}} \end{equation}
将此式子与 $ f_i^{-1} = \sum_{j=1}^M c_{i j} p_j^{u} $ ($ f_i $ 使得 $ \sum_j p_{ij}^b=1 $) 进行联合优化即为WHAM核心思想,。一般设置 $ f_i $ 为 $ 1 $ , 求得无偏估计 $ p_j^{u} $, 以此更新 $ f_i $ , 然后重复这个过程直到收敛。
本文由 KoN 创作,如果您觉得本文不错,请随意赞赏
采用 知识共享署名4.0 国际许可协议进行许可
本站文章除注明转载/出处外,均为本站原创或翻译,转载前请务必署名
原文链接:/archives/wham---weightedhistogramanalysismethod
最后更新:2023-11-27 18:43:36
Update your browser to view this website correctly. Update my browser now