Unbiased Rare Event Sampling in Spatial Stochastic Systems Biology Models Using a Weighted Ensemble of Trajectories

Molecular simulations reveal intricate coupling between agonist-bound β-adrenergic receptors and G protein

G protein-coupled receptors (GPCRs) and G proteins transmit signals from hormones and neurotransmitters across cell membranes, initiating downstream signaling and modulating cellular behavior. Using advanced computer modeling and simulation, we identified atomistic-level structural, dynamic, and energetic mechanisms of norepinephrine (NE) and stimulatory G protein (Gs) interactions with β-adrenergic receptors (βARs), crucial GPCRs for heart function regulation and major drug targets. Our analysis revealed distinct binding behaviors of NE within β1AR and β2AR despite identical orthosteric binding pockets. β2AR had an additional binding site, explaining variations in NE binding affinities. Simulations showed significant differences in NE dissociation pathways and receptor interactions with the Gs. β1AR binds Gs more strongly, while β2AR induces greater conformational changes in the α subunit of Gs. Furthermore, GTP and GDP binding to Gs may disrupt coupling between NE and βAR, as well as between βAR and Gs. These findings may aid in designing precise βAR-targeted drugs.

Weighted ensemble: Recent mathematical developments

Weighted ensemble (WE) is an enhanced sampling method based on periodically replicating and pruning trajectories generated in parallel. WE has grown increasingly popular for computational biochemistry problems due, in part, to improved hardware and accessible software implementations. Algorithmic and analytical improvements have played an important role, and progress has accelerated in recent years. Here, we discuss and elaborate on the WE method from a mathematical perspective, highlighting recent results that enhance the computational efficiency. The mathematical theory reveals a new strategy for optimizing trajectory management that approaches the best possible variance while generalizing to systems of arbitrary dimension.

A Suite of Tutorials for the WESTPA 2.0 Rare-Events Sampling Software [Article v2.0]

The weighted ensemble (WE) strategy has been demonstrated to be highly efficient in generating pathways and rate constants for rare events such as protein folding and protein binding using atomistic molecular dynamics simulations. Here we present two sets of tutorials instructing users in the best practices for preparing, carrying out, and analyzing WE simulations for various applications using the WESTPA software. The first set of more basic tutorials describes a range of simulation types, from a molecular association process in explicit solvent to more complex processes such as host-guest association, peptide conformational sampling, and protein folding. The second set ecompasses six advanced tutorials instructing users in the best practices of using key new features and plugins/extensions of the WESTPA 2.0 software package, which consists of major upgrades for larger systems and/or slower processes. The advanced tutorials demonstrate the use of the following key features: (i) a generalized resampler module for the creation of "binless" schemes, (ii) a minimal adaptive binning scheme for more efficient surmounting of free energy barriers, (iii) streamlined handling of large simulation datasets using an HDF5 framework, (iv) two different schemes for more efficient rate-constant estimation, (v) a Python API for simplified analysis of WE simulations, and (vi) plugins/extensions for Markovian Weighted Ensemble Milestoning and WE rule-based modeling for systems biology models. Applications of the advanced tutorials include atomistic and non-spatial models, and consist of complex processes such as protein folding and the membrane permeability of a drug-like molecule. Users are expected to already have significant experience with running conventional molecular dynamics or systems biology simulations.

Identifying a Feasible Transition Pathway between Two Conformational States for a Protein

Proteins usually need to transit between different conformational states to fulfill their biological functions. In the mechanistic study of such transition processes by molecular dynamics simulations, identification of the minimum free energy path (MFEP) can substantially reduce the sampling space, thus enabling rigorous thermodynamic evaluation of the process. Conventionally, the MFEP is derived by iterative local optimization from an initial path, which is typically generated by simple brute force techniques like the targeted molecular dynamics (tMD). Therefore, the quality of the initial path determines the successfulness of MFEP estimation. In this work, we propose a method to improve derivation of the initial path. Through iterative relaxation-biasing simulations in a bidirectional manner, this method can construct a feasible transition pathway connecting two known states for a protein. Evaluation on small, fast-folding proteins against long equilibrium trajectories supports the good sampling efficiency of our method. When applied to larger proteins including the catalytic domain of human c-Src kinase as well as the converter domain of myosin VI, the paths generated by our method deviate significantly from those computed with the generic tMD approach. More importantly, free energy profiles and intermediate states obtained from our paths exhibit remarkable improvements over those from tMD paths with respect to both physical rationality and consistency with a priori knowledge.

Sampling rare trajectories using stochastic bridges

The numerical quantification of the statistics of rare events in stochastic processes is a challenging computational problem. We present a sampling method that constructs an ensemble of stochastic trajectories that are constrained to have fixed start and end points (so-called stochastic bridges). We then show that by carefully choosing a set of such bridges and assigning an appropriate statistical weight to each bridge, one can focus more processing power on the rare events of a target stochastic process while faithfully preserving the statistics of these rare trajectories. Further, we also compare the stochastic bridges we produce to the Wentzel-Kramers-Brillouin (WKB) optimal paths of the target process, derived in the limit of low noise. We see that the generated paths, encoding the full statistics of the process, collapse onto the WKB optimal path as the level of noise is reduced. We propose that the method can also be used to judge the accuracy of the WKB approximation at finite levels of noise.

Simulation of receptor triggering by kinetic segregation shows role of oligomers and close-contacts

The activation of T cells, key players of the immune system, involves local evacuation of phosphatase CD45 from a region of the T cell's surface, segregating it from the T cell receptor. What drives this evacuation? In the presence of antigen, what ensures evacuation happens in the subsecond timescales necessary to initiate signaling? In the absence of antigen, what mechanisms ensure that evacuation does not happen spontaneously, which could cause signaling errors? Phenomena known to influence spatial organization of CD45 or similar surface molecules include diffusive motion in the lipid bilayer, oligomerization reactions, and mechanical compression against a nearby surface, such as that of the cell presenting the antigen. Computer simulations can investigate hypothesized spatiotemporal mechanisms of T cell signaling. The challenge to computational studies of evacuation is that the base process, spontaneous evacuation by simple diffusion, is in the extreme rare event limit, meaning direct stochastic simulation is unfeasible. Here, we combine particle-based spatial stochastic simulation with the weighted ensemble method for rare events to compute the mean first passage time for cell surface availability by surface reorganization of CD45. We confirm mathematical estimates that, at physiological concentrations, spontaneous evacuation is extremely rare, roughly 300 years. We find that dimerization decreases the time required for evacuation. A weak bimolecular interaction (dissociation constant estimate 460 μM) is sufficient for an order of magnitude reduction of spontaneous evacuation times, and oligomerization to hexamers reduces times to below 1 s. This introduces a mechanism whereby externally induced CD45 oligomerization could significantly modify T cell function. For large regions of close contact, such as those induced by large microvilli, molecular size and compressibility imply a nonzero reentry probability of 60%, decreasing evacuation times. Simulations show that these reduced evacuation times are still unrealistically long (even with a fourfold variation centered around previous estimates of parameters), suggesting that a yet-to-be-described mechanism, besides compressional exclusion at a close contact, drives evacuation.

WESTPA 2.0: High-Performance Upgrades for Weighted Ensemble Simulations and Analysis of Longer-Timescale Applications

The weighted ensemble (WE) family of methods is one of several statistical mechanics-based path sampling strategies that can provide estimates of key observables (rate constants and pathways) using a fraction of the time required by direct simulation methods such as molecular dynamics or discrete-state stochastic algorithms. WE methods oversee numerous parallel trajectories using intermittent overhead operations at fixed time intervals, enabling facile interoperability with any dynamics engine. Here, we report on the major upgrades to the WESTPA software package, an open-source, high-performance framework that implements both basic and recently developed WE methods. These upgrades offer substantial improvements over traditional WE methods. The key features of the new WESTPA 2.0 software enhance the efficiency and ease of use: an adaptive binning scheme for more efficient surmounting of large free energy barriers, streamlined handling of large simulation data sets, exponentially improved analysis of kinetics, and developer-friendly tools for creating new WE methods, including a Python API and resampler module for implementing both binned and "binless" WE strategies.

A gentle introduction to the non-equilibrium physics of trajectories: Theory, algorithms, and biomolecular applications

Despite the importance of non-equilibrium statistical mechanics in modern physics and related fields, the topic is often omitted from undergraduate and core-graduate curricula. Key aspects of non-equilibrium physics, however, can be understood with a minimum of formalism based on a rigorous trajectory picture. The fundamental object is the ensemble of trajectories, a set of independent time-evolving systems, which easily can be visualized or simulated (e.g., for protein folding) and which can be analyzed rigorously in analogy to an ensemble of static system configurations. The trajectory picture provides a straightforward basis for understanding first-passage times, "mechanisms" in complex systems, and fundamental constraints on the apparent reversibility of complex processes. Trajectories make concrete the physics underlying the diffusion and Fokker-Planck partial differential equations. Last but not least, trajectory ensembles underpin some of the most important algorithms that have provided significant advances in biomolecular studies of protein conformational and binding processes.

A first-passage approach to diffusion-influenced reversible binding and its insights into nanoscale signaling at the presynapse

Synaptic transmission between neurons is governed by a cascade of stochastic calcium ion reaction–diffusion events within nerve terminals leading to vesicular release of neurotransmitter. Since experimental measurements of such systems are challenging due to their nanometer and sub-millisecond scale, numerical simulations remain the principal tool for studying calcium-dependent neurotransmitter release driven by electrical impulses, despite the limitations of time-consuming calculations. In this paper, we develop an analytical solution to rapidly explore dynamical stochastic reaction–diffusion problems based on first-passage times. This is the first analytical model that accounts simultaneously for relevant statistical features of calcium ion diffusion, buffering, and its binding/unbinding reaction with a calcium sensor for synaptic vesicle fusion. In particular, unbinding kinetics are shown to have a major impact on submillisecond sensor occupancy probability and therefore cannot be neglected. Using Monte Carlo simulations we validated our analytical solution for instantaneous calcium influx and that through voltage-gated calcium channels. We present a fast and rigorous analytical tool that permits a systematic exploration of the influence of various biophysical parameters on molecular interactions within cells, and which can serve as a building block for more general cell signaling simulators.

Quantifying the roles of space and stochasticity in computer simulations for cell biology and cellular biochemistry

Most of the fascinating phenomena studied in cell biology emerge from interactions among highly organized multimolecular structures embedded into complex and frequently dynamic cellular morphologies. For the exploration of such systems, computer simulation has proved to be an invaluable tool, and many researchers in this field have developed sophisticated computational models for application to specific cell biological questions. However, it is often difficult to reconcile conflicting computational results that use different approaches to describe the same phenomenon. To address this issue systematically, we have defined a series of computational test cases ranging from very simple to moderately complex, varying key features of dimensionality, reaction type, reaction speed, crowding, and cell size. We then quantified how explicit spatial and/or stochastic implementations alter outcomes, even when all methods use the same reaction network, rates, and concentrations. For simple cases, we generally find minor differences in solutions of the same problem. However, we observe increasing discordance as the effects of localization, dimensionality reduction, and irreversible enzymatic reactions are combined. We discuss the strengths and limitations of commonly used computational approaches for exploring cell biological questions and provide a framework for decision making by researchers developing new models. As computational power and speed continue to increase at a remarkable rate, the dream of a fully comprehensive computational model of a living cell may be drawing closer to reality, but our analysis demonstrates that it will be crucial to evaluate the accuracy of such models critically and systematically.

Wepy: A Flexible Software Framework for Simulating Rare Events with Weighted Ensemble Resampling

Here, we introduce the open-source software framework wepy (https://github.com/ADicksonLab/wepy) which is a toolkit for running and analyzing weighted ensemble (WE) simulations. The wepy toolkit is in pure Python and as such is highly portable and extensible, making it an excellent platform to develop and use new WE resampling algorithms such as WExplore, REVO, and others while leveraging the entire Python ecosystem. In addition, wepy simplifies WE-specific analyses by defining out-of-core tree-like data structures using the cross-platform HDF5 file format. In this paper, we discuss the motivations and challenges for simulating rare events in biomolecular systems. As has previously been shown, high-dimensional WE resampling algorithms such as WExplore and REVO have been successful at these tasks, especially for rare events that are difficult to describe by one or two collective variables. We explain in detail how wepy facilitates implementation of these algorithms, as well as aids in analyzing the unique structure of WE simulation results. To explain how wepy and WE work in general, we describe the mathematical formalism of WE, an overview of the architecture of wepy, and provide code examples of how to construct, run, and analyze simulation results for a protein-ligand system (T4 Lysozyme in an implicit solvent). This paper is written with a variety of readers in mind, including (1) those curious about how to leverage WE rare-event simulations for their domain, (2) current WE users who want to begin using new high-dimensional resamplers such as WExplore and REVO, and (3) expert users who would like to prototype or implement their own algorithms that can be easily adopted by others.

Accelerated Estimation of Long-Timescale Kinetics from Weighted Ensemble Simulation via Non-Markovian "Microbin" Analysis

The weighted ensemble (WE) simulation strategy provides unbiased sampling of nonequilibrium processes, such as molecular folding or binding, but the extraction of rate constants relies on characterizing steady-state behavior. Unfortunately, WE simulations of sufficiently complex systems will not relax to steady state on observed simulation times. Here, we show that a postsimulation clustering of molecular configurations into "microbins" using methods developed in the Markov State Model (MSM) community can yield unbiased kinetics from WE data before steady-state convergence of the WE simulation itself. Because WE trajectories are directional and not equilibrium distributed, the history-augmented MSM (haMSM) formulation can be used, which yields the mean first-passage time (MFPT) without bias for arbitrarily small lag times. Accurate kinetics can be obtained while bypassing the often prohibitive convergence requirements of the nonequilibrium weighted ensemble. We validate the method in a simple diffusive process on a two-dimensional (2D) random energy landscape and then analyze atomistic protein folding simulations using WE molecular dynamics. We report significant progress toward the unbiased estimation of protein folding times and pathways, though key challenges remain.

Automatic error control during forward flux sampling of rare events in master equation models

Enhanced sampling methods, such as forward flux sampling (FFS), have great capacity for accelerating stochastic simulations of nonequilibrium biochemical systems involving rare events. However, the description of the tradeoffs between simulation efficiency and error in FFS remains incomplete. We present a novel and mathematically rigorous analysis of the errors in FFS that, for the first time, covers the contribution of every phase of the simulation. We derive a closed form expression for the optimally efficient count of samples to take in each FFS phase in terms of a fixed constraint on sampling error. We introduce a new method, forward flux pilot sampling (FFPilot), that is designed to take full advantage of our optimizing equation without prior information or assumptions about the phase weights and costs along the transition path. In simulations of both single and multidimensional gene regulatory networks, FFPilot is able to completely control sampling error. We then discuss how memory effects can introduce additional error when relaxation along the transition path is slow. This extra error can be traced to correlations between the FFS phases and can be controlled by monitoring the covariance between them. Finally, we show that, in sets of simulations with matched error, FFPilot is on the order of tens-to-hundreds of times faster than direct sampling and noticeably more efficient than previous FFS methods.

Potential based, spatial simulation of dynamically nested particles

BACKGROUND

To study cell biological phenomena which depend on diffusion, active transport processes, or the locations of species, modeling and simulation studies need to take space into account. To describe the system as a collection of discrete objects moving and interacting in continuous space, various particle-based reaction diffusion simulators for cell-biological system have been developed. So far the focus has been on particles as solid spheres or points. However, spatial dynamics might happen at different organizational levels, such as proteins, vesicles or cells with interrelated dynamics which requires spatial approaches that take this multi-levelness of cell biological systems into account.

RESULTS

Based on the perception of particles forming hollow spheres, ML-Force contributes to the family of particle-based simulation approaches: in addition to excluded volumes and forces, it also supports compartmental dynamics and relating dynamics between different organizational levels explicitly. Thereby, compartmental dynamics, e.g., particles entering and leaving other particles, and bimolecular reactions are modeled using pair-wise potentials (forces) and the Langevin equation. In addition, forces that act independently of other particles can be applied to direct the movement of particles. Attributes and the possibility to define arbitrary functions on particles, their attributes and content, to determine the results and kinetics of reactions add to the expressiveness of ML-Force. Its implementation comprises a rudimentary rule-based embedded domain-specific modeling language for specifying models and a simulator for executing models continuously. Applications inspired by cell biological models from literature, such as vesicle transport or yeast growth, show the value of the realized features. They facilitate capturing more complex spatial dynamics, such as the fission of compartments or the directed movement of particles, and enable the integration of non-spatial intra-compartmental dynamics as stochastic events.

CONCLUSIONS

By handling all dynamics based on potentials (forces) and the Langevin equation, compartmental dynamics, such as dynamic nesting, fusion and fission of compartmental structures are handled continuously and are seamlessly integrated with traditional particle-based reaction-diffusion dynamics within the cell. Thereby, attributes and arbitrary functions allow to flexibly describe diverse spatial phenomena, and relate dynamics across organizational levels. Also they prove crucial in modeling intra-cellular or intra-compartmental dynamics in a non-spatial manner, and, thus, to abstract from spatial dynamics, on demand which increases the range of multi-compartmental processes that can be captured.

A Suite of Tutorials for the WESTPA Rare-Events Sampling Software [Article v1.0]

The weighted ensemble (WE) strategy has been demonstrated to be highly efficient in generating pathways and rate constants for rare events such as protein folding and protein binding using atomistic molecular dynamics simulations. Here we present five tutorials instructing users in the best practices for preparing, carrying out, and analyzing WE simulations for various applications using the WESTPA software. Users are expected to already have significant experience with running standard molecular dynamics simulations using the underlying dynamics engine of interest (e.g. Amber, Gromacs, OpenMM). The tutorials range from a molecular association process in explicit solvent to more complex processes such as host-guest association, peptide conformational sampling, and protein folding.

Computational Estimation of Microsecond to Second Atomistic Folding Times

Despite the development of massively parallel computing hardware including inexpensive graphics processing units (GPUs), it has remained infeasible to simulate the folding of atomistic proteins at room temperature using conventional molecular dynamics (MD) beyond the microsecond scale. Here, we report the folding of atomistic, implicitly solvated protein systems with folding times τ ranging from ∼10 μs to ∼100 ms using the weighted ensemble (WE) strategy in combination with GPU computing. Starting from an initial structure or set of structures, WE organizes an ensemble of GPU-accelerated MD trajectory segments via intermittent pruning and replication events to generate statistically unbiased estimates of rate constants for rare events such as folding; no biasing forces are used. Although the variance among atomistic WE folding runs is significant, multiple independent runs are used to reduce and quantify statistical uncertainty. Folding times are estimated directly from WE probability flux and from history-augmented Markov analysis of the WE data. Three systems were examined: NTL9 at low solvent viscosity (yielding τf = 0.8-9 μs), NTL9 at water-like viscosity (τf = 0.2-2 ms), and Protein G at low viscosity (τf = 3-200 ms). In all cases, the folding time, uncertainty, and ensemble properties could be estimated from WE simulation; for Protein G, this characterization required significantly less overall computing than would be required to observe a single folding event with conventional MD simulations. Our results suggest that the use and calibration of force fields and solvent models for precise estimation of kinetic quantities is becoming feasible.

MCell-R: A Particle-Resolution Network-Free Spatial Modeling Framework

Spatial heterogeneity can have dramatic effects on the biochemical networks that drive cell regulation and decision-making. For this reason, a number of methods have been developed to model spatial heterogeneity and incorporated into widely used modeling platforms. Unfortunately, the standard approaches for specifying and simulating chemical reaction networks become untenable when dealing with multistate, multicomponent systems that are characterized by combinatorial complexity. To address this issue, we developed MCell-R, a framework that extends the particle-based spatial Monte Carlo simulator, MCell, with the rule-based model specification and simulation capabilities provided by BioNetGen and NFsim. The BioNetGen syntax enables the specification of biomolecules as structured objects whose components can have different internal states that represent such features as covalent modification and conformation and which can bind components of other molecules to form molecular complexes. The network-free simulation algorithm used by NFsim enables efficient simulation of rule-based models even when the size of the network implied by the biochemical rules is too large to enumerate explicitly, which frequently occurs in detailed models of biochemical signaling. The result is a framework that can efficiently simulate systems characterized by combinatorial complexity at the level of spatially resolved individual molecules over biologically relevant time and length scales.

Sampling rare events in stochastic reaction-diffusion systems within trajectory looping

In bistable reaction-diffusion systems, transitions between stable states typically occur on timescales orders of magnitude longer than the chemical equilibration time. Estimation of transition rates within explicit Brownian dynamics simulations is computationally prohibitively costly. We present a method that exploits a single trajectory, generated by a prior simulation of diffusive motions of molecules, to sample chemical kinetic processes on timescales several orders of magnitude longer than the duration of the diffusive trajectory. In this approach, we "loop" the diffusive trajectory by transferring chemical states of the molecules from the last to the first time step of the trajectory. Trajectory looping can be applied to enhance sampling of rare events in biochemical systems in which the number of reacting molecules is constant, as in cellular signal transduction pathways. As an example, we consider a bistable system of autophosphorylating kinases, for which we calculate state-to-state transition rates and traveling wave velocities. We provide an open-source implementation of the method.

Multifidelity Analysis for Predicting Rare Events in Stochastic Computational Models of Complex Biological Systems

Rare events such as genetic mutations or cell-cell interactions are important contributors to dynamics in complex biological systems, eg, in drug-resistant infections. Computational approaches can help analyze rare events that are difficult to study experimentally. However, analyzing the frequency and dynamics of rare events in computational models can also be challenging due to high computational resource demands, especially for high-fidelity stochastic computational models. To facilitate analysis of rare events in complex biological systems, we present a multifidelity analysis approach that uses medium-fidelity analysis (Monte Carlo simulations) and/or low-fidelity analysis (Markov chain models) to analyze high-fidelity stochastic model results. Medium-fidelity analysis can produce large numbers of possible rare event trajectories for a single high-fidelity model simulation. This allows prediction of both rare event dynamics and probability distributions at much lower frequencies than high-fidelity models. Low-fidelity analysis can calculate probability distributions for rare events over time for any frequency by updating the probabilities of the rare event state space after each discrete event of the high-fidelity model. To validate the approach, we apply multifidelity analysis to a high-fidelity model of tuberculosis disease. We validate the method against high-fidelity model results and illustrate the application of multifidelity analysis in predicting rare event trajectories, performing sensitivity analyses and extrapolating predictions to very low frequencies in complex systems. We believe that our approach will complement ongoing efforts to enable accurate prediction of rare event dynamics in high-fidelity computational models.

Rare-event sampling of epigenetic landscapes and phenotype transitions

Stochastic simulation has been a powerful tool for studying the dynamics of gene regulatory networks, particularly in terms of understanding how cell-phenotype stability and fate-transitions are impacted by noisy gene expression. However, gene networks often have dynamics characterized by multiple attractors. Stochastic simulation is often inefficient for such systems, because most of the simulation time is spent waiting for rare, barrier-crossing events to occur. We present a rare-event simulation-based method for computing epigenetic landscapes and phenotype-transitions in metastable gene networks. Our computational pipeline was inspired by studies of metastability and barrier-crossing in protein folding, and provides an automated means of computing and visualizing essential stationary and dynamic information that is generally inaccessible to conventional simulation. Applied to a network model of pluripotency in Embryonic Stem Cells, our simulations revealed rare phenotypes and approximately Markovian transitions among phenotype-states, occurring with a broad range of timescales. The relative probabilities of phenotypes and the transition paths linking pluripotency and differentiation are sensitive to global kinetic parameters governing transcription factor-DNA binding kinetics. Our approach significantly expands the capability of stochastic simulation to investigate gene regulatory network dynamics, which may help guide rational cell reprogramming strategies. Our approach is also generalizable to other types of molecular networks and stochastic dynamics frameworks.

ML-Space: Hybrid Spatial Gillespie and Particle Simulation of Multi-Level Rule-Based Models in Cell Biology

Spatio-temporal dynamics of cellular processes can be simulated at different levels of detail, from (deterministic) partial differential equations via the spatial Stochastic Simulation algorithm to tracking Brownian trajectories of individual particles. We present a spatial simulation approach for multi-level rule-based models, which includes dynamically hierarchically nested cellular compartments and entities. Our approach ML-Space combines discrete compartmental dynamics, stochastic spatial approaches in discrete space, and particles moving in continuous space. The rule-based specification language of ML-Space supports concise and compact descriptions of models and to adapt the spatial resolution of models easily.

Weighted Ensemble Simulation: Review of Methodology, Applications, and Software

The weighted ensemble (WE) methodology orchestrates quasi-independent parallel simulations run with intermittent communication that can enhance sampling of rare events such as protein conformational changes, folding, and binding. The WE strategy can achieve superlinear scaling-the unbiased estimation of key observables such as rate constants and equilibrium state populations to greater precision than would be possible with ordinary parallel simulation. WE software can be used to control any dynamics engine, such as standard molecular dynamics and cell-modeling packages. This article reviews the theoretical basis of WE and goes on to describe successful applications to a number of complex biological processes-protein conformational transitions, (un)binding, and assembly processes, as well as cell-scale processes in systems biology. We furthermore discuss the challenges that need to be overcome in the next phase of WE methodological development. Overall, the combined advances in WE methodology and software have enabled the simulation of long-timescale processes that would otherwise not be practical on typical computing resources using standard simulation.

Path-sampling strategies for simulating rare events in biomolecular systems

Despite more than three decades of effort with molecular dynamics simulations, long-timescale (ms and beyond) biologically relevant phenomena remain out of reach in most systems of interest. This is largely because important transitions, such as conformational changes and (un)binding events, tend to be rare for conventional simulations (<10μs). That is, conventional simulations will predominantly dwell in metastable states instead of making large transitions in complex biomolecular energy landscapes. In contrast, path sampling approaches focus computing effort specifically on transitions of interest. Such approaches have been in use for nearly 20 years in biomolecular systems and enabled the generation of pathways and calculation of rate constants for ms processes, including large protein conformational changes, protein folding, and protein (un)binding.

Efficient Atomistic Simulation of Pathways and Calculation of Rate Constants for a Protein-Peptide Binding Process: Application to the MDM2 Protein and an Intrinsically Disordered p53 Peptide

The characterization of protein binding processes - with all of the key conformational changes - has been a grand challenge in the field of biophysics. Here, we have used the weighted ensemble path sampling strategy to orchestrate molecular dynamics simulations, yielding atomistic views of protein-peptide binding pathways involving the MDM2 oncoprotein and an intrinsically disordered p53 peptide. A total of 182 independent, continuous binding pathways were generated, yielding a kon that is in good agreement with experiment. These pathways were generated in 15 days using 3500 cores of a supercomputer, substantially faster than would be possible with "brute force" simulations. Many of these pathways involve the anchoring of p53 residue F19 into the MDM2 binding cleft when forming the metastable encounter complex, indicating that F19 may be a kinetically important residue. Our study demonstrates that it is now practical to generate pathways and calculate rate constants for protein binding processes using atomistic simulation on typical computing resources.