Simulation of an Ensemble of Conformational Transitions in a United-Residue Model of Calmodulin

Experiment and Simulation Reveal Residue Details for How Target Binding Tunes Calmodulin's Calcium-Binding Properties

We aim to elucidate the molecular mechanism of the reciprocal relation of calmodulin's (CaM) target binding and its affinity for calcium ions (Ca²⁺), which is central to decoding CaM-dependent Ca²⁺ signaling in a cell. We employed stopped-flow experiments and coarse-grained molecular simulations that learn the coordination chemistry of Ca²⁺ in CaM from first-principle calculations. The associative memories as part of the coarse-grained force fields built on known protein structures further influence CaM's selection of its polymorphic target peptides in the simulations. We modeled the peptides from the Ca²⁺/CaM-binding domain of Ca²⁺/CaM-dependent kinase II (CaMKII), CaMKIIp (293-310) and selected distinctive mutations at the N-terminus. Our stopped-flow experiments have shown that the CaM's affinity for Ca²⁺ in the bound complex of Ca²⁺/CaM/CaMKIIp decreased significantly when Ca²⁺/CaM bound to the mutant peptide (296-AAA-298) compared to that bound to the wild-type peptide (296-RRK-298). The coarse-grained molecular simulations revealed that the 296-AAA-298 mutant peptide destabilized the structures of Ca²⁺-binding loops at the C-domain of CaM (c-CaM) due to both loss of electrostatic interactions and differences in polymorphic structures. We have leveraged a powerful coarse-grained approach to advance a residue-level understanding of the reciprocal relation in CaM, that could not be possibly achieved by other computational approaches.

Coarse-Grained Modeling and Molecular Dynamics Simulations of Ca2+-Calmodulin

Calmodulin (CaM) is a calcium-binding protein that transduces signals to downstream proteins through target binding upon calcium binding in a time-dependent manner. Understanding the target binding process that tunes CaM’s affinity for the calcium ions (Ca²⁺), or vice versa, may provide insight into how Ca²⁺-CaM selects its target binding proteins. However, modeling of Ca²⁺-CaM in molecular simulations is challenging because of the gross structural changes in its central linker regions while the two lobes are relatively rigid due to tight binding of the Ca²⁺ to the calcium-binding loops where the loop forms a pentagonal bipyramidal coordination geometry with Ca²⁺. This feature that underlies the reciprocal relation between Ca²⁺ binding and target binding of CaM, however, has yet to be considered in the structural modeling. Here, we presented a coarse-grained model based on the Associative memory, Water mediated, Structure, and Energy Model (AWSEM) protein force field, to investigate the salient features of CaM. Particularly, we optimized the force field of CaM and that of Ca²⁺ ions by using its coordination chemistry in the calcium-binding loops to match with experimental observations. We presented a “community model” of CaM that is capable of sampling various conformations of CaM, incorporating various calcium-binding states, and carrying the memory of binding with various targets, which sets the foundation of the reciprocal relation of target binding and Ca²⁺ binding in future studies.

Overcoming the Bottleneck of the Enzymatic Cycle by Steric Frustration

The enormous catalytic power of natural enzymes relies on the ability to overcome the bottleneck event in the enzymatic cycle, yet the underlying physical mechanisms are not fully understood. Here, by performing molecular simulations of the whole enzymatic cycle for a model multisubstrate enzyme with a dynamic energy landscape model, we show that multisubstrate enzymes can utilize steric frustration to facilitate the rate-limiting product-release step. During the enzymatic cycles, the bottleneck product is actively squeezed out by the binding of a new substrate at the neighboring site through the population of a substrate-product cobound complex, in which the binding pockets are frustrated due to steric incompatibility. Such steric frustration thereby enables an active mechanism of product release driven by substrate-binding energy, facilitating the enzymatic cycle.

Middle-way flexible docking: Pose prediction using mixed-resolution Monte Carlo in estrogen receptor α

There is a vast gulf between the two primary strategies for simulating protein-ligand interactions. Docking methods significantly limit or eliminate protein flexibility to gain great speed at the price of uncontrolled inaccuracy, whereas fully flexible atomistic molecular dynamics simulations are expensive and often suffer from limited sampling. We have developed a flexible docking approach geared especially for highly flexible or poorly resolved targets based on mixed-resolution Monte Carlo (MRMC), which is intended to offer a balance among speed, protein flexibility, and sampling power. The binding region of the protein is treated with a standard atomistic force field, while the remainder of the protein is modeled at the residue level with a Gō model that permits protein flexibility while saving computational cost. Implicit solvation is used. Here we assess three facets of the MRMC approach with implications for other docking studies: (i) the role of receptor flexibility in cross-docking pose prediction; (ii) the use of non-equilibrium candidate Monte Carlo (NCMC) and (iii) the use of pose-clustering in scoring. We examine 61 co-crystallized ligands of estrogen receptor α, an important cancer target known for its flexibility. We also compare the performance of the MRMC approach with Autodock smina. Adding protein flexibility, not surprisingly, leads to significantly lower total energies and stronger interactions between protein and ligand, but notably we document the important role of backbone flexibility in the improvement. The improved backbone flexibility also leads to improved performance relative to smina. Somewhat unexpectedly, our implementation of NCMC leads to only modestly improved sampling of ligand poses. Overall, the addition of protein flexibility improves the performance of docking, as measured by energy-ranked poses, but we do not find significant improvements based on cluster information or the use of NCMC. We discuss possible improvements for the model including alternative coarse-grained force fields, improvements to the treatment of solvation, and adding additional types of NCMC moves.

Quantifying the Intrinsic Conformation Energy Landscape Topography of Proteins with Large-Scale Open-Closed Transition

Large-scale conformational changes of proteins, including the open-closed transitions, are crucial for a variety of protein functions. These open-closed transitions are often associated with ligand binding. However, the understandings of the underlying mechanisms of the conformational changes within proteins during the open-closed transitions are still challenging at present. In this study, we quantified the intrinsic underlying conformational energy landscapes of five different proteins with large-scale open-closed transitions. This is realized by exploring the underlying density of states and the intrinsic conformational energy landscape topography measure Λ. Λ is a dimensionless ratio of conformational energy gap δE versus conformational energy roughness δE and configurational entropy S or size of the intrinsic conformational energy landscape. By quantifying the Λ of intrinsic open-closed conformational (Λ^oc) and intrinsic global folding (Λ^global) energy landscapes, we show that both intrinsic open-closed conformation energy and entropy landscapes are funneled toward the closed state. Furthermore, our results indicate the strong correlations between Λ and thermodynamics (conformational state transition temperature against trapping temperature) as well as between Λ and kinetics (open-closed kinetic time) of these proteins. This shows that the intrinsic conformational landscape topography determines both the conformational thermodynamic stability and kinetic speed of the conformational dynamics. Our investigations provide important insights for understanding the fundamental mechanisms of the protein conformational dynamics in a physical and global way.

Fast Protein Translation Can Promote Co- and Posttranslational Folding of Misfolding-Prone Proteins

Chemical kinetic modeling has previously been used to predict that fast-translating codons can enhance cotranslational protein folding by helping to avoid misfolded intermediates. Consistent with this prediction, protein aggregation in yeast and worms was observed to increase when translation was globally slowed down, possibly due to increased cotranslational misfolding. Observation of similar behavior in molecular simulations would confirm predictions from the simpler chemical kinetic model and provide a molecular perspective on cotranslational folding, misfolding, and the impact of translation speed on these processes. All-atom simulations cannot reach the timescales relevant to protein synthesis, and most conventional structure-based coarse-grained models do not allow for nonnative structure formation. Here, we introduce a protocol to incorporate misfolding using the functional forms of publicly available force fields. With this model we create two artificial proteins that are capable of undergoing structural transitions between a native and a misfolded conformation and simulate their synthesis by the ribosome. Consistent with the chemical kinetic predictions, we find that rapid synthesis of misfolding-prone nascent-chain segments increases the fraction of folded proteins by kinetically partitioning more molecules through on-pathway intermediates, decreasing the likelihood of sampling misfolded conformations. Novel to this study, to our knowledge, we observe that differences in protein dynamics, arising from different translation-elongation schedules, can persist long after the nascent protein has been released from the ribosome, and that a sufficient level of energetic frustration is needed for fast-translating codons to be beneficial for folding. These results provide further evidence that fast-translating codons can be as biologically important as pause sites in coordinating cotranslational folding.

Tabulation as a high-resolution alternative to coarse-graining protein interactions: Initial application to virus capsid subunits

Traditional coarse-graining based on a reduced number of interaction sites often entails a significant sacrifice of chemical accuracy. As an alternative, we present a method for simulating large systems composed of interacting macromolecules using an energy tabulation strategy previously devised for small rigid molecules or molecular fragments [S. Lettieri and D. M. Zuckerman, J. Comput. Chem. 33, 268-275 (2012); J. Spiriti and D. M. Zuckerman, J. Chem. Theory Comput. 10, 5161-5177 (2014)]. We treat proteins as rigid and construct distance and orientation-dependent tables of the interaction energy between them. Arbitrarily detailed interactions may be incorporated into the tables, but as a proof-of-principle, we tabulate a simple α-carbon Gō-like model for interactions between dimeric subunits of the hepatitis B viral capsid. This model is significantly more structurally realistic than previous models used in capsid assembly studies. We are able to increase the speed of Monte Carlo simulations by a factor of up to 6700 compared to simulations without tables, with only minimal further loss in accuracy. To obtain further enhancement of sampling, we combine tabulation with the weighted ensemble (WE) method, in which multiple parallel simulations are occasionally replicated or pruned in order to sample targeted regions of a reaction coordinate space. In the initial study reported here, WE is able to yield pathways of the final ∼25% of the assembly process.

Demonstrated Convergence of the Equilibrium Ensemble for a Fast United-Residue Protein Model

Because of the time-scale limitations of all-atom simulation of proteins, there has been substantial interest in coarse-grained approaches. Some methods, like "resolution exchange" (Lyman, E.; Ytreberg, F. M.; Zuckerman, D. M. Phys. Rev. Lett. 2006, 96, 028105-1-4), can accelerate canonical all-atom sampling but require properly distributed coarse ensembles. We therefore demonstrate that full sampling can indeed be achieved in a sufficiently simplified protein model, as verified by a recently developed convergence analysis. The model accounts for protein backbone geometry, in that rigid peptide planes rotate according to atomistically defined dihedral angles, but there are only two degrees of freedom (φ and ψ dihedrals) per residue. Our convergence analysis indicates that small proteins (up to 89 residues in our tests) can be simulated for more than 50 "structural decorrelation times" in less than a week on a single processor. We show that the fluctuation behavior is reasonable, and we discuss applications, limitations, and extensions of the model.

Domain Motion Enhanced (DoME) Model for Efficient Conformational Sampling of Multidomain Proteins

Large conformational changes of multidomain proteins are difficult to simulate using all-atom molecular dynamics (MD) due to the slow time scale. We show that a simple modification of the structure-based coarse-grained (CG) model enables a stable and efficient MD simulation of those proteins. "Motion Tree", a tree diagram that describes conformational changes between two structures in a protein, provides information on rigid structural units (domains) and the magnitudes of domain motions. In our new CG model, which we call the DoME (domain motion enhanced) model, interdomain interactions are defined as being inversely proportional to the magnitude of the domain motions in the diagram, whereas intradomain interactions are kept constant. We applied the DoME model in combination with the Go model to simulations of adenylate kinase (AdK). The results of the DoME-Go simulation are consistent with an all-atom MD simulation for 10 μs as well as known experimental data. Unlike the conventional Go model, the DoME-Go model yields stable simulation trajectories against temperature changes and conformational transitions are easily sampled despite domain rigidity. Evidently, identification of domains and their interfaces is useful approach for CG modeling of multidomain proteins.

Tunable Coarse Graining for Monte Carlo Simulations of Proteins via Smoothed Energy Tables: Direct and Exchange Simulations

Many commonly used coarse-grained models for proteins are based on simplified interaction sites and consequently may suffer from significant limitations, such as the inability to properly model protein secondary structure without the addition of restraints. Recent work on a benzene fluid (Lettieri S.; Zuckerman D. M.J. Comput. Chem.2012, 33, 268-275) suggested an alternative strategy of tabulating and smoothing fully atomistic orientation-dependent interactions among rigid molecules or fragments. Here we report our initial efforts to apply this approach to the polar and covalent interactions intrinsic to polypeptides. We divide proteins into nearly rigid fragments, construct distance and orientation-dependent tables of the atomistic interaction energies between those fragments, and apply potential energy smoothing techniques to those tables. The amount of smoothing can be adjusted to give coarse-grained models that range from the underlying atomistic force field all the way to a bead-like coarse-grained model. For a moderate amount of smoothing, the method is able to preserve about 70-90% of the α-helical structure while providing a factor of 3-10 improvement in sampling per unit computation time (depending on how sampling is measured). For a greater amount of smoothing, multiple folding-unfolding transitions of the peptide were observed, along with a factor of 10-100 improvement in sampling per unit computation time, although the time spent in the unfolded state was increased compared with less smoothed simulations. For a β hairpin, secondary structure is also preserved, albeit for a narrower range of the smoothing parameter and, consequently, for a more modest improvement in sampling. We have also applied the new method in a "resolution exchange" setting, in which each replica runs a Monte Carlo simulation with a different degree of smoothing. We obtain exchange rates that compare favorably to our previous efforts at resolution exchange (Lyman E.; Zuckerman D. M.J. Chem. Theory Comput.2006, 2, 656-666).

Aggregation and network formation in self-assembly of protein (H3.1) by a coarse-grained Monte Carlo simulation

Multi-scale aggregation to network formation of interacting proteins (H3.1) are examined by a knowledge-based coarse-grained Monte Carlo simulation as a function of temperature and the number of protein chains, i.e., the concentration of the protein. Self-assembly of corresponding homo-polymers of constitutive residues (Cys, Thr, and Glu) with extreme residue-residue interactions, i.e., attractive (Cys-Cys), neutral (Thr-Thr), and repulsive (Glu-Glu), are also studied for comparison with the native protein. Visual inspections show contrast and similarity in morphological evolutions of protein assembly, aggregation of small aggregates to a ramified network from low to high temperature with the aggregation of a Cys-polymer, and an entangled network of Glu and Thr polymers. Variations in mobility profiles of residues with the concentration of the protein suggest that the segmental characteristic of proteins is altered considerably by the self-assembly from that in its isolated state. The global motion of proteins and Cys polymer chains is enhanced by their interacting network at the low temperature where isolated chains remain quasi-static. Transition from globular to random coil transition, evidenced by the sharp variation in the radius of gyration, of an isolated protein is smeared due to self-assembly of interacting networks of many proteins. Scaling of the structure factor S(q) with the wave vector q provides estimates of effective dimension D of the mass distribution at multiple length scales in self-assembly. Crossover from solid aggregates (D ∼ 3) at low temperature to a ramified fibrous network (D ∼ 2) at high temperature is observed for the protein H3.1 and Cys polymers in contrast to little changes in mass distribution (D ∼ 1.6) of fibrous Glu- and Thr-chain configurations.

Energy landscape views for interplays among folding, binding, and allostery of calmodulin domains

Ligand binding modulates the energy landscape of proteins, thus altering their folding and allosteric conformational dynamics. To investigate such interplay, calmodulin has been a model protein. Despite much attention, fully resolved mechanisms of calmodulin folding/binding have not been elucidated. Here, by constructing a computational model that can integrate folding, binding, and allosteric motions, we studied in-depth folding of isolated calmodulin domains coupled with binding of two calcium ions and associated allosteric conformational changes. First, mechanically pulled simulations revealed coexistence of three distinct conformational states: the unfolded, the closed, and the open states, which is in accord with and augments structural understanding of recent single-molecule experiments. Second, near the denaturation temperature, we found the same three conformational states as well as three distinct binding states: zero, one, and two calcium ion bound states, leading to as many as nine states. Third, in terms of the nine-state representation, we found multiroute folding/binding pathways and shifts in their probabilities with the calcium concentration. At a lower calcium concentration, "combined spontaneous folding and induced fit" occurs, whereas at a higher concentration, "binding-induced folding" dominates. Even without calcium binding, we observed that the folding pathway of calmodulin domains can be modulated by the presence of metastable states. Finally, full-length calmodulin also exhibited an intriguing coupling between two domains when applying tension.

Discriminating binding mechanisms of an intrinsically disordered protein via a multi-state coarse-grained model

Many proteins undergo a conformational transition upon binding to their cognate binding partner, with intrinsically disordered proteins (IDPs) providing an extreme example in which a folding transition occurs. However, it is often not clear whether this occurs via an "induced fit" or "conformational selection" mechanism, or via some intermediate scenario. In the first case, transient encounters with the binding partner favour transitions to the bound structure before the two proteins dissociate, while in the second the bound structure must be selected from a subset of unbound structures which are in the correct state for binding, because transient encounters of the incorrect conformation with the binding partner are most likely to result in dissociation. A particularly interesting situation involves those intrinsically disordered proteins which can bind to different binding partners in different conformations. We have devised a multi-state coarse-grained simulation model which is able to capture the binding of IDPs in alternate conformations, and by applying it to the binding of nuclear coactivator binding domain (NCBD) to either ACTR or IRF-3 we are able to determine the binding mechanism. By all measures, the binding of NCBD to either binding partner appears to occur via an induced fit mechanism. Nonetheless, we also show how a scenario closer to conformational selection could arise by choosing an alternative non-binding structure for NCBD.

Perspective: Coarse-grained models for biomolecular systems

By focusing on essential features, while averaging over less important details, coarse-grained (CG) models provide significant computational and conceptual advantages with respect to more detailed models. Consequently, despite dramatic advances in computational methodologies and resources, CG models enjoy surging popularity and are becoming increasingly equal partners to atomically detailed models. This perspective surveys the rapidly developing landscape of CG models for biomolecular systems. In particular, this review seeks to provide a balanced, coherent, and unified presentation of several distinct approaches for developing CG models, including top-down, network-based, native-centric, knowledge-based, and bottom-up modeling strategies. The review summarizes their basic philosophies, theoretical foundations, typical applications, and recent developments. Additionally, the review identifies fundamental inter-relationships among the diverse approaches and discusses outstanding challenges in the field. When carefully applied and assessed, current CG models provide highly efficient means for investigating the biological consequences of basic physicochemical principles. Moreover, rigorous bottom-up approaches hold great promise for further improving the accuracy and scope of CG models for biomolecular systems.

Protein recognition and selection through conformational and mutually induced fit

Protein-protein interactions drive most every biological process, but in many instances the domains mediating recognition are disordered. How specificity in binding is attained in the absence of defined structure contrasts with well-established experimental and theoretical work describing ligand binding to protein. The signaling protein calmodulin presents a unique opportunity to investigate mechanisms for target recognition given that it interacts with several hundred different targets. By advancing coarse-grained computer simulations and experimental techniques, mechanistic insights were gained in defining the pathways leading to recognition and in how target selectivity can be achieved at the molecular level. A model requiring mutually induced conformational changes in both calmodulin and target proteins was necessary and broadly informs how proteins can achieve both high affinity and high specificity.

Entropic mechanism of allosteric communication in conformational transitions of dihydrofolate reductase

The mechanism of allosteric conformational transitions of Escherichia coli dihydrofolate reductase (DHFR) is investigated theoretically by applying a newly developed coarse-grained model. Functional forms of interaction potentials in the model depend on the local structural environments around those interactions to represent the many-residue effects due to atomic packing in each local region, and hence, this model is called "the chameleon model". The chameleon model consistently describes the free-energy landscape of two conformational transitions in the catalytic cycle of DHFR, which we call conformational transition 1 (CT1) and conformational transition 2 (CT2); CT1 is accompanied by the hydride transfer reaction, and CT2 is accompanied by the product ligand release. The transition state of CT1 is entropically stabilized by the disordering of loops at the peripheral regions of the protein, which enhances the positively correlated fluctuations at the center part of the protein, showing that the allosteric communication between distant regions through the central region is intrinsically associated with the entropic stabilization of the transition state. The transition state of CT2 is entropically stabilized through the mechanism that enhances the breathing motion of two domains, showing that the difference in the distribution of interactions brings about the difference in the transition mechanism between CT1 and CT2. The chameleon model opens a way to consistently describe the dynamical energy landscape of enzymatic reactions.

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.

Biomolecular dynamics: order-disorder transitions and energy landscapes

While the energy landscape theory of protein folding is now a widely accepted view for understanding how relatively weak molecular interactions lead to rapid and cooperative protein folding, such a framework must be extended to describe the large-scale functional motions observed in molecular machines. In this review, we discuss (1) the development of the energy landscape theory of biomolecular folding, (2) recent advances toward establishing a consistent understanding of folding and function and (3) emerging themes in the functional motions of enzymes, biomolecular motors and other biomolecular machines. Recent theoretical, computational and experimental lines of investigation have provided a very dynamic picture of biomolecular motion. In contrast to earlier ideas, where molecular machines were thought to function similarly to macroscopic machines, with rigid components that move along a few degrees of freedom in a deterministic fashion, biomolecular complexes are only marginally stable. Since the stabilizing contribution of each atomic interaction is on the order of the thermal fluctuations in solution, the rigid body description of molecular function must be revisited. An emerging theme is that functional motions encompass order-disorder transitions and structural flexibility provides significant contributions to the free energy. In this review, we describe the biological importance of order-disorder transitions and discuss the statistical-mechanical foundation of theoretical approaches that can characterize such transitions.

Elastic Network Models are Robust to Variations in Formalism

Understanding the functions of biomolecules requires insight not only from structures, but from dynamics as well. Often, the most interesting processes occur on time scales too slow for exploration by conventional molecular dynamics (MD) simulations. For this reason, alternative computational methods such as elastic network models (ENMs) have become increasingly popular. These simple, coarse-grained models represent molecules as beads connected by harmonic springs; the system's motions are solved analytically by normal mode analysis. In the past few years, many different formalisms for performing ENM calculations have emerged, and several have been optimized using all-atom MD simulations. In contrast to other studies, we have compared the various formalisms in a systematic, quantitative way. In this study, we optimize many ENM functional forms using a uniform dataset containing only long (> 1 μs) all-atom MD simulations. Our results show that all models once optimized produce spring constants for immediate neighboring residues that are orders of magnitude stiffer than more distal contacts. In addition, the statistical significance of ENM performance varied with model resolution. We also show that fitting long trajectories does not improve ENM performance due to a problem inherent in all network models tested: they underestimate the relative importance of the most concerted motions. Finally, we characterize ENMs' resilience by tessellating the parameter space to show that broad ranges of parameters produce similar quality predictions. Taken together our data reveals that choice of spring function and parameters are not vital to performance of a network model and that simple parameters can by derived "by hand" when no data is available for fitting, thus illustrating the robustness of these models.

Exploration of multi-state conformational dynamics and underlying global functional landscape of maltose binding protein

An increasing number of biological machines have been revealed to have more than two macroscopic states. Quantifying the underlying multiple-basin functional landscape is essential for understanding their functions. However, the present models seem to be insufficient to describe such multiple-state systems. To meet this challenge, we have developed a coarse grained triple-basin structure-based model with implicit ligand. Based on our model, the constructed functional landscape is sufficiently sampled by the brute-force molecular dynamics simulation. We explored maltose-binding protein (MBP) which undergoes large-scale domain motion between open, apo-closed (partially closed) and holo-closed (fully closed) states responding to ligand binding. We revealed an underlying mechanism whereby major induced fit and minor population shift pathways co-exist by quantitative flux analysis. We found that the hinge regions play an important role in the functional dynamics as well as that increases in its flexibility promote population shifts. This finding provides a theoretical explanation of the mechanistic discrepancies in PBP protein family. We also found a functional “backtracking” behavior that favors conformational change. We further explored the underlying folding landscape in response to ligand binding. Consistent with earlier experimental findings, the presence of ligand increases the cooperativity and stability of MBP. This work provides the first study to explore the folding dynamics and functional dynamics under the same theoretical framework using our triple-basin functional model.

A central goal of biology is to understand the function of the organism and its constituent parts at each of its scales of complexity. Function at the molecular level is often realized by changes in conformation. Unfortunately, experimental explorations of global motions critical for functional conformational changes are still challenging. In the present work, we developed a coarse grained triple-well structure-based model to explore the underlying functional landscape of maltose-binding protein (MBP). By quantitative flux analysis, we uncover the underlying mechanism by which the major induced fit and minor population shift pathways co-exist. Though we have previously lent credence to the assertion that dynamical equilibrium between open and minor closed conformations exist for all the free PBPs, the generality of this rule is still a matter of open debate. We found that the hinge flexibility is favorable to population shift mechanism. This finding provides a theoretical explanation of the mechanism discrepancies in PBP protein family. We also simulated the folding dynamics using this functional multi-basin model which successfully reproduced earlier protein melting experiment. This represents an exciting opportunity to characterize the interplay between folding and function, which is a long-standing question in the community. The theoretical approach employed in this study is general and can be applied to other systems.

αC helix as a switch in the conformational transition of Src/CDK-like kinase domains

One mechanism of regulating the catalytic activity of protein kinases is through conformational transitions. Despite great diversity in the structural changes involved in the transitions, a certain set of changes within the kinase domain (KD) has been observed for many kinases including Src and CDK2. We investigated this conformational transition computationally to identify the topological features that are energetically critical to the transition. Results from both molecular dynamics sampling and transition path optimization highlight the displacement of the αC helix as the major energy barrier, mediating the switch of the KD between the active and down-regulated states. The critical role of the αC helix is noteworthy by providing a rationale for a number of activation and deactivation mechanisms known to occur in cells. We find that kinases with the αC helix displacement exist throughout the kinome, suggesting that this feature may have emerged early in evolution.

Massive conformation change in the prion protein: Using dual-basin structure-based models to find misfolding pathways

We employ all-atom structure-based models with a force field with multiple energetic basins for the C-terminal (residues 166-226) of the mammalian prion protein. One basin represents the known alpha-helical (αH) structure while the other represents the same residues in a left-handed beta-helical (LHBH) conformation. The LHBH structure has been proposed to help describe one class of in vitro grown fibrils, as well as possibly self-templating the conversion of normal cellular prion protein to the infectious form. Yet, it is unclear how the protein may make this global rearrangement. Our results demonstrate that the conformation changes are not strongly limited by large-scale geometry modification and that there may exist an overall preference for the LHBH conformation. Furthermore, our model presents novel intermediate trapping conformations with twisted LHBH structure.

The two-pathway model of the biological catch-bond as a limit of the allosteric model

Catch-binding is a counterintuitive phenomenon in which the lifetime of a receptor/ligand bond increases when a force is applied to break the bond. Several mechanisms have been proposed to rationalize catch-binding. In the two-pathway model, the force drives the system away from its native dissociation pathway into an alternative pathway involving a higher energy barrier. Here, we analyze an allosteric model suggesting that a force applied to the complex alters the distribution of receptor conformations, and as a result, induces changes in the ligand-binding site. The model assumes explicitly that the allosteric transitions govern the properties of the ligand site. We demonstrate that the dynamics of the ligand is described by two relaxation times, one of which arises from the allosteric site. Therefore, we argue that one can characterize the allosteric transitions by studying the receptor/ligand binding. We show that the allosteric description reduces to the two-pathway model in the limit when the allosteric transitions are faster than the bond dissociation. The formal results are illustrated with two systems, P-selectin/PSGL-1 and FimH/mannose, subjected to both constant and time-dependent forces. The report advances our understanding of catch-binding by combining alternative physical models into a unified description and makes the problem more tractable for the bond mechanics community.

Conformational flexibility and the mechanisms of allosteric transitions in topologically similar proteins

Conformational flexibility plays a central role in allosteric transition of proteins. In this paper, we extend the analysis of our previous study [S. Tripathi and J. J. Portman, Proc. Natl. Acad. Sci. U.S.A. 106, 2104 (2009)] to investigate how relatively minor structural changes of the meta-stable states can significantly influence the conformational flexibility and allosteric transition mechanism. We use the allosteric transitions of the domains of calmodulin as an example system to highlight the relationship between the transition mechanism and the inter-residue contacts present in the meta-stable states. In particular, we focus on the origin of transient local unfolding (cracking), a mechanism that can lower free energy barriers of allosteric transitions, in terms of the inter-residue contacts of the meta-stable states and the pattern of local strain that develops during the transition. We find that the magnitude of the local strain in the protein is not the sole factor determining whether a region will ultimately crack during the transition. These results emphasize that the residue interactions found exclusively in one of the two meta-stable states is the key in understanding the mechanism of allosteric conformational change.

Accelerating molecular Monte Carlo simulations using distance and orientation-dependent energy tables: tuning from atomistic accuracy to smoothed "coarse-grained" models

Typically, the most time consuming part of any atomistic molecular simulation is the repeated calculation of distances, energies, and forces between pairs of atoms. However, many molecules contain nearly rigid multi-atom groups such as rings and other conjugated moieties, whose rigidity can be exploited to significantly speed-up computations. The availability of GB-scale random-access memory (RAM) offers the possibility of tabulation (precalculation) of distance- and orientation-dependent interactions among such rigid molecular bodies. Here, we perform an investigation of this energy tabulation approach for a fluid of atomistic-but rigid-benzene molecules at standard temperature and density. In particular, using O(1) GB of RAM, we construct an energy look-up table, which encompasses the full range of allowed relative positions and orientations between a pair of whole molecules. We obtain a hardware-dependent speed-up of a factor of 24-50 as compared to an ordinary ("exact") Monte Carlo simulation and find excellent agreement between energetic and structural properties. Second, we examine the somewhat reduced fidelity of results obtained using energy tables based on much less memory use. Third, the energy table serves as a convenient platform to explore potential energy smoothing techniques, akin to coarse-graining. Simulations with smoothed tables exhibit near atomistic accuracy while increasing diffusivity. The combined speed-up in sampling from tabulation and smoothing exceeds a factor of 100. For future applications, greater speed-ups can be expected for larger rigid groups, such as those found in biomolecules.

Simulations of the alternating access mechanism of the sodium symporter Mhp1

Sodium coupled cotransporters of the five-helix inverted repeat (5HIR) superfamily use an alternating access mechanism to transport a myriad of small molecules across the cell membrane. One of the primary steps in this mechanism is the conformational transition from a state poised to bind extracellular substrates to a state that is competent to deliver substrate to the cytoplasm. Here, we construct a coarse-grained model of the 5HIR benzylhydantoin transporter Mhp1 that incorporates experimental structures of the outward- and inward-open states to investigate the mechanism of this conformational change. Using the weighted ensemble path-sampling method, we rigorously sample the outward- to inward-facing transition path ensemble. The transition path ensemble reveals a heterogeneous set of pathways connecting the two states and identifies two modes of transport: one consistent with a strict alternating access mechanism and another where decoupling of the inner and outer gates causes the transient formation of a continuous permeation pathway through the transporter. We also show that the conformational switch between the outward- and inward-open states results from rigid body motions of the hash motif relative to the substrate bundle, supporting the rocking bundle hypothesis. Finally, our methodology provides the groundwork for more chemically detailed investigations of the alternating mechanism.

The effect of macromolecular crowding, ionic strength and calcium binding on calmodulin dynamics

The flexibility in the structure of calmodulin (CaM) allows its binding to over 300 target proteins in the cell. To investigate the structure-function relationship of CaM, we combined methods of computer simulation and experiments based on circular dichroism (CD) to investigate the structural characteristics of CaM that influence its target recognition in crowded cell-like conditions. We developed a unique multiscale solution of charges computed from quantum chemistry, together with protein reconstruction, coarse-grained molecular simulations, and statistical physics, to represent the charge distribution in the transition from apoCaM to holoCaM upon calcium binding. Computationally, we found that increased levels of macromolecular crowding, in addition to calcium binding and ionic strength typical of that found inside cells, can impact the conformation, helicity and the EF hand orientation of CaM. Because EF hand orientation impacts the affinity of calcium binding and the specificity of CaM's target selection, our results may provide unique insight into understanding the promiscuous behavior of calmodulin in target selection inside cells.

Proteins are workhorses for driving biological functions inside cells. Calmodulin (CaM) is a protein that can carry cellular signals by triggered conformational changes due to calcium binding that alters target binding. Interestingly, CaM is able to bind over 300 targets. One of the challenges in characterizing CaM's ability to bind multiple targets lies in that CaM is a flexible protein and its structure is easily modulated by the physicochemical changes in its surroundings, particularly inside a complex cellular milieu. In order to determine structure-function relationships of CaM, we employed a combined approach of experiments, computer simulations and statistical physics in the investigation of the effect of calcium-binding, salt concentration, and macromolecular crowding on CaM. The results revealed unique folding energy landscapes of CaM in the absence and presence of calcium ions and the structural implications of CaM are interpreted under cell-like conditions. Further, a large conformational change in CaM in response to environmental impacts, dictates the packing of local helices that may be critical to its function of target binding and recognition among vast target selections.

FRET-FCS detection of intralobe dynamics in calmodulin

Fluorescence correlation spectroscopy (FCS) can be coupled with Förster resonance energy transfer (FRET) to detect intramolecular dynamics of proteins on the microsecond time scale. Here we describe application of FRET-FCS to detect fluctuations within the N-terminal and C-terminal domains of the Ca(2+)-signaling protein calmodulin. Intramolecular fluctuations were resolved by global fitting of the two fluorescence autocorrelation functions (green-green and red-red) together with the two cross-correlation functions (green-red and red-green). To match the Förster radius for FRET to the dimensions of the N-terminal and C-terminal domains, a near-infrared acceptor fluorophore (Atto 740) was coupled with a green-emitting donor (Alexa Fluor 488). Fluctuations were detected in both domains on the time scale of 30 to 40 μs. In the N-terminal domain, the amplitude of the fluctuations was dependent on occupancy of Ca(2+) binding sites. A high amplitude of dynamics in apo-calmodulin (in the absence of Ca(2+)) was nearly abolished at a high Ca(2+) concentration. For the C-terminal domain, the dynamic amplitude changed little with Ca(2+) concentration. The Ca(2+) dependence of dynamics for the N-terminal domain suggests that the fluctuations detected by FCS in the N-terminal domain are coupled to the opening and closing of the EF-hand Ca(2+)-binding loops.

Equilibrium sampling in biomolecular simulations

Equilibrium sampling of biomolecules remains an unmet challenge after more than 30 years of atomistic simulation. Efforts to enhance sampling capability, which are reviewed here, range from the development of new algorithms to parallelization to novel uses of hardware. Special focus is placed on classifying algorithms--most of which are underpinned by a few key ideas--in order to understand their fundamental strengths and limitations. Although algorithms have proliferated, progress resulting from novel hardware use appears to be more clear-cut than from algorithms alone, due partly to the lack of widely used sampling measures.

The origin of nonmonotonic complex behavior and the effects of nonnative interactions on the diffusive properties of protein folding

We present a method for calculating the configurational-dependent diffusion coefficient of a globular protein as a function of the global folding process. Using a coarse-grained structure-based model, we determined the diffusion coefficient, in reaction coordinate space, as a function of the fraction of native contacts formed Q for the cold shock protein (TmCSP). We find nonmonotonic behavior for the diffusion coefficient, with high values for the folded and unfolded ensembles and a lower range of values in the transition state ensemble. We also characterized the folding landscape associated with an energetically frustrated variant of the model. We find that a low-level of frustration can actually stabilize the native ensemble and increase the associated diffusion coefficient. These findings can be understood from a mechanistic standpoint, in that the transition state ensemble has a more homogeneous structural content when frustration is present. Additionally, these findings are consistent with earlier calculations based on lattice models of protein folding and more recent single-molecule fluorescence measurements.

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

We performed "weighted ensemble" path-sampling simulations of adenylate kinase, using several semi-atomistic protein models. The models have an all-atom backbone with various levels of residue interactions. The primary result is that full statistically rigorous path sampling required only a few weeks of single-processor computing time with these models, indicating the addition of further chemical detail should be readily feasible. Our semi-atomistic path ensembles are consistent with previous biophysical findings: the presence of two distinct pathways, identification of intermediates, and symmetry of forward and reverse pathways.

Regulation of catch binding by allosteric transitions

An allosteric model is used to describe changes in lifetimes of biological receptor-ligand bonds subjected to an external force. Force-induced transitions between the two states of the allosteric site lead to changes in the receptor conformation. The ligand bound to the receptor fluctuates between two different potentials formed by the two receptor conformations. The effect of the force on the receptor-ligand interaction potential is described by the Bell mechanism. The probability of detecting the ligand in the bound state is found to depend on the relaxation times of both ligand and allosteric sites. An analytic expression for the bond lifetime is derived as a function of force. The formal theoretical results are used to explain the anomalous force and time dependences of the integrin-fibronectin bond lifetimes measured by atomic force microscopy (Kong, F.; et al J. Cell Biol. 2009, 185, 1275-1284). The analytic expression and model parameters describe very well all anomalous dependences identified in the experiments.

Automated sampling assessment for molecular simulations using the effective sample size

To quantify the progress in the development of algorithms and forcefields used in molecular simulations, a general method for the assessment of the sampling quality is needed. Statistical mechanics principles suggest the populations of physical states characterize equilibrium sampling in a fundamental way. We therefore develop an approach for analyzing the variances in state populations, which quantifies the degree of sampling in terms of the effective sample size (ESS). The ESS estimates the number of statistically independent configurations contained in a simulated ensemble. The method is applicable to both traditional dynamics simulations as well as more modern (e.g., multi-canonical) approaches. Our procedure is tested in a variety of systems from toy models to atomistic protein simulations. We also introduce a simple automated procedure to obtain approximate physical states from dynamic trajectories: this allows sample-size estimation in systems for which physical states are not known in advance.

Lymphotactin: how a protein can adopt two folds

Metamorphic proteins such as lymphotactin are a notable exception of the empirical principle that structured natural proteins possess a unique three-dimensional structure. In particular, the human chemokine lymphotactin protein exists in two distinct conformations (one monomeric and one dimeric) under physiological conditions. In this work, we use a C(alpha) Go model to show how this very peculiar behavior can be reproduced. From the study of the thermodynamics and of the kinetics, we characterize the interconversion mechanism. In particular, this takes place through the docking of the two chains living in a third monomeric, partially unfolded, state which shows a residual structure involving a set of local contacts common to the two native conformations. The main feature of two fold proteins appears to be the sharing of a common set of local contacts between the two distinct folds as confirmed by the study of two designed two fold proteins. Metamorphic proteins may be more common than expected.

Gating transition of pentameric ligand-gated ion channels

Pentameric ligand-gated ion channels are an important family of membrane proteins and play key roles in physiological processes, including signal transduction at chemical synapses. Here, we study the conformational changes associated with the opening and closing of the channel pore. Based on recent crystal structures of two prokaryotic members of the family in open and closed states, respectively, mixed elastic network models are constructed for the transmembrane domain. To explore the conformational changes in the gating transition, a coarse-grained transition path is computed that smoothly connects the closed and open conformations of the channel. We find that the conformational transition involves no major rotations of the transmembrane helices, and is instead characterized by a concerted tilting of helices M2 and M3. In addition, helix M2 changes its bending state, which results in an early closure of the pore during the open-to-closed transition.

General library-based Monte Carlo technique enables equilibrium sampling of semi-atomistic protein models

We introduce "library-based Monte Carlo" (LBMC) simulation, which performs Boltzmann sampling of molecular systems based on precalculated statistical libraries of molecular-fragment configurations, energies, and interactions. The library for each fragment can be Boltzmann distributed and thus account for all correlations internal to the fragment. LBMC can be applied to both atomistic and coarse-grained models, as we demonstrate in this "proof-of-principle" report. We first verify the approach in a toy model and in implicitly solvated all-atom polyalanine systems. We next study five proteins, up to 309 residues in size. On the basis of atomistic equilibrium libraries of peptide-plane configurations, the proteins are modeled with fully atomistic backbones and simplified Go-like interactions among residues. We show that full equilibrium sampling can be obtained in days to weeks on a single processor, suggesting that more accurate models are well within reach. For the future, LBMC provides a convenient platform for constructing adjustable or mixed-resolution models: the configurations of all atoms can be stored at no run-time cost, while an arbitrary subset of interactions is "turned on".

Combining an Elastic Network With a Coarse-Grained Molecular Force Field: Structure, Dynamics, and Intermolecular Recognition

Structure-based and physics-based coarse-grained molecular force fields have become attractive approaches to gain mechanistic insight into the function of large biomolecular assemblies. Here, we study how both approaches can be combined into a single representation, that we term ELNEDIN. In this representation an elastic network is used as a structural scaffold to describe and maintain the overall shape of a protein and a physics-based coarse-grained model (MARTINI-2.1) is used to describe both inter- and intramolecular interactions in the system. The results show that when used in molecular dynamics simulations ELNEDIN models can be built so that the resulting structural and dynamical properties of a protein, including its collective motions, are comparable to those obtained using atomistic protein models. We then evaluate the behavior of such models in (1) long, microsecond time-scale, simulations, (2) the modeling of very large macromolecular assemblies, a viral capsid, and (3) the study of a protein-protein association process, the reassembly of the ROP homodimer. The results for this series of tests indicate that ELNEDIN models allow microsecond time-scale molecular dynamics simulations to be carried out readily, that large biological entities such as the viral capsid of the cowpea mosaic virus can be stably modeled as assemblies of independent ELNEDIN models, and that ELNEDIN models show significant promise for modeling protein-protein association processes.

An all-atom structure-based potential for proteins: bridging minimal models with all-atom empirical forcefields

Protein dynamics take place on many time and length scales. Coarse-grained structure-based (Go) models utilize the funneled energy landscape theory of protein folding to provide an understanding of both long time and long length scale dynamics. All-atom empirical forcefields with explicit solvent can elucidate our understanding of short time dynamics with high energetic and structural resolution. Thus, structure-based models with atomic details included can be used to bridge our understanding between these two approaches. We report on the robustness of folding mechanisms in one such all-atom model. Results for the B domain of Protein A, the SH3 domain of C-Src Kinase, and Chymotrypsin Inhibitor 2 are reported. The interplay between side chain packing and backbone folding is explored. We also compare this model to a C(alpha) structure-based model and an all-atom empirical forcefield. Key findings include: (1) backbone collapse is accompanied by partial side chain packing in a cooperative transition and residual side chain packing occurs gradually with decreasing temperature, (2) folding mechanisms are robust to variations of the energetic parameters, (3) protein folding free-energy barriers can be manipulated through parametric modifications, (4) the global folding mechanisms in a C(alpha) model and the all-atom model agree, although differences can be attributed to energetic heterogeneity in the all-atom model, and (5) proline residues have significant effects on folding mechanisms, independent of isomerization effects. Because this structure-based model has atomic resolution, this work lays the foundation for future studies to probe the contributions of specific energetic factors on protein folding and function.