Markov chain Monte Carlo (MCMC)¶
openmmtools
provides an extensible Markov chain Monte Carlo simulation framework.
This module provides a framework for equilibrium sampling from a given thermodynamic state of a biomolecule using a Markov chain Monte Carlo scheme.
It currently offer supports for
Langevin dynamics (assumed to be free of integration error; use at your own risk]),
hybrid Monte Carlo,
generalized hybrid Monte Carlo, and
Monte Carlo barostat moves,
which can be combined through the SequenceMove
and WeightedMove
classes.
By default, the MCMCMoves
use the fastest OpenMM platform available and a shared global ContextCache
that minimizes the number of OpenMM Context
objects that must be maintained at once.
The examples below show how to configure these aspects.
Note
To use the ContextCache
on the CUDA platform, the NVIDIA driver must be set to shared
mode to allow the process to create multiple GPU contexts.
Using the MCMC framework requires importing ThermodynamicState
and SamplerState
from openmmtools.states
:
from simtk import unit
from openmmtools import testsystems, cache
from openmmtools.states import ThermodynamicState, SamplerState
Create the initial state (thermodynamic and microscopic) for an alanine dipeptide system in vacuum.
test = testsystems.AlanineDipeptideVacuum()
thermodynamic_state = ThermodynamicState(system=test.system, temperature=298*unit.kelvin)
sampler_state = SamplerState(positions=test.positions)
Create an MCMC move to perform at every iteration of the simulation, and initialize a sampler instance.
ghmc_move = GHMCMove(timestep=1.0*unit.femtosecond, n_steps=50)
langevin_move = LangevinDynamicsMove(n_steps=10)
sampler = MCMCSampler(thermodynamic_state, sampler_state, move=ghmc_move)
You can combine them to form a sequence of moves
sequence_move = SequenceMove([ghmc_move, langevin_move])
sampler = MCMCSampler(thermodynamic_state, sampler_state, move=sequence_move)
or create a move that selects one of them at random with given probability at each iteration.
weighted_move = WeightedMove([(ghmc_move, 0.5), (langevin_move, 0.5)])
sampler = MCMCSampler(thermodynamic_state, sampler_state, move=weighted_move)
By default the MCMCMove
use a global ContextCache that creates Context
on the
fastest available OpenMM platform. You can configure the default platform to use
before starting the simulation
reference_platform = openmm.Platform.getPlatformByName('Reference')
cache.global_context_cache.platform = reference_platform
cache.global_context_cache.time_to_live = 10 # number of read/write operations
Minimize and run the simulation for few iterations.
sampler.minimize()
sampler.run(n_iterations=2)
If you don’t want to use a global cache, you can create local ones.
local_cache1 = cache.ContextCache(capacity=5, time_to_live=50)
local_cache2 = cache.ContextCache(platform=reference_platform, capacity=1)
sequence_move = SequenceMove([HMCMove(), LangevinDynamicsMove()], context_cache=local_cache1)
ghmc_move = GHMCMove(context_cache=local_cache2)
If you don’t want to cache Context
at all but create one every time, you can use
the DummyCache
.
dummy_cache = cache.DummyContextCache(platform=reference_platform)
ghmc_move = GHMCMove(context_cache=dummy_cache)
This book by Jun Liu is an excellent overview of Markov chain Monte Carlo:
Jun S. Liu. Monte Carlo Strategies in Scientific Computing. Springer, 2008.
MCMC samplers¶
An MCMC sampler driver is provided that can either utilize a programmed sequence of moves or draw from a weighted set of moves.
Basic Markov chain Monte Carlo sampler. |
|
A sequence of MCMC moves. |
|
Pick an MCMC move out of set with given probability at each iteration. |
MCMC move types¶
A number of MCMC component move types that can be arranged into groups or subclassed are provided.
Markov chain Monte Carlo (MCMC) move abstraction. |
|
A general MCMC move that applies an integrator. |
|
A base class for metropolized moves. |
|
An MCMCMove that propagate the system with an integrator. |
|
Langevin dynamics segment as a (pseudo) Monte Carlo move. |
|
Langevin dynamics segment with custom splitting of the operators and optional Metropolized Monte Carlo validation. |
|
Generalized hybrid Monte Carlo (GHMC) Markov chain Monte Carlo move. |
|
Hybrid Monte Carlo dynamics. |
|
Monte Carlo barostat move. |
|
A metropolized move that randomly displace a subset of atoms. |
|
A metropolized move that randomly rotate a subset of atoms. |