^{1,a)}, Robert H. Swendsen

^{2}and Daniel M. Zuckerman

^{3,b)}

### Abstract

We present a detailed comparison of computational efficiency and precision for several free energy difference methods. The analysis includes both equilibrium and nonequilibrium approaches, and distinguishes between unidirectional and bidirectional methodologies. We are primarily interested in comparing two recently proposed approaches, adaptive integration, and single-ensemble path sampling to more established methodologies. As test cases, we study relative solvation free energies of large changes to the size or charge of a Lennard-Jones particle in explicit water. The results show that, for the systems used in this study, both adaptive integration and path sampling offer unique advantages over the more traditional approaches. Specifically, adaptive integration is found to provide very precise long-simulation estimates as compared to other methods used in this report, while also offering rapid estimation of . The results demonstrate that the adaptive integration approach is the best overall method for the systems studied here. The single-ensemble path sampling approach is found to be superior to ordinary Jarzynski averaging for the unidirectional, “fast-growth” nonequilibrium case. Closer examination of the path sampling approach on a two-dimensional system suggests it may be the overall method of choice when conformational sampling barriers are high. However, it appears that the free energy landscapes for the systems used in this study have rather modest configurational sampling barriers.

The authors would like to thank Ron White and Hagai Meirovitch for valuable discussions, and also Manuel Athènes and Gilles Adjanor for helpful comments regarding the manuscript. Funding for this research was provided by the Department of Computational Biology and the Department of Environmental and Occupational Health at the University of Pittsburgh, and the National Institutes of Health under Fellowship No. GM073517 [to one of the authors (F.M.Y.)] and Grant Nos. ES007318, GM073517, and CA078039.

I. INTRODUCTION

II. EQUILIBRIUM FREE ENERGY CALCULATION

A. Thermodynamic integration

B. Adaptive integration

C. Free energyperturbation

D. Equilibrium Bennett estimation

III. NONEQUILIBRIUM FREE ENERGY ESTIMATION

A. Jarzynski averaging

B. Bennett averaging of Jarzynski work values

C. Single-ensemble path sampling

D. Bennett averaging of path sampled work values

IV. SIMULATION DETAILS

A. Thermodynamic integration calculations

B. Adaptive integration calculations

C. Free energyperturbation and equilibrium Bennett calculations

D. Jarzynski estimate calculations

E. Single-ensemble path sampling calculations

V. RESULTS AND DISCUSSION

A. Growing a Lennard-Jones particle

1. Fast-growth unidirectional data

B. Charging a Lennard-Jones particle

C. A second look at a two-dimensional model

VI. CONCLUSIONS

### Key Topics

- Free energy
- 69.0
- Brownian dynamics
- 14.0
- Entropy
- 8.0
- Solvents
- 7.0
- Thermal analysis
- 5.0

## Figures

The slope of the free energy as a function of for changing the Lennard-Jones size of a neutral particle in a box of explicit water. Results for both TI and AIM methods are shown for dynamics steps. The data show the averages (data points) and standard deviations (error bars) from 16 independent simulations for each method. The figure demonstrates that AIM has the ability to sample the path more efficiently, thus producing a much smoother and more precise profile compared to TI. Thus, AIM is preferred over TI for computing the potential of mean force for this system. In addition, the smoothness of the profile suggests that the switching function of Eq. (25) used in this report is adequate.

The slope of the free energy as a function of for changing the Lennard-Jones size of a neutral particle in a box of explicit water. Results for both TI and AIM methods are shown for dynamics steps. The data show the averages (data points) and standard deviations (error bars) from 16 independent simulations for each method. The figure demonstrates that AIM has the ability to sample the path more efficiently, thus producing a much smoother and more precise profile compared to TI. Thus, AIM is preferred over TI for computing the potential of mean force for this system. In addition, the smoothness of the profile suggests that the switching function of Eq. (25) used in this report is adequate.

(a) “Fast-growth” unidirectional free energy difference estimates obtained for changing the Lennard-Jones size of a neutral particle in a box of explicit water. Results are shown for both SEPS and Jarz methods as a function of the number of dynamics steps used in the simulation. For both methods, fast-growth work values were generated by simulating roughly 2000 dynamics steps per path, which is ten times shorter than optimal. The solid horizontal line represents the best estimate of the free energy difference based on averaging all results shown in Table I at dynamics steps. The averages (data points) and standard deviations (errorbars) are from 16 independent simulations. (b) Histograms of the work values used to generate the free energy estimates for both the SEPS and Jarz methods. The plots demonstrate the potential usefulness of using path sampling over regular Jarzynski averaging. Specifically, if the work values are fast growth and unidirectional, then SEPS is able to bias the work values in such a way to improve the free energy estimate. Note that for all the SEPS data shown, the first 50 work values are thrown away for equilibration, as described in Sec. IV E.

(a) “Fast-growth” unidirectional free energy difference estimates obtained for changing the Lennard-Jones size of a neutral particle in a box of explicit water. Results are shown for both SEPS and Jarz methods as a function of the number of dynamics steps used in the simulation. For both methods, fast-growth work values were generated by simulating roughly 2000 dynamics steps per path, which is ten times shorter than optimal. The solid horizontal line represents the best estimate of the free energy difference based on averaging all results shown in Table I at dynamics steps. The averages (data points) and standard deviations (errorbars) are from 16 independent simulations. (b) Histograms of the work values used to generate the free energy estimates for both the SEPS and Jarz methods. The plots demonstrate the potential usefulness of using path sampling over regular Jarzynski averaging. Specifically, if the work values are fast growth and unidirectional, then SEPS is able to bias the work values in such a way to improve the free energy estimate. Note that for all the SEPS data shown, the first 50 work values are thrown away for equilibration, as described in Sec. IV E.

The slope of the free energy as a function of for a changing the charge of a Lennard-Jones particle in a box of explicit water from to . Results for both TI and AIM methods are shown for dynamics steps. The data show the averages (data points) and standard deviations (error bars) from 16 independent simulations for each method. The differences between TI and AIM are too small to resolve on the plot shown; however, it should be noted that the average uncertainty in the data for AIM is and for TI is , suggesting that AIM has the ability to produce a more precise profile compared to TI. Thus, AIM is preferred over TI for computing the potential of mean force for this system. The smoothness of the profile also suggests that the switching function of Eq. (24) used in this report is adequate.

The slope of the free energy as a function of for a changing the charge of a Lennard-Jones particle in a box of explicit water from to . Results for both TI and AIM methods are shown for dynamics steps. The data show the averages (data points) and standard deviations (error bars) from 16 independent simulations for each method. The differences between TI and AIM are too small to resolve on the plot shown; however, it should be noted that the average uncertainty in the data for AIM is and for TI is , suggesting that AIM has the ability to produce a more precise profile compared to TI. Thus, AIM is preferred over TI for computing the potential of mean force for this system. The smoothness of the profile also suggests that the switching function of Eq. (24) used in this report is adequate.

## Tables

Free energy difference estimates in units of kcal/mol obtained for changing the Lennard-Jones size of a neutral particle in a box of explicit water. Results are shown for various methods described in the text as a function of the number of dynamics steps used in the simulation. Table entries are the mean estimates from 16 independent simulations with the standard deviation shown in parentheses. For single-ensemble path sampling (SEPS and BSEPS) and Jarzynski methods (Jarz and BJarz), only the most efficient results are shown. The table shows that in the limit of long simulation times ( dynamics steps) all methods produce average estimates that roughly agree. The table also shows that AIM provides the most precise long-simulation estimate.

Free energy difference estimates in units of kcal/mol obtained for changing the Lennard-Jones size of a neutral particle in a box of explicit water. Results are shown for various methods described in the text as a function of the number of dynamics steps used in the simulation. Table entries are the mean estimates from 16 independent simulations with the standard deviation shown in parentheses. For single-ensemble path sampling (SEPS and BSEPS) and Jarzynski methods (Jarz and BJarz), only the most efficient results are shown. The table shows that in the limit of long simulation times ( dynamics steps) all methods produce average estimates that roughly agree. The table also shows that AIM provides the most precise long-simulation estimate.

Number of dynamics steps necessary to be within a specified tolerance of the correct result , average estimate at dynamics steps for all methods, for growing a Lennard-Jones particle in explicit solvent. The first column is the method used to obtain the estimate. The second column is the number of dynamics steps needed to estimate within of with an uncertainty less than . The third column is the number of dynamics steps needed to obtain an estimate within with an uncertainty less than .

Number of dynamics steps necessary to be within a specified tolerance of the correct result , average estimate at dynamics steps for all methods, for growing a Lennard-Jones particle in explicit solvent. The first column is the method used to obtain the estimate. The second column is the number of dynamics steps needed to estimate within of with an uncertainty less than . The third column is the number of dynamics steps needed to obtain an estimate within with an uncertainty less than .

Free energy difference estimates in units of kcal/mol obtained for changing the charge of a Lennard-Jones particle from to in a box of explicit water. Results are the averages from 16 independent simulations for various methods described in the text as a function of the number of dynamics steps used in the simulation. The standard deviation is shown in parentheses. For single-ensemble path sampling (SEPS and BSEPS) and Jarzynski methods (Jarz and BJarz), only the most efficient results are shown. The table shows that in the limit of long simulation times ( dynamics steps) all methods produce average estimates that roughly agree. The table also shows that AIM and BJarz approaches provide the most precise long-simulation estimate.

Free energy difference estimates in units of kcal/mol obtained for changing the charge of a Lennard-Jones particle from to in a box of explicit water. Results are the averages from 16 independent simulations for various methods described in the text as a function of the number of dynamics steps used in the simulation. The standard deviation is shown in parentheses. For single-ensemble path sampling (SEPS and BSEPS) and Jarzynski methods (Jarz and BJarz), only the most efficient results are shown. The table shows that in the limit of long simulation times ( dynamics steps) all methods produce average estimates that roughly agree. The table also shows that AIM and BJarz approaches provide the most precise long-simulation estimate.

Number of dynamics steps necessary to be within a specified tolerance of the correct result , average estimate at dynamics steps for all methods, for charging a Lennard-Jones particle in explicit solvent. The first column is the method used to obtain the estimate. The second column is the number of dynamics steps needed to estimate within of with an uncertainty less than . The third column is the number of dynamics steps needed to obtain an estimate within with an uncertainty less than .

Number of dynamics steps necessary to be within a specified tolerance of the correct result , average estimate at dynamics steps for all methods, for charging a Lennard-Jones particle in explicit solvent. The first column is the method used to obtain the estimate. The second column is the number of dynamics steps needed to estimate within of with an uncertainty less than . The third column is the number of dynamics steps needed to obtain an estimate within with an uncertainty less than .

Number of dynamics steps necessary to be within of the analytical result for with a or less standard deviation for the two-dimensional model in Ref. 10. The first column is the barrier height of the potential energy surface in units. The second and third columns are the total numbers of dynamics steps using SEPS with, respectively, 200 work values and 20 000 work values. The fourth column is the total number of dynamics steps using TI with using 51 equally spaced values of . For both TI and SEPS, half of the generated data was thrown away for equilibration.

Number of dynamics steps necessary to be within of the analytical result for with a or less standard deviation for the two-dimensional model in Ref. 10. The first column is the barrier height of the potential energy surface in units. The second and third columns are the total numbers of dynamics steps using SEPS with, respectively, 200 work values and 20 000 work values. The fourth column is the total number of dynamics steps using TI with using 51 equally spaced values of . For both TI and SEPS, half of the generated data was thrown away for equilibration.

Article metrics loading...

Full text loading...

Commenting has been disabled for this content