Heterogeneous path ensembles for conformational transitions in semi-atomistic models of adenylate kinase

Worth the Weight: Sub-Pocket EXplorer (SubPEx), a Weighted Ensemble Method to Enhance Binding-Pocket Conformational Sampling

Structure-based virtual screening (VS) is an effective method for identifying potential small-molecule ligands, but traditional VS approaches consider only a single binding-pocket conformation. Consequently, they struggle to identify ligands that bind to alternate conformations. Ensemble docking helps address this issue by incorporating multiple conformations into the docking process, but it depends on methods that can thoroughly explore pocket flexibility. We here introduce Sub-Pocket EXplorer (SubPEx), an approach that uses weighted ensemble (WE) path sampling to accelerate binding-pocket sampling. As proof of principle, we apply SubPEx to three proteins relevant to drug discovery: heat shock protein 90, influenza neuraminidase, and yeast hexokinase 2. SubPEx is available free of charge without registration under the terms of the open-source MIT license: http://durrantlab.com/subpex/.

Worth the weight: Sub-Pocket EXplorer (SubPEx), a weighted-ensemble method to enhance binding-pocket conformational sampling

Structure-based virtual screening (VS) is an effective method for identifying potential small-molecule ligands, but traditional VS approaches consider only a single binding-pocket conformation. Consequently, they struggle to identify ligands that bind to alternate conformations. Ensemble docking helps address this issue by incorporating multiple conformations into the docking process, but it depends on methods that can thoroughly explore pocket flexibility. We here introduce Sub-Pocket EXplorer (SubPEx), an approach that uses weighted ensemble (WE) path sampling to accelerate binding-pocket sampling. As proof of principle, we apply SubPEx to three proteins relevant to drug discovery: heat shock protein 90, influenza neuraminidase, and yeast hexokinase 2. SubPEx is available free of charge without registration under the terms of the open-source MIT license: http://durrantlab.com/subpex/.

Simulation of ligand dissociation kinetics from the protein kinase PYK2

Early-stage drug discovery projects often focus on equilibrium binding affinity to the target alongside selectivity and other pharmaceutical properties. The kinetics of drug binding are ignored but can have significant influence on drug efficacy. Therefore, increasing attention has been paid on evaluating drug-binding kinetics early in a drug discovery process. Simulating drug-binding kinetics at the atomic level is challenging for the long time scale involved. Here, we used the transition-based reweighting analysis method (TRAM) with the Markov state model to study the dissociation of a ligand from the protein kinase PYK2. TRAM combines biased and unbiased simulations to reduce computational costs. This work used the umbrella sampling technique for the biased simulations. Although using the potential of mean force from umbrella sampling simulations with the transition-state theory over-estimated the dissociation rate by three orders of magnitude, TRAM gave a dissociation rate within an order of magnitude of the experimental value.

Understanding enzyme behavior in a crowded scenario through modulation in activity, conformation and dynamics

Macromolecular crowding, inside the physiological interior, modulates the energy landscape of biological macromolecules in multiple ways. Amongst these, enzymes occupy a special place and hence understanding the function of the same in the crowded interior is of utmost importance. In this study, we have investigated the manner in which the multidomain enzyme, AK3L1 (PDB ID: 1ZD8), an isoform of adenylate kinase, has its features affected in presence of commonly used crowders (PEG 8, Dextran 40, Dextran 70, and Ficoll 70). Michaelis Menten plots reveal that the crowders in general enhance the activity of the enzyme, with the K_m and V_max values showing significant variations. Ficoll 70, induced the maximum activity for AK3L1 at 100 g/L, beyond which the activity reduced. Ensemble FRET studies were performed to provide insights into the relative domain (LID and CORE) displacements in presence of the crowders. Solvation studies reveal that the protein matrix surrounding the probe CPM (7-diethylamino-3-(4-maleimido-phenyl)-4-methylcoumarin) gets restricted in presence of the crowders, with Ficoll 70 providing the maximum rigidity, the same being linked to the decrease in the activity of the enzyme. Through our multipronged approach, we have observed a distinct correlation between domain displacement, enzyme activity and associated dynamics. Thus, keeping in mind the complex nature of enzyme activity and the surrounding bath of dense soup that the biological entity remains immersed in, indeed more such approaches need to be undertaken to have a better grasp of the "enzymes in the crowd".

Coarse-Grained Modeling of Multiple Pathways in Conformational Transitions of Multi-Domain Proteins

Large-scale conformational transitions in multi-domain proteins are often essential for their functions. To investigate the transitions, it is necessary to explore multiple potential pathways, which involve different intermediate structures. Here, we present a multi-basin (MB) coarse-grained (CG) structure-based Go̅ model for describing transitions in proteins with more than two moving domains. This model is an extension of our dual-basin Go̅ model in which system-dependent parameters are determined systematically using the multistate Bennett acceptance ratio method. In the MB Go̅ model for multi-domain proteins, we assume that intermediate structures may have partial inter-domain native contacts. This approach allows us to search multiple transition pathways that involve distinct intermediate structures using the CG molecular dynamics (MD) simulations. We apply this scheme to an enzyme, adenylate kinase (AdK), which has three major domains and can move along two different pathways. Using the optimized mixing parameters for each pathway, AdK shows frequent transitions between the Open, Closed, and the intermediate basins and samples a wide variety of conformations within each basin. The explored multiple transition pathways could be compared with experimental data and examined in more detail by atomistic MD simulations.

Qualitative Prediction of Ligand Dissociation Kinetics from Focal Adhesion Kinase Using Steered Molecular Dynamics

Most early-stage drug discovery projects focus on equilibrium binding affinity to the target alongside selectivity and other pharmaceutical properties. Since many approved drugs have nonequilibrium binding characteristics, there has been increasing interest in optimizing binding kinetics early in the drug discovery process. As focal adhesion kinase (FAK) is an important drug target, we examine whether steered molecular dynamics (SMD) can be useful for identifying drug candidates with the desired drug-binding kinetics. In simulating the dissociation of 14 ligands from FAK, we find an empirical power–law relationship between the simulated time needed for ligand unbinding and the experimental rate constant for dissociation, with a strong correlation depending on the SMD force used. To improve predictions, we further develop regression models connecting experimental dissociation rate with various structural and energetic quantities derived from the simulations. These models can be used to predict dissociation rates from FAK for related compounds.

Escape of a Small Molecule from Inside T4 Lysozyme by Multiple Pathways

The T4 lysozyme L99A mutant is often used as a model system to study small-molecule binding to proteins, but pathways for ligand entry and exit from the buried binding site and the associated protein conformational changes have not been fully resolved. Here, molecular dynamics simulations were employed to model benzene exit from its binding cavity using the weighted ensemble (WE) approach to enhance sampling of low-probability unbinding trajectories. Independent WE simulations revealed four pathways for benzene exit, which correspond to transient tunnels spontaneously formed in previous simulations of apo T4 lysozyme. Thus, benzene unbinding occurs through multiple pathways partially created by intrinsic protein structural fluctuations. Motions of several α-helices and side chains were involved in ligand escape from metastable microstates. WE simulations also provided preliminary estimates of rate constants for each exit pathway. These results complement previous works and provide a semiquantitative characterization of pathway heterogeneity for binding of small molecules to proteins.

Multiple Pathways and Time Scales for Conformational Transitions in apo-Adenylate Kinase

The open/close transition in adenylate kinase (AK) is regarded as a representative example for large-scale conformational transition in proteins, yet its mechanism remains unclear despite numerous experimental and computational studies. Using extensive (∼50 μs) explicit solvent atomistic simulations and Markov state analysis, we shed new lights on the mechanism of this transition in the apo form of AK. The closed basin of apo AK features an open NMP domain while the LID domain closes and rotates toward it. Therefore, although the computed structural properties of the closed ensemble are consistent with previously reported FRET and PRE measurements, our simulations suggest that NMP closure is likely to follow AMP binding, in contrast to the previous interpretation of FRET and PRE data that the apo state was able to sample the fully closed conformation for "ligand selection". The closed state ensemble is found to be kinetically heterogeneous; multiple pathways and time scales are associated with the open/close transition, providing new clues to the disparate time scales observed in different experiments. Besides interdomain interactions, a novel mutual information analysis identifies specific intradomain interactions that correlate strongly to transition kinetics, supporting observations from previous chimera experiments. While our results underscore the role of internal domain properties in determining the kinetics of open/close transition in apo AK, no evidence is observed for any significant degree of local unfolding during the transition. These observations about AK have general implications to our view of conformational states, transition pathways, and time scales of conformational changes in proteins. The key features and time scales of observed transition pathways are robust and similar from simulations using two popular fixed charge force fields.

Substrate Binding Specifically Modulates Domain Arrangements in Adenylate Kinase

The enzyme adenylate kinase (ADK) features two substrate binding domains that undergo large-scale motions during catalysis. In the apo state, the enzyme preferentially adopts a globally open state with accessible binding sites. Binding of two substrate molecules (AMP + ATP or ADP + ADP) results in a closed domain conformation, allowing efficient phosphoryl-transfer catalysis. We employed molecular dynamics simulations to systematically investigate how the individual domain motions are modulated by the binding of substrates. Two-dimensional free-energy landscapes were calculated along the opening of the two flexible lid domains for apo and holo ADK as well as for all single natural substrates bound to one of the two binding sites of ADK. The simulations reveal a strong dependence of the conformational ensembles on type and binding position of the bound substrates and a nonsymmetric behavior of the lid domains. Altogether, the ensembles suggest that, upon initial substrate binding to the corresponding lid site, the opposing lid is maintained open and accessible for subsequent substrate binding. In contrast, ATP binding to the AMP-lid induces global domain closing, preventing further substrate binding to the ATP-lid site. This might constitute a mechanism by which the enzyme avoids the formation of a stable but enzymatically unproductive state.

Principles and Overview of Sampling Methods for Modeling Macromolecular Structure and Dynamics

Investigation of macromolecular structure and dynamics is fundamental to understanding how macromolecules carry out their functions in the cell. Significant advances have been made toward this end in silico, with a growing number of computational methods proposed yearly to study and simulate various aspects of macromolecular structure and dynamics. This review aims to provide an overview of recent advances, focusing primarily on methods proposed for exploring the structure space of macromolecules in isolation and in assemblies for the purpose of characterizing equilibrium structure and dynamics. In addition to surveying recent applications that showcase current capabilities of computational methods, this review highlights state-of-the-art algorithmic techniques proposed to overcome challenges posed in silico by the disparate spatial and time scales accessed by dynamic macromolecules. This review is not meant to be exhaustive, as such an endeavor is impossible, but rather aims to balance breadth and depth of strategies for modeling macromolecular structure and dynamics for a broad audience of novices and experts.

A comparison of weighted ensemble and Markov state model methodologies

Computation of reaction rates and elucidation of reaction mechanisms are two of the main goals of molecular dynamics (MD) and related simulation methods. Since it is time consuming to study reaction mechanisms over long time scales using brute force MD simulations, two ensemble methods, Markov State Models (MSMs) and Weighted Ensemble (WE), have been proposed to accelerate the procedure. Both approaches require clustering of microscopic configurations into networks of "macro-states" for different purposes. MSMs model a discretization of the original dynamics on the macro-states. Accuracy of the model significantly relies on the boundaries of macro-states. On the other hand, WE uses macro-states to formulate a resampling procedure that kills and splits MD simulations for achieving better efficiency of sampling. Comparing to MSMs, accuracy of WE rate predictions is less sensitive to the definition of macro-states. Rigorous numerical experiments using alanine dipeptide and penta-alanine support our analyses. It is shown that MSMs introduce significant biases in the computation of reaction rates, which depend on the boundaries of macro-states, and Accelerated Weighted Ensemble (AWE), a formulation of weighted ensemble that uses the notion of colors to compute fluxes, has reliable flux estimation on varying definitions of macro-states. Our results suggest that whereas MSMs provide a good idea of the metastable sets and visualization of overall dynamics, AWE provides reliable rate estimations requiring less efforts on defining macro-states on the high dimensional conformational space.

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

The long-term goal of connecting scales in biological simulation can be facilitated by scale-agnostic methods. We demonstrate that the weighted ensemble (WE) strategy, initially developed for molecular simulations, applies effectively to spatially resolved cell-scale simulations. The WE approach runs an ensemble of parallel trajectories with assigned weights and uses a statistical resampling strategy of replicating and pruning trajectories to focus computational effort on difficult-to-sample regions. The method can also generate unbiased estimates of non-equilibrium and equilibrium observables, sometimes with significantly less aggregate computing time than would be possible using standard parallelization. Here, we use WE to orchestrate particle-based kinetic Monte Carlo simulations, which include spatial geometry (e.g., of organelles, plasma membrane) and biochemical interactions among mobile molecular species. We study a series of models exhibiting spatial, temporal and biochemical complexity and show that although WE has important limitations, it can achieve performance significantly exceeding standard parallel simulation—by orders of magnitude for some observables.

Stochastic simulations (simulations where randomness plays a role) of even simple biological systems are often so computationally intensive that it is impossible, in practice, to simulate them exhaustively and gather good statistics about the likelihood of different outcomes. The difficulty is compounded for the observation of rare events in these simulations; unfortunately, rare events, such as state transitions and barrier crossings, are often those of particular interest. Using the weighted ensemble (WE) method, we are able to enhance the characterization of rare events in cell biology simulations, but in such a way that the statistics for these events remain unbiased. The histogram of outcomes that WE produces has the same shape as a naive one, but the resolution of events in the tails of the histogram is greatly improved. This improved resolution in rare event statistics can be used to infer unbiased estimates of long timescale dynamics from short simulations, and we show that using a weighted ensemble can result in a reduction in total simulation time needed to sample certain events of interest in spatial, stochastic models of biological systems.

Estimating first-passage time distributions from weighted ensemble simulations and non-Markovian analyses

First-passage times (FPTs) are widely used to characterize stochastic processes such as chemical reactions, protein folding, diffusion processes or triggering a stock option. In previous work (Suarez et al., JCTC 2014;10:2658-2667), we demonstrated a non-Markovian analysis approach that, with a sufficient subset of history information, yields unbiased mean first-passage times from weighted-ensemble (WE) simulations. The estimation of the distribution of the first-passage times is, however, a more ambitious goal since it cannot be obtained by direct observation in WE trajectories. Likewise, a large number of events would be required to make a good estimation of the distribution from a regular "brute force" simulation. Here, we show how the previously developed non-Markovian analysis can generate approximate, but highly accurate, FPT distributions from WE data. The analysis can also be applied to any other unbiased trajectories, such as from standard molecular dynamics simulations. The present study employs a range of systems with independent verification of the distributions to demonstrate the success and limitations of the approach. By comparison to a standard Markov analysis, the non-Markovian approach is less sensitive to the user-defined discretization of configuration space.

The Atomistic Mechanism of Conformational Transition of Adenylate Kinase Investigated by Lorentzian Structure-Based Potential

We present a new all-atom structure-based method to study protein conformational transitions using Lorentzian attractive interactions based on native structures. The variability of each native contact is estimated based on evolutionary information using a machine learning method. To test the validity of this approach, we have investigated the conformational transition of adenylate kinase (ADK). The intrinsic boundedness of the Lorentzian attractive interactions facilitated frequent conformational transitions, and consequently we were able to observe more than 1000 structural interconversions between the open and closed states of ADK out of a total of 6 μs MD simulations. ADK has three domains: the nucleoside monophosphate (NMP) binding domain, the LID-domain, and the CORE domain, which catalyze the interconversion between ATP and ADP. We identified two transition states: a more frequent LID-closed-NMP-open (TS1) state and a less frequent LID-open-NMP-closed (TS2) state. The transition was found to be symmetric in both directions via TS1. We also obtained an off-pathway metastable state that was previously observed with physics-based all-atom simulations but not with coarse-grained models. In the metastable state, the LID domain was slightly twisted and formed contacts with the NMP domain. Our model correctly identified a total of 14 out of the top 16 residues with highest fluctuation by NMR experiment, thus showing excellent agreement with experimental NMR relaxation data and overwhelmingly better results than existing models.

Features of large hinge-bending conformational transitions. Prediction of closed structure from open state

We performed a detailed analysis of conformational transition pathways for a set of 10 proteins, which undergo large hinge-bending-type motions with 4-12 Å RMSD (root mean-square distance) between open and closed crystal structures. Anisotropic network model-Monte Carlo (ANM-MC) algorithm generates a targeted pathway between two conformations, where the collective modes from the ANM are used for deformation at each iteration and the conformational energy of the deformed structure is minimized via an MC algorithm. The target structure was approached successfully with an RMSD of 0.9-4.1 Å when a relatively low cutoff radius of 10 Å was used in ANM. Even though one predominant mode (first or second) directed the open-to-closed conformational transition, changes in the dominant mode character were observed for most cases along the transition. By imposing radius of gyration constraint during mode selection, it was possible to predict the closed structure for eight out of 10 proteins (with initial 4.1-7.1 Å and final 1.7-2.9 Å RMSD to target). Deforming along a single mode leads to most successful predictions. Based on the previously reported free energy surface of adenylate kinase, deformations along the first mode produced an energetically favorable path, which was interestingly facilitated by a change in mode shape (resembling second and third modes) at key points. Pathway intermediates are provided in our database of conformational transitions (http://safir.prc.boun.edu.tr/anmmc/method/1).

WESTPA: an interoperable, highly scalable software package for weighted ensemble simulation and analysis

The weighted ensemble (WE) path sampling approach orchestrates an ensemble of parallel calculations with intermittent communication to enhance the sampling of rare events, such as molecular associations or conformational changes in proteins or peptides. Trajectories are replicated and pruned in a way that focuses computational effort on underexplored regions of configuration space while maintaining rigorous kinetics. To enable the simulation of rare events at any scale (e.g., atomistic, cellular), we have developed an open-source, interoperable, and highly scalable software package for the execution and analysis of WE simulations: WESTPA (The Weighted Ensemble Simulation Toolkit with Parallelization and Analysis). WESTPA scales to thousands of CPU cores and includes a suite of analysis tools that have been implemented in a massively parallel fashion. The software has been designed to interface conveniently with any dynamics engine and has already been used with a variety of molecular dynamics (e.g., GROMACS, NAMD, OpenMM, AMBER) and cell-modeling packages (e.g., BioNetGen, MCell). WESTPA has been in production use for over a year, and its utility has been demonstrated for a broad set of problems, ranging from atomically detailed host–guest associations to nonspatial chemical kinetics of cellular signaling networks. The following describes the design and features of WESTPA, including the facilities it provides for running WE simulations and storing and analyzing WE simulation data, as well as examples of input and output.

AWE-WQ: fast-forwarding molecular dynamics using the accelerated weighted ensemble

A limitation of traditional molecular dynamics (MD) is that reaction rates are difficult to compute. This is due to the rarity of observing transitions between metastable states since high energy barriers trap the system in these states. Recently the weighted ensemble (WE) family of methods have emerged which can flexibly and efficiently sample conformational space without being trapped and allow calculation of unbiased rates. However, while WE can sample correctly and efficiently, a scalable implementation applicable to interesting biomolecular systems is not available. We provide here a GPLv2 implementation called AWE-WQ of a WE algorithm using the master/worker distributed computing WorkQueue (WQ) framework. AWE-WQ is scalable to thousands of nodes and supports dynamic allocation of computer resources, heterogeneous resource usage (such as central processing units (CPU) and graphical processing units (GPUs) concurrently), seamless heterogeneous cluster usage (i.e., campus grids and cloud providers), and support for arbitrary MD codes such as GROMACS, while ensuring that all statistics are unbiased. We applied AWE-WQ to a 34 residue protein which simulated 1.5 ms over 8 months with peak aggregate performance of 1000 ns/h. Comparison was done with a 200 μs simulation collected on a GPU over a similar timespan. The folding and unfolded rates were of comparable accuracy.

Conformational Transition Pathways of Epidermal Growth Factor Receptor Kinase Domain from Multiple Molecular Dynamics Simulations and Bayesian Clustering

The epidermal growth factor receptor (EGFR) is aberrantly activated in various cancer cells and an important target for cancer treatment. Deep understanding of EGFR conformational changes between the active and inactive states is of pharmaceutical interest. Here we present a strategy combining multiply targeted molecular dynamics simulations, unbiased molecular dynamics simulations, and Bayesian clustering to investigate transition pathways during the activation/inactivation process of EGFR kinase domain. Two distinct pathways between the active and inactive forms are designed, explored, and compared. Based on Bayesian clustering and rough two-dimensional free energy surfaces, the energy-favorable pathway is recognized, though DFG-flip happens in both pathways. In addition, another pathway with different intermediate states appears in our simulations. Comparison of distinct pathways also indicates that disruption of the Lys745-Glu762 interaction is critically important in DFG-flip while movement of the A-loop significantly facilitates the conformational change. Our simulations yield new insights into EGFR conformational transitions. Moreover, our results verify that this approach is valid and efficient in sampling of protein conformational changes and comparison of distinct pathways.

Global transitions of proteins explored by a multiscale hybrid methodology: application to adenylate kinase

Efficient and accurate mapping of transition pathways is a challenging problem in allosteric proteins. We propose here a to our knowledge new methodology called collective molecular dynamics (coMD). coMD takes advantage of the collective modes of motions encoded by the fold, simultaneously evaluating the interactions and energetics via a full-atomic MD simulation protocol. The basic approach is to deform the structure collectively along the modes predicted by the anisotropic network model, upon selecting them via a Monte Carlo/Metropolis algorithm from among the complete pool of all accessible modes. Application to adenylate kinase, an allosteric enzyme composed of three domains, CORE, LID, and NMP, shows that both open-to-closed and closed-to-open transitions are readily sampled by coMD, with large-scale motions of the LID dominating. An energy-barrier crossing occurs during the NMP movements. The energy barrier originates from a switch between the salt bridges K136-D118 at the LID-CORE interface and K57-E170 and D33-R156 at the CORE-NMP and LID-NMP interfaces, respectively. Despite its simplicity and computing efficiency, coMD yields ensembles of transition pathways in close accord with detailed full atomic simulations, lending support to its utility as a multiscale hybrid method for efficiently exploring the allosteric transitions of multidomain or multimeric proteins.

Simultaneous Computation of Dynamical and Equilibrium Information Using a Weighted Ensemble of Trajectories

Equilibrium formally can be represented as an ensemble of uncoupled systems undergoing unbiased dynamics in which detailed balance is maintained. Many nonequilibrium processes can be described by suitable subsets of the equilibrium ensemble. Here, we employ the “weighted ensemble” (WE) simulation protocol [Huber and Kim, Biophys. J.1996, 70, 97–110] to generate equilibrium trajectory ensembles and extract nonequilibrium subsets for computing kinetic quantities. States do not need to be chosen in advance. The procedure formally allows estimation of kinetic rates between arbitrary states chosen after the simulation, along with their equilibrium populations. We also describe a related history-dependent matrix procedure for estimating equilibrium and nonequilibrium observables when phase space has been divided into arbitrary non-Markovian regions, whether in WE or ordinary simulation. In this proof-of-principle study, these methods are successfully applied and validated on two molecular systems: explicitly solvated methane association and the implicitly solvated Ala4 peptide. We comment on challenges remaining in WE calculations.

Learning Kinetic Distance Metrics for Markov State Models of Protein Conformational Dynamics

Statistical modeling of long timescale dynamics with Markov state models (MSMs) has been shown to be an effective strategy for building quantitative and qualitative insight into protein folding processes. Existing methodologies, however, rely on geometric clustering using distance metrics such as root mean square deviation (RMSD), assuming that geometric similarity provides an adequate basis for the kinetic partitioning of phase space. Here, inspired by advances in the machine learning community, we introduce a new approach for learning a distance metric explicitly constructed to model kinetic similarity. This approach enables the construction of models, especially in the regime of high anisotropy in the diffusion constant, with fewer states than was previously possible. Application of this technique to the analysis of two ultralong molecular dynamics simulations of the FiP35 WW domain identifies discrete near-native relaxation dynamics in the millisecond regime that were not resolved in previous analyses.

Simulating rare events using a weighted ensemble-based string method

We introduce an extension to the weighted ensemble (WE) path sampling method to restrict sampling to a one-dimensional path through a high dimensional phase space. Our method, which is based on the finite-temperature string method, permits efficient sampling of both equilibrium and non-equilibrium systems. Sampling obtained from the WE method guides the adaptive refinement of a Voronoi tessellation of order parameter space, whose generating points, upon convergence, coincide with the principle reaction pathway. We demonstrate the application of this method to several simple, two-dimensional models of driven Brownian motion and to the conformational change of the nitrogen regulatory protein C receiver domain using an elastic network model. The simplicity of the two-dimensional models allows us to directly compare the efficiency of the WE method to conventional brute force simulations and other path sampling algorithms, while the example of protein conformational change demonstrates how the method can be used to efficiently study transitions in the space of many collective variables.

Allosteric activation transitions in enzymes and biomolecular motors: insights from atomistic and coarse-grained simulations

The chemical step in enzymes is usually preceded by a kinetically distinct activation step that involves large-scale conformational transitions. In "simple" enzymes this step corresponds to the closure of the active site; in more complex enzymes, such as biomolecular motors, the activation step is more complex and may involve interactions with other biomolecules. These activation transitions are essential to the function of enzymes and perturbations in the scale and/or rate of these transitions are implicated in various serious human diseases; incorporating key flexibilities into engineered enzymes is also considered a major remaining challenge in rational enzyme design. Therefore it is important to understand the underlying mechanism of these transitions. This is a significant challenge to both experimental and computational studies because of the allosteric and multi-scale nature of such transitions. Using our recent studies of two enzyme systems, myosin and adenylate kinase (AK), we discuss how atomistic and coarse-grained simulations can be used to provide insights into the mechanism of activation transitions in realistic systems. Collectively, the results suggest that although many allosteric transitions can be viewed as domain displacements mediated by flexible hinges, there are additional complexities and various deviations. For example, although our studies do not find any evidence for "cracking" in AK, our results do underline the contribution of intra-domain properties (e.g., dihedral flexibility) to the rate of the transition. The study of mechanochemical coupling in myosin highlights that local changes important to chemistry require stabilization from more extensive structural changes; in this sense, more global structural transitions are needed to activate the chemistry in the active site. These discussions further emphasize the importance of better understanding factors that control the degree of co-operativity for allosteric transitions, again hinting at the intimate connection between protein stability and functional flexibility. Finally, a number of topics of considerable future interest are briefly discussed.

Exploring the Dynamic Functional Landscape of Adenylate Kinase Modulated by Substrates

Adenylate kinase (ADK) has been explored widely, through both experimental and theoretical studies. However, still less is known about how the functional dynamics of ADK is modulated explicitly by its natural substrates. Here, we report a quantitative study of the dynamic energy landscape for ADK responding to the substrate binding by integrating both experimental investigations and theoretical modeling. We make theoretical predictions which are in remarkable agreement with the single molecule experiments on the substrate-bound complex. With our combined models of ADK in its apo form, in the presence of AMP or ATP, and in complex with both substrates, we specifically address the following key questions: (1) Are there intermediate state(s) during their catalytic cycle and if so how many? (2) How many pathways are there along the open-to-closed transitions and what are their corresponding weights? (3) How do substrates influence the pathway weights and the stability of the intermediates? (4) Which lid's motion is rate-limiting along the turnover cycle, the NMP or the LID domain? Our models predict two major parallel stepwise pathways and two on-pathway intermediates which are denoted as IN (NMP domain open while LID domain closed) and IL (LID domain open and NMP domain closed), respectively. Further investigation of temperature effects suggests that the IN pathway is dominant at room temperature, but the IL pathway is dominant at the optimal temperature. This leads us to propose that the IL pathway is more dominant by entropy and IN pathway by enthalpy. Remarkably, our results show that even with maximum concentrations of natural substrates, ADK still fluctuates between multiple functional states, reflecting an intrinsic capability of large-scale conformational fluctuations which may be essential to its biological function. The results based on the dual-ligands model provide the theoretical validation of random bisubstrate biproducts (Bi-Bi) mechanism for the enzymatic reaction of ADK. Additionally, the pathway flux analysis strongly suggests that the motion of the NMP domain is the rate-determining step for the conformational cycle (opening and closing).

An adaptive weighted ensemble procedure for efficient computation of free energies and first passage rates

We introduce an adaptive weighted-ensemble procedure (aWEP) for efficient and accurate evaluation of first-passage rates between states for two-state systems. The basic idea that distinguishes aWEP from conventional weighted-ensemble (WE) methodology is the division of the configuration space into smaller regions and equilibration of the trajectories within each region upon adaptive partitioning of the regions themselves into small grids. The equilibrated conditional∕transition probabilities between each pair of regions lead to the determination of populations of the regions and the first-passage times between regions, which in turn are combined to evaluate the first passage times for the forward and backward transitions between the two states. The application of the procedure to a non-trivial coarse-grained model of a 70-residue calcium binding domain of calmodulin is shown to efficiently yield information on the equilibrium probabilities of the two states as well as their first passage times. Notably, the new procedure is significantly more efficient than the canonical implementation of the WE procedure, and this improvement becomes even more significant at low temperatures.

Tunable, mixed-resolution modeling using library-based Monte Carlo and graphics processing units

Building on our recently introduced library-based Monte Carlo (LBMC) approach, we describe a flexible protocol for mixed coarse-grained (CG)/all-atom (AA) simulation of proteins and ligands. In the present implementation of LBMC, protein side chain configurations are pre-calculated and stored in libraries, while bonded interactions along the backbone are treated explicitly. Because the AA side chain coordinates are maintained at minimal run-time cost, arbitrary sites and interaction terms can be turned on to create mixed-resolution models. For example, an AA region of interest such as a binding site can be coupled to a CG model for the rest of the protein. We have additionally developed a hybrid implementation of the generalized Born/surface area (GBSA) implicit solvent model suitable for mixed-resolution models, which in turn was ported to a graphics processing unit (GPU) for faster calculation. The new software was applied to study two systems: (i) the behavior of spin labels on the B1 domain of protein G (GB1) and (ii) docking of randomly initialized estradiol configurations to the ligand binding domain of the estrogen receptor (ERα). The performance of the GPU version of the code was also benchmarked in a number of additional systems.

Flow-Dependent Unfolding and Refolding of an RNA by Nonequilibrium Umbrella Sampling

Nonequilibrium experiments of single biomolecules such as force-induced unfolding reveal details about a few degrees of freedom of a complex system. Molecular dynamics simulations can provide complementary information, but exploration of the space of possible configurations is often hindered by large barriers in phase space that separate metastable regions. To solve this problem, enhanced sampling methods have been developed that divide a phase space into regions and integrate trajectory segments in each region. These methods boost the probability of passage over barriers and facilitate parallelization since integration of the trajectory segments does not require communication, aside from their initialization and termination. Here, we present a parallel version of an enhanced sampling method suitable for systems driven far from equilibrium: nonequilibrium umbrella sampling (NEUS). We apply this method to a coarse-grained model of a 262-nucleotide RNA molecule that unfolds and refolds in an explicit flow field modeled with stochastic rotation dynamics. Using NEUS, we are able to observe extremely rare unfolding events that have mean first passage times as long as 45 s (1.1 × 10(15) dynamics steps). We examine the unfolding process for a range of flow rates of the medium, and we describe two competing pathways in which different intramolecular contacts are broken.

Beyond microscopic reversibility: Are observable non-equilibrium processes precisely reversible?

Although the principle of microscopic reversibility has been studied for many decades, there remain ambiguities in its application to non-equilibrium processes of importance to chemistry, physics and biology. Examples include whether protein unfolding should follow the same pathways and in the same proportions as folding, and whether unbinding should likewise mirror binding. Using continuum-space calculations which extend previous kinetic analyses, we demonstrate formally that the precise symmetry of forward and reverse processes is expected only under certain special conditions. Approximate symmetry will be exhibited under a separate set of conditions. Exact, approximate, and broken symmetry scenarios are verified in several ways: using numerical calculations on toy and molecular systems; using exact calculations on kinetic models of induced fit in protein-ligand binding; and based on reported experimental results. The analysis highlights intrinsic challenges and ambiguities in the design and analysis of both experiments and simulations.