Determination of the potential of mean force via bidirectional non equilibrium techniques

The Jarzynski identity is seemingly a better route than the CT to evaluate the full potential of mean force along $ F(z)$ in the interval $ [z_0,z_t]$, with $ 0<t<\tau$. However the exponential averages in Eq. 8.5 is known to be strongly biased, i.e. it contains a systematic error[149] that grows with decreasing number of non equilibrium experiments. This can be qualitatively explained with the fact that, for dissipative fast non equilibrium experiments, the forward work distribution $ P(W)$ has its maximum where the exponential factor $ e^{-\beta W}$ is negligibly small, so that the size of the integrand $ P(W)e^{-\beta W}$ is de facto controlled by the left tail of the $ P(W)$ distribution.[65] An unfortunate consequence of this, is that the PMF calculated through the JI becomes more and more biased as the reaction $ z$ coordinate is advanced, since the accumulated dissipation work shift the maximum of the $ P(W)$ distribution

The CT is far more precise than the JI to evaluate free energy differences. Shirts and Pande[66] have restated the CT theorem showing that the maximum likelihood estimate (MLE) of the free energy difference exactly correspond to the so-called Bennett acceptance ratio[118]8.4. The MLE restatement of the CT is the following

$\displaystyle \sum_{i=1}^{n_F} \frac{1}{ 1 + \frac{n_F}{n_R} ~ {\rm e}^{ \beta ...
...1}{ 1 + \frac{n_R}{n_F} ~ {\rm e}^{ \beta (
W[{\mathbf R}_i] + \Delta F ) } } -$     (8.6)

where the $ n_F$, $ n_r$ are the number of forward and backward non equilibrium experiments and $ W[{\mathbf F}_i]$ $ W[{\mathbf R}_i]$ indicate the outcome of i-th forward and backward work measurement. This equation has only one solution for $ \Delta F$, i.e. the MLE. As such, however, the Crooks theorem allows, through the MLE estimate based on bidirectional work measurements, to compute the free energy difference $ \Delta F$ between the end points (i.e. between thermodynamic states at fixed and given reaction coordinates $ z=z_0$ and $ z=z_t$). In principle, to reconstruct the full PMF along the reaction coordinate $ z$, in the spirit of thermodynamics integration, One should provide a series of equilibrium ensembles of configurations at intermediate values of $ z_t$. Here, we briefly sketch out a methodology for reconstructing the full PMF in the segment $ [z_0,z_\tau]$ doing only the two work measurements from $ z_0$ to $ z_\tau$ and back. We first rewrite the Crooks equation, Eq. 8.2, as follows

$\displaystyle \rho_F(\Gamma) = \rho_R (\hat \Gamma) e^{\beta ( W - \Delta F ) },$ (8.7)

where $ \rho_F$, $ \rho_b$ are the probability to observe a particular trajectory $ \Gamma$ in the forward and reverse process, respectively and $ \hat \Gamma$ indicate the time trajectory taken with inverted time schedule. Eq. 8.7 trivially implies that
$\displaystyle <{\cal F}>_F$ $\displaystyle =$ $\displaystyle <\hat {\cal F} e^{\beta (W-\Delta F)}>_R$ (8.8)
$\displaystyle <{\cal F}>_R$ $\displaystyle =$ $\displaystyle <\hat {\cal F} e^{-\beta (W-\Delta F)}>_F$ (8.9)

where $ {\cal F} = {\cal F}(\Gamma)$, $ \hat {\cal F} = {\cal F}(\hat
\Gamma)$ is an arbitrary functional of the trajectory $ \Gamma$ and of its inverted time schedule counterpart $ \hat \Gamma$. Using Eq. 8.7, we thus can combine the direct estimate of $ \rho_F(\Gamma)$ with the indirect estimate of the same quantity obtained from $ \rho_R (\hat \Gamma)$. This latter, according to Eq. 8.7, must be unbiased with the weight factor corresponding to the exponential of the dissipated work in the forward measurement. If the direct and indirect (Eq. 8.7) estimates are done with $ n_F$ forward measurements and $ n_R$ reverse measurements, respectively, the optimal (minimum variance) combination of these two estimates of $ \rho_F(\Gamma)$ is done according to the WHAM formula[54]
$\displaystyle \rho_F(\Gamma)$ $\displaystyle =$ $\displaystyle \frac {n_F
\rho_F(\Gamma) + n_R \rho_R(\hat \Gamma)} {n_F + n_R
e^{-\beta(W-\Delta F)}}.$ (8.10)

Here W is the work done in the full $ \Gamma$ path from the end point at $ t=0$ to the end point at $ t=\tau$. We now calculate the average of the trajectory functional $ e^{-\beta W_0^{t}}$ at intermediate times $ 0<t<\tau$, using the optimized above density. Taking the average of this functional over forward ($ \Gamma$) and reverse ( $ \hat \Gamma)$ work measurements, exploiting the Jarzynski identity 8.5 in the form $ <e^{-\beta W_0^{t}}> = e^{-\beta (F(z=z_t) - F(z=z_0))} $ , using the fact that $ W$ is odd under time reversal and that $ W_0^t[\hat \Gamma]=
-W^\tau_{(\tau-t)}[\Gamma]$, we obtain the following estimate for the free energy at intermediate $ t$, with $ 0<t<\tau$:


$\displaystyle e^{-\beta(F_t
-F_0) }$ $\displaystyle =$ $\displaystyle \left \langle \frac{n_F e^{-\beta W_0^t} }
{n_F+n_Re^{-\beta(W-\D...
...R
e^{\beta W_{(\tau-t)}^\tau} } {n_F+n_Re^{\beta(W+\Delta F)}} \right \rangle_R$ (8.11)

This equation, due to Minh and Adib[125], allows to reconstruct the entire potential of mean force $ F_t-F_0$ along the reaction coordinate spanned during the bidirectional non equilibrium experiments of duration $ \tau $, no matter how fast the driven processes are done. Note that $ \Delta W=F_\tau-F_0$ and $ W$ in Eq. 8.11 are the forward free energy difference and work relative the end points, respectively.

For fast pulling experiments, i.e. when the dissipated work is large, it can be shown[150], that Eq. 8.11 reduces to

$\displaystyle e^{-\beta(F_t -F_0)} = < e^{-\beta W_0^t} >_F + e^{-\beta \Delta F } < e^{-\beta W_{\tau)}^{(t)}>_R }$ (8.12)

In both Eq. 8.12 and Eq. 8.11 one needs to know the free energy difference between the end points $ \Delta F$. An unbiased estimate of $ \Delta F$ is easily available through the Bennett acceptance ratio, Eq. 8.6.

procacci 2021-12-29