This package is for use with 'Wepy', and can be used to calculate (un)binding free energies from pulling trajectories of Lennard-Jones Pairs. This code was used to generate the data in N. M. Roussey & A. Dickson, 2020.
The theory for this project is based in nonequilibrium stastical mechanics, specifically the Jarzynski Nonequilibrium Work Relation AKA 'The Jarzynski Equality'. This is an equation that relates the work done by a nonequilibrium process on a system to the equilibrium free energy difference between the inital and final states. This equation is of interest because it allows us to perform nonequilibrium work (along some system-relevant lambda from lambda_0 to lambda_1) to generate equilibrium free energy differences. The JE is as follows: <exp(-beta W)> = exp(-beta* Delta F), where <...> indicate the ensemble average, W is work, and Delta F is the Helmholtz free energy. This is significant because these nonequilibrium simualtions can be very fast and potentially less computationly expensive than traditional equilibrium simualtions used to determine the same value.
The final free energy determined is dominated by low work / high importance trajectories, and we are interested in determining if enhancing the sampling of low work trajectories may allow us to find an accurate value of Delta F without the necessary number of trajectories otherwise needed for convergence. We have aimed to develop a new method of calculating free energy differences between states of a system utilizing nonequilibrium processes and means of analyzing these nonequilibrium processes. This software conatins a new "importance resampler" as well as a "history-dependent" version of the novelty-REVO resampler. A Wepy compatible Diffusion Monte Carlo engine is also provided. Free energy surfaces can be generated for these simulations with the provided analysis code.
Nicole Roussey(1), Samuel Lotz(1), & Alex Dickson(1)(2)
(1). The Department of Biochemistry and Molecular Biology, Michigan State University, East Lansing, Michigan
(2). The Department of Computational Mathemetaics, Science, & Engineering, Michigan State University, East Lansing, Michigan
To use this software:
First clone the necessary git repos.
git clone https://gitlab.com/ADicksonLab/wepy.git
install wepy with "pip install -e ."
git clone https://gitlab.com/nmroussey/lj_activity.git
conda install -c omnia -c conda-forge openmm
conda install -c conda-forge mdtraj
'Wepy' will install most dependencies.
Then to test the example scripts (containing values for 500 cycle simulations):
python lj_noresampler.py 500 100 50 test_noresampler.wepy.h5 10 #args are num_cycles num_steps num_walkers output_file_name epsilon
python lj_revo.py 500 100 50 test_revo.wepy.h5 10 char_dist* #args are num_cycles num_steps num_walkers output_file_name epsilon *characteristic distance (see: Donyapour, 2019)
*For epsilon = 10 kcal/mol, char_dist is 0.46396780
python lj_importance.py 500 100 50 test_importance.wepy.h5 10 1.0 #args are num_cycles num_steps num_walkers output_file_name epsilon amplification
python lj_diffmc.py 500 100 50 test_diffmc.wepy.h5 10 # args are num_cycles num_steps num_walkers output_file_name epsilon
The analysis provided here is equation 8 from Hummer and Szabo, 2001. This creates a free energy surface from the probabilities and work values saved during a simulation.
To run this analysis code:
-
run the REVO example above
-
cd analysis/
-
python hz_equation8.py 10 50 ../outputs/resampler/10_eps/test_revo.wepy.h5 #args are epsilon num_walkers path/to/file
A plot of the target free energy surface for epsilon 10 as well as the free energy surface generated from the simulation data will appear. The diffmc version of this analysis is available in the diff_mc folder and is run with the same inputs.
Software Packages
'Wepy'
'OpenMM'
'mdtraj'
'NumPy'
'h5py'
Papers
'Nonequilibrium Equality for Free Energy Differences' Jarzynski, C., PRL, 1997
Analysis from:
'Free energy reconstruction from nonequilibrium single-molecule pulling experiments' Hummer, G., & Szabo, A., PNAS, 2001