BioSimSpace.FreeEnergy

The FreeEnergy package contains tools to configure, run, and analyse relative free energy simulations.

Free-energy perturbation simulations require a System containing a merged molecule that can be perturbed between two molecular end states by use of an alchemical potential. To create merged molecules, please use the BioSimSpace.Align package.

Relative free-energy calculations require the simulation of two perturbations, typically referred to as legs. A potential of mean force (PMF) is computed for each leg, which can then be used to computed the relative free-energy difference. For generality and flexibility, BioSimSpace decouples the two legs, allowing the use of difference molecular simulation engines, protocols, and for legs to be re-used in different calculations.

Simulations are typically used to compute solvation (currently hydration only) or binding free-energies. In the examples that follow, merged refers to a perturbable molecule created by merging two ligands, ligA and ligB, merged_sol refers to the same perturbable molecule in solvent, and complex_sol is a solvated protein-ligand complex containing the same perturbable molecule. We assume that each molecule/system has been appropriately minimised and equlibrated.

To setup, run, and analyse a binding free-energy calculation:

import BioSimSpace as BSS

...

# Create two a protocol for the two legs of a binding free-energy simulation.
# Use more lambda windows for the "bound" leg.
protocol_bound = BSS.Protocol.FreeEnergy(num_lam=20)
protocol_free  = BSS.Protocol.FreeEnergy(num_lam=12)

# Setup the perturbations for each leg, using the SOMD engine. This will
# create all of the input files and simulation processes that are required.
fep_bound = BSS.FreeEnergy.Relative(complex_sol,
                                    protocol_bound,
                                    engine="somd",
                                    work_dir="ligA_ligB/bound")
fep_free  = BSS.FreeEnergy.Relative(merged_sol,
                                    protocol_bound,
                                    engine="somd",
                                    work_dir="ligA_ligB/free")

# Run all simulations for each leg. Note that the lambda windows are run
# sequentially, so this is a sub-optimal way of executing the simulation
# if you have access to HPC resources.

# Bound leg.
fep_bound.run()
fep_bound.wait()

# Free leg.
fep_free.run()
fep_free.wait()

# Analyse the simulation data from each leg, returning the PMF and overlap
# matrix.
pmf_bound, overlap_bound = fep_bound.analyse()
pmf_free,  overlap_free  = fep_free.analyse()

# Compute the relative free-energy difference.
free_nrg_binding = BSS.FreeEnergy.Relative.difference(pmf_bound, pmf_free)

Similarly, for a solvation free-energy calculation:

# Here we are assuming that we are using the same ligands, so will re-use
# the free leg from the previous example.

# Setup the perturbation for the vacuum leg using a default protocol.
fep_vacuum = BSS.FreeEnergy.Relative(merged.toSystem(),
                                     engine="somd",
                                     work_dir="ligA_ligB/vacuum")

# Run the simulations for the perturbation.
fep_vacuum.run()
fep_vacuum.wait()

# Analyse the simulation data.
pmf_vacuum, overlap_vacuum = fep_vacuum.analyse()

# Compute the relative free-energy difference.
free_nrg_solvation = BSS.FreeEnergy.Relative.difference(pmf_free, pmf_vacuum)

Since it is usually preferable to run simulations intensive simulation such as these on external HPC resources, the BioSimSpace.FreeEnergy package also provides support for only creating the input files that are needed by passing the setup_only=True argument. This saves the overhead of creating Process objects. The input files can then be copied to a remote server, with the indivual simulations curated in a job submission script. (We don’t yet provide support for configuring and writing submission scripts for you.)

To just setup the vacuum leg input files:

# Setup the input for the vacuum leg. No processes are created so the .run()
# method won't do anything.
fep_vacuum = BSS.FreeEnergy.Relative(merged.toSystem(),
                                     engine="somd",
                                     work_dir="ligA_ligB/vacuum",
                                     setup_only=True)

It is also possible to analyse existing simulation output directly by passing the path to a working directory to FreeEnergy.Relative.analyse:

pmf_vacuum, overlap_vacuum = BSS.FreeEnergy.Relative.analyse("ligA_ligB/vacuum")

Classes

Relative(system[, protocol, work_dir, …])

Class for configuring and running relative free-energy perturbation simulations.

Functions

engines()	List the supported molecular dynamics engines for running free-energy perturbation simulations.
getData([name, file_link, work_dir])	Return a link to a zip file containing the data files required for post-simulation analysis.