Exploring transition pathway and free-energy profile of large-scale protein conformational change by combining normal mode analysis and umbrella sampling molecular dynamics

A Mode Evolution Metric to Extract Reaction Coordinates for Biomolecular Conformational Transitions

The complex, multidimensional energy landscape of biomolecules makes the extraction of suitable, nonintuitive collective variables (CVs) that describe their conformational transitions challenging. At present, dimensionality reduction approaches and machine learning (ML) schemes are employed to obtain CVs from molecular dynamics (MD)/Monte Carlo (MC) trajectories or structural databanks for biomolecules. However, minimum sampling conditions to generate reliable CVs that accurately describe the underlying energy landscape remain unclear. Here, we address this issue by developing a Mode evolution Metric (MeM) to extract CVs that can pinpoint new states and describe local transitions in the vicinity of a reference minimum from nonequilibrated MD/MC trajectories. We present a general mathematical formulation of MeM for both statistical dimensionality reduction and machine learning approaches. Application of MeM to MC trajectories of model potential energy landscapes and MD trajectories of solvated alanine dipeptide reveals that the principal components which locate new states in the vicinity of a reference minimum emerge well before the trajectories locally equilibrate between the associated states. Finally, we demonstrate a possible application of MeM in designing efficient biased sampling schemes to construct accurate energy landscape slices that link transitions between states. MeM can help speed up the search for new minima around a biomolecular conformational state and enable the accurate estimation of thermodynamics for states lying on the energy landscape and the description of associated transitions.

A theoretical study on the activity and selectivity of IDO/TDO inhibitors

Indoleamine 2,3-dioxygenase 1 (IDO) is a tryptophan (Trp) metabolic enzyme along the kynurenine (NFK) pathway. Under pathological conditions, IDO overexpressed by tumor cells causes depletion of tryptophan and the accumulation of metabolic products, which inhibit the local immune response and form immune escape. Therefore, the suppression of IDO activity is one of the strategies for tumor immunotherapy, and drug design for this target has been the focus of research for more than two decades. Apart from IDO, tryptophan dioxygenase (TDO) of the same family can also catalyze the same biochemical reaction in the human body, but it has different tissue distribution and substrate selectivity from IDO. Based on the principle of drug design with high potency and low cross-reactivity to specific targets, in this subject, the activity and selectivity of IDO and TDO toward small molecular inhibitors were studied from the perspective of thermodynamics and kinetics. The aim was to elucidate the structural requirements for achieving favorable biological activity and selectivity of IDO and TDO inhibitors. Specifically, the interactions of inhibitors from eight families with IDO and TDO were initially investigated through molecular docking and molecular dynamics simulations, and the thermodynamic data for binding of inhibitors were predicted by the molecular mechanics/generalized Born surface area (MM/GBSA) method. Secondly, we explored the free energy landscape of JKloops, the kinetic control element of IDO/TDO, using temperature replica exchange molecular dynamics (T-REMD) simulations and elucidated the connection between the rules of IDO/TDO conformational changes and the inhibitor selectivity mechanism. Furthermore, the binding and dissociation processes of the C1 inhibitor (NLG919) were simulated by the adaptive steering molecular dynamics (ASMD) method, which not only addressed the possible stable, metastable, and transition states for C1 inhibitor-IDO/TDO interactions, but also accurately predicted kinetic data for C1 inhibitor binding and dissociation. In conclusion, we have constructed a complete process from enzyme (IDO/TDO) conformational activation to inhibitor binding/dissociation and used the thermodynamic and kinetic data of each link as clues to verify the control mechanism of IDO/TDO on inhibitor selectivity. This is of great significance for us to understand the design principles of tumor immunotherapy drugs and to avoid drug resistance of immunotherapy drugs.

Energetic and structural insights behind calcium induced conformational transition in calmodulin

Calmodulin (CaM) is a key signaling protein that triggers several cellular and physiological processes inside the cell. Upon binding with calcium ion, CaM undergoes large scale conformational transition from a closed state to an open state that facilitates its interaction with various target protein and regulates their activity. This work explores the origin of the energetic and structural variation of the wild type and mutated CaM and explores the molecular origin for the structural differences between them. We first calculated the sequential calcium binding energy to CaM using the PDLD/S-LRA/β approach. This study  shows a very good correlation with experimental calcium binding energies. Next we calculated the calcium binding energies to the wild type CaM and several mutated CaM systems which were reported experimentally. On the structural aspect, it has been reported experimentally that certain mutation (Q41L-K75I) in calcium bound CaM leads to complete conformational transition from an open to a closed state. By using equilibrium molecular dynamics simulation, free energy calculation and contact frequency map analysis, we have shown that the formation of a cluster of long-range hydrophobic contacts, initiated by the Q41L-K75I CaM variant is the driving force behind its closing motion. This study unravels the energetics and structural aspects behind calcium ion induced conformational changes in wild type CaM and its variant.

Tackling Hysteresis in Conformational Sampling: How to Be Forgetful with MEMENTO

The structure of proteins has long been recognized to hold the key to understanding and engineering their function, and rapid advances in structural biology and protein structure prediction are now supplying researchers with an ever-increasing wealth of structural information. Most of the time, however, structures can only be determined in free energy minima, one at a time. While conformational flexibility may thus be inferred from static end-state structures, their interconversion mechanisms─a central ambition of structural biology─are often beyond the scope of direct experimentation. Given the dynamical nature of the processes in question, many studies have attempted to explore conformational transitions using molecular dynamics (MD). However, ensuring proper convergence and reversibility in the predicted transitions is extremely challenging. In particular, a commonly used technique to map out a path from a starting to a target conformation called steered MD (SMD) can suffer from starting-state dependence (hysteresis) when combined with techniques such as umbrella sampling (US) to compute the free energy profile of a transition. Here, we study this problem in detail on conformational changes of increasing complexity. We also present a new, history-independent approach that we term "MEMENTO" (Morphing End states by Modelling Ensembles with iNdependent TOpologies) to generate paths that alleviate hysteresis in the construction of conformational free energy profiles. MEMENTO utilizes template-based structure modelling to restore physically reasonable protein conformations based on coordinate interpolation (morphing) as an ensemble of plausible intermediates, from which a smooth path is picked. We compare SMD and MEMENTO on well-characterized test cases (the toy peptide deca-alanine and the enzyme adenylate kinase) before discussing its use in more complicated systems (the kinase P38α and the bacterial leucine transporter LeuT). Our work shows that for all but the simplest systems SMD paths should not in general be used to seed umbrella sampling or related techniques, unless the paths are validated by consistent results from biased runs in opposite directions. MEMENTO, on the other hand, performs well as a flexible tool to generate intermediate structures for umbrella sampling. We also demonstrate that extended end-state sampling combined with MEMENTO can aid the discovery of collective variables on a case-by-case basis.

Local Normal Mode Analysis for Fast Loop Conformational Sampling

We propose and validate a novel method to efficiently explore local protein loop conformations based on a new formalism for constrained normal mode analysis (NMA) in internal coordinates. The manifold of possible loop configurations imposed by the position and orientation of the fixed loop ends is reduced to an orthogonal set of motions (or modes) encoding concerted rotations of all the backbone dihedral angles. We validate the sampling power on a set of protein loops with highly variable experimental structures and demonstrate that our approach can efficiently explore the conformational space of closed loops. We also show an acceptable resemblance of the ensembles around equilibrium conformations generated by long molecular simulations and constrained NMA on a set of exposed and diverse loops. In comparison with other methods, the main advantage is the lack of restrictions on the number of dihedrals that can be altered simultaneously. Furthermore, the method is computationally efficient since it only requires the diagonalization of a tiny matrix, and the modes of motions are energetically contextualized by the elastic network model, which includes both the loop and the neighboring residues.

A Free-Energy Landscape Analysis of Calmodulin Obtained from an NMR Data-Utilized Multi-Scale Divide-and-Conquer Molecular Dynamics Simulation

Calmodulin (CaM) is a multifunctional calcium-binding protein, which regulates a variety of biochemical processes. CaM acts through its conformational changes and complex formation with its target enzymes. CaM consists of two globular domains (N-lobe and C-lobe) linked by an extended linker region. Upon calcium binding, the N-lobe and C-lobe undergo local conformational changes, followed by a major conformational change of the entire CaM to wrap the target enzyme. However, the regulation mechanisms, such as allosteric interactions, which regulate the large structural changes, are still unclear. In order to investigate the series of structural changes, the free-energy landscape of CaM was obtained by multi-scale divide-and-conquer molecular dynamics (MSDC-MD). The resultant free-energy landscape (FEL) shows that the Ca²⁺ bound CaM (holo-CaM) would take an experimentally famous elongated structure, which can be formed in the early stage of structural change, by breaking the inter-domain interactions. The FEL also shows that important interactions complete the structural change from the elongated structure to the ring-like structure. In addition, the FEL might give a guiding principle to predict mutational sites in CaM. In this study, it was demonstrated that the movement process of macroscopic variables on the FEL may be diffusive to some extent, and then, the MSDC-MD is suitable to the parallel computation.

Activation of Abl1 Kinase Explored Using Well-Tempered Metadynamics Simulations on an Essential Dynamics Sampled Path

Well-tempered metadynamics (wT-metaD) simulations using path collective variables (CVs) have been successfully applied in recent years to explore conformational transitions in protein kinases and other biomolecular systems. While this methodology has the advantage of describing the transitions with a limited number of predefined path CVs, it requires as an input a reference path connecting the initial and target states of the system. It is desirable to automate the path generation using approaches that do not rely on the choice of geometric CVs to describe the transition of interest. To this end, we developed an approach that couples essential dynamics sampling with wT-metaD simulations. We used this newly developed procedure to explore the activation mechanism of Abl1 kinase and compute the associated free energy barriers. Through these simulations, we identified a three-step mechanism for the activation that involved two metastable intermediates that possessed a partially open activation loop and differed primarily in the “in” or “out” conformation of the aspartate residue of the DFG motif. One of these states is similar to a conformation that was detected in previous spectroscopic studies of Abl1 kinase, albeit its mechanistic role in the activation was hitherto not well understood. The present study establishes its intermediary role in the activation and predicts a rate-determining free energy barrier of 13.8 kcal/mol that is in good agreement with previous experimental and computational estimates. Overall, our study demonstrates the usability of essential dynamics sampling as a path CV in wT-metaD to conveniently study conformational transitions and accurately calculate the associated barriers.

Multiscale computational study of ligand binding pathways: Case of p38 MAP kinase and its inhibitors

Protein kinases are one of the most important drug targets in the past 10 years. Understanding the inhibitor association processes will profoundly impact new binder designs with preferred binding kinetics. However, after more than a decade of effort, a complete atomistic-level study of kinase inhibitor binding pathways is still lacking. As all kinases share a similar scaffold, we used p38 kinase as a model system to investigate the conformational dynamics and free energy transition of inhibitor binding toward kinases. Two major kinase conformations, Asp-Phe-Gly (DFG)-in and DFG-out, and three types of inhibitors, type I, II, and III, were thoroughly investigated in this work. We performed Brownian dynamics simulations and up to 340 μs Gaussian-accelerated molecular dynamics simulations to capture the inhibitor binding paths and a series of conformational transitions of the p38 kinase from its apo to inhibitor-bound form. Eighteen successful binding trajectories, including all types of inhibitors, are reported herein. Our simulations suggest a mechanism of inhibitor recruitment, a faster ligand association step to a pre-existing DFG-in/DFG-out p38 protein, followed by a slower molecular rearrangement step to adjust the protein-ligand conformation followed by a shift in the energy landscape to reach the final bound state. The ligand association processes also reflect the energetic favor of type I and type II/III inhibitor binding through ATP and allosteric channels, respectively. These different binding routes are directly responsible for the fast (type I binders) and slow (type II/III binders) kinetics of different types of p38 inhibitors. Our findings also echo the recent study of p38 inhibitor dissociation, implying that ligand unbinding could undergo a reverse path of binding, and both processes share similar metastates. This study deepens the understanding of molecular and energetic features of kinase inhibitor-binding processes and will inspire future drug development from a kinetic point of view.

Activation Pathways and Free Energy Landscapes of the SARS-CoV-2 Spike Protein

Severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2) uses a spike protein (S-protein) to recognize the receptor protein ACE2 of human cells and initiate infection, during which the conformational transition of the S-protein from inactive (down) state to active (up) state is one of the key molecular events determining the infectivity but the underlying mechanism remains poorly understood. In this work, we investigated the activation pathways and free energy landscape of the S-protein of SARS-CoV-2 and compared with those of the closely related counterpart SARS-CoV using molecular dynamics simulations. Our results revealed a large difference between the activation pathways of the two S-proteins. The transition from inactive to an active state for the S-protein of SARS-CoV-2 is more cooperative, involving simultaneous disruptions of several key interfacial hydrogen bonds, and the transition encounters a much higher free energy barrier. In addition, the conformational equilibrium of the SARS-CoV-2 S-protein is more biased to the inactive state compared to that of the SARS-CoV S-protein, suggesting the transient feature of the active state before binding to the receptor protein of the host cell. The key interactions contributing to the difference of the activation pathways and free energy landscapes were discussed. The results provide insights into the molecular mechanism involved in the initial stage of the SARS-CoV-2 infection.

Combining Virtual Reality Visualization with Ensemble Molecular Dynamics to Study Complex Protein Conformational Changes

Molecular dynamics (MD) simulations are increasingly used to elucidate relationships between protein structure, dynamics, and their biological function. Currently, it is extremely challenging to perform MD simulations of large-scale structural rearrangements in proteins that occur on millisecond timescales or beyond, as this requires very significant computational resources, or the use of cumbersome "collective variable" enhanced sampling protocols. Here, we describe a framework that combines ensemble MD simulations and virtual reality visualization (eMD-VR) to enable users to interactively generate realistic descriptions of large amplitude, millisecond timescale protein conformational changes in proteins. Detailed tests demonstrate that eMD-VR substantially decreases the computational cost of folding simulations of a WW domain, without the need to define collective variables a priori. We further show that eMD-VR generated pathways can be combined with Markov state models to describe the thermodynamics and kinetics of large-scale loop motions in the enzyme cyclophilin A. Our results suggest eMD-VR is a powerful tool for exploring protein energy landscapes in bioengineering efforts.

Increasing the Sampling Efficiency of Protein Conformational Change by Combining a Modified Replica Exchange Molecular Dynamics and Normal Mode Analysis

Understanding conformational change at an atomic level is significant when determining a protein functional mechanism. Replica exchange molecular dynamics (REMD) is a widely used enhanced sampling method to explore protein conformational space. However, REMD with an explicit solvent model requires huge computational resources, immensely limiting its application. In this study, a variation of parallel tempering metadynamics (PTMetaD) with the omission of solvent-solvent interactions in exchange attempts and the use of low-frequency modes calculated by normal-mode analysis (NMA) as collective variables (CVs), namely ossPTMetaD, is proposed with the aim to accelerate MD simulations simultaneously in temperature and geometrical spaces. For testing the performance of ossPTMetaD, five protein systems with diverse biological functions and motion patterns were selected, including large-scale domain motion (AdK), flap movement (HIV-1 protease and BACE1), and DFG-motif flip in kinases (p38α and c-Abl). The simulation results showed that ossPTMetaD requires much fewer numbers of replicas than temperature REMD (T-REMD) with a reduction of ∼70% to achieve a similar exchange ratio. Although it does not obey the detailed balance condition, ossPTMetaD provides consistent results with T-REMD and experimental data. The high accessibility of the large conformational change of protein systems by ossPTMetaD, especially in simulating the very challenging DFG-motif flip of protein kinases, demonstrated its high efficiency and robustness in the characterization of the large-scale protein conformational change pathway and associated free energy profile.

Computational Insights into the Conformational Accessibility and Binding Strength of SARS-CoV-2 Spike Protein to Human Angiotensin-Converting Enzyme 2

The spike protein of SARS-CoV-2 (CoV-2-S) mediates the virus entry into human cells. Experimental studies have shown the stronger binding affinity of the RBD (receptor binding domain) of CoV-2-S to angiotensin-converting enzyme 2 (ACE2) as compared to that of SARS-CoV spike (CoV-S). However, a similar or weaker binding affinity of CoV-2-S compared to that of CoV-S is observed if entire spikes are used in the bioassay. To explore the underlying mechanism, we calculated the binding affinities of the RBDs to ACE2 and simulated the transitions between ACE2-inaccessible and -accessible conformations. We found that the ACE2-accessible angle of CoV-2-S is 52.2° and that the ACE2 binding strength of CoV-2-S RBD is much stronger than that of CoV-S RBD. However, CoV-2-S has much less of an ACE2-accessible conformation and is much more difficult to shift from ACE2-inaccessible to -accessible than CoV-S, making the binding affinity of the entire protein decrease. Further analysis revealed key interactional residues for strong binding and five potential ligand-binding pockets for drug research.

Exploring Conformational Change of Adenylate Kinase by Replica Exchange Molecular Dynamic Simulation

Replica exchange molecular dynamics (REMD) simulation is a popular enhanced sampling method that is widely used for exploring the atomic mechanism of protein conformational change. However, the requirement of huge computational resources for REMD, especially with the explicit solvent model, largely limits its application. In this study, the availability and efficiency of a variant of velocity-scaling REMD (vsREMD) was assessed with adenylate kinase as an example. Although vsREMD achieved results consistent with those from conventional REMD and experimental studies, the number of replicas required for vsREMD (30) was much less than that for conventional REMD (80) to achieve a similar acceptance rate (∼0.2), demonstrating high efficiency of vsREMD to characterize the protein conformational change and associated free-energy profile. Thus, vsREMD is a highly efficient approach for studying the large-scale conformational change of protein systems.

Large-Scale Biomolecular Conformational Transitions Explored by a Combined Elastic Network Model and Enhanced Sampling Molecular Dynamics

Biomolecules often undergo large-scale conformational transitions when carrying out their functions. However, it is still challenging for conventional molecular dynamics simulations to provide adequate structural dynamics information to interpret associated mechanisms. Here, we present a combined elastic network model and enhanced sampling-based strategy (iterANM-IaMD) by adopting iterANM to construct initial conformation space and enhanced sampling IaMD to explore the free energy landscape along specific large-scale conformational transitions. We applied this strategy to three functionally and structurally distinct proteins (adenylate kinase, calmodulin, and p38α kinase), which undergo striking conformational change upon ligand binding. The simulation results for both free and ligand-bound proteins show qualitative and quantitative agreement with existing studies, suggesting iterANM-IaMD as an accurate and efficient tool to investigate structural dynamics involved in complicated biological processes. Our work also provides insights into the relationship between the dynamics and functionality of biomolecules.

Activation and Inactivation of the FLT3 Kinase: Pathway Intermediates and the Free Energy of Transition

The aberrant expression of kinases is often associated with pathologies such as cancer and autoimmune diseases. Like other types of enzymes, kinases can adopt active and inactive states, where a shift toward more stable active state often leads to disease. Dozens of kinase inhibitors are, therefore, used as drugs. Most of these bind to either the inactive or active state. In this work, we study the transitions between these two states in FLT3, an important drug target in leukemias. Kinases are composed of two lobes (N- and C-terminal lobes) with the catalytic site in-between. Through projection of the largest motions obtained through molecular dynamics (MD) simulations, we show that each of the end-states (active or inactive) already possess the ability for transition as the two lobes rotate which initiates the transition. A targeted simulation approach known as essential dynamics sampling (EDS) was used to speed up the transition between the two protein states. Coupling the EDS to implicit-solvent MD was performed to estimate the free energy barriers of the transitions. The activation energies were found in good agreement with previous estimates obtained for other kinases. Finally, we identified FLT3 intermediates that assumed configurations that resemble that of the c-Src nonreceptor tyrosine kinase. The intermediates show better binding to the drug ponatinib than c-Src and the inactive state of FLT3. This suggests that targeting intermediate states can be used to explain the drug-binding patterns of kinases and for rational drug design.

T7 RNA Polymerase Discriminates Correct and Incorrect Nucleoside Triphosphates by Free Energy

RNA polymerase (RNAP) is the primary machine responsible for transcription. Its ability to distinguish between correct (cognate) and incorrect (noncognate) nucleoside triphosphates (NTPs) is important for fidelity control in transcription. In this work, we investigated the substrate selection mechanism of T7 RNAP from the perspective of energetics. The dissociation free energies were determined for matched and unmatched base pairs in the preinsertion complex using the umbrella sampling method. A clear hydrogen-bond-rupture peak is observed in the potential of mean force curve for a matched base pair, whereas no such peaks are present in the position of mean force profiles for unmatched ones. The free-energy barrier could prevent correct substrates from being separated from the active site. Therefore, when NTPs diffuse into the active site, correct ones will stay for chemistry once they establish effective base pairing contacts with the template nucleotide, whereas incorrect ones will be withdrawn from the active site and rejected back to solution. This result provides an important energy evidence for the substrate selection mechanism of RNAP. Then we elucidated energetics and molecular details for correct NTP binding to the active site of the insertion complex. Our observations reveal that strong interactions act on the triphosphate of NTP to constrain its movement, whereas relatively weak interactions serve to position the base in the correct conformation. Triple interactions, hydrophobic contacts from residues M635 and Y639, base stacking from the 3' RNA terminal nucleotide, and base pairing from the template nucleotide act together to position the NTP base in a catalytically competent conformation. At last, we observed that incorrect NTPs cannot be as well-stabilized as the correct one in the active site when they are misincorporated in the insertion site. It is expected that our work can be helpful for comprehensively understanding details of this basic step in genetic transcription.

How calcium ion binding induces the conformational transition of the calmodulin N-terminal domain—an atomic level characterization

The Ca²⁺binding and triggering conformation transition of nCaM were detected in unbiased molecular dynamics simulations.

Different structures of berberine and five other protoberberine alkaloids that affect P-glycoprotein-mediated efflux capacity

Berberine, berberrubine, thalifendine, demethyleneberberine, jatrorrhizine, and columbamine are six natural protoberberine alkaloid (PA) compounds that display extensive pharmacological properties and share the same protoberberine molecular skeleton with only slight substitution differences. The oral delivery of most PAs is hindered by their poor bioavailability, which is largely caused by P-glycoprotein (P-gp)-mediated drug efflux. Meanwhile, P-gp undergoes large-scale conformational changes (from an inward-facing to an outward-facing state) when transporting substrates, and these changes might strongly affect the P-gp-binding specificity. To confirm whether these six compounds are substrates of P-gp, to investigate the differences in efflux capacity caused by their trivial structural differences and to reveal the key to increasing their binding affinity to P-gp, we conducted a series of in vivo, in vitro, and in silico assays. Here, we first confirmed that all six compounds were substrates of P-gp by comparing the drug concentrations in wild-type and P-gp-knockout mice in vivo. The efflux capacity (net efflux) ranked as berberrubine > berberine > columbamine ~ jatrorrhizine > thalifendine > demethyleneberberine based on in vitro transport studies in Caco-2 monolayers. Using molecular dynamics simulation and molecular docking techniques, we determined the transport pathways of the six compounds and their binding affinities to P-gp. The results suggested that at the early binding stage, different hydrophobic and electrostatic interactions collectively differentiate the binding affinities of the compounds to P-gp, whereas electrostatic interactions are the main determinant at the late release stage. In addition to hydrophobic interactions, hydrogen bonds play an important role in discriminating the binding affinities.

Exploring molecular mechanism of allosteric inhibitor to relieve drug resistance of multiple mutations in HIV-1 protease by enhanced conformational sampling

Recently, allosteric regulations of HIV-1 protease (PR) are suggested as a promising approach to relieve drug resistance of mutations toward inhibitors targeting the active site of PR. Replica-exchange molecular dynamics (REMD) simulations and normal mode analysis (NMA) are integrated to enhance conformational sampling of PR. Molecular mechanics generalized Born surface area (MM-GBSA) method was applied to calculate binding free energies of three inhibitors APV, DRV, and NIT to the wild-type (WT) and multidrug resistance (MDR) PRs. The results suggest that binding free energies of APV and DRV are decreased in the MDR PR relative to the WT PR, suggesting drug resistance of mutations on these two inhibitors. However, the binding ability of the allosteric inhibitor NIT is not impaired in the MDR PR. In addition, internal dynamics analysis based on REMD simulations proves that mutations hardly produce obvious effect on the conformation of the MDR PR in comparison to the WT PR. Scanning of hydrophobic contacts and hydrogen bond contacts of inhibitors with residues of PRs on the concatenated trajectories of REMD demonstrates that mutations change the symmetric interaction networks of APV and DRV with PR, but do not generate obvious influence on the asymmetric interaction network of NIT with PR. In summary, allosteric inhibitor NIT can adapt the MDR PR better than those inhibitors toward the active site of PR, thus allosteric inhibitors of PR may be a possible channel to overcome drug resistance of PR.

Zinc ion-induced conformational changes in new Delphi metallo-β-lactamase 1 probed by molecular dynamics simulations and umbrella sampling

The hydrolysis of a β-lactam core ring caused by new Delphi metallo-β-lactamase 1 (NDM-1) with the help of two zinc cofactors induces significant resistance toward β-lactam antibiotics. Molecular dynamics (MD) simulations and the umbrella sampling method are integrated to study the conformational change mechanism of NDM-1 mediated by zinc ion binding. The statistical analyses of interaction contacts of the antibiotic ampicillin (AMP) with residues based on MD trajectories suggest that two Zn ions are essential for maintaining the binding of AMP with NDM-1. Umbrella sampling simulations further reveal that double-Zn coordination exerts strong restriction on the motions of loop L10 relative to loops L3 and L4. Principal component (PC) analysis also demonstrates that zinc ion binding totally inhibits the motion extent of NDM-1 and changes internal motion modes in NDM-1. We expect that the current study can provide significant dynamical information involving conformational changes of NDM-1 for the development of efficient inhibitors to decrease drug resistance of NDM-1 toward antibiotics.

Computational Exploration of Conformational Transitions in Protein Drug Targets

Protein drug targets vary from highly structured to completely disordered; either way dynamics governs function. Hence, understanding the dynamical aspects of how protein targets function can enable improved interventions with drug molecules. Computational approaches offer highly detailed structural models of protein dynamics which are becoming more predictive as model quality and sampling power improve. However, the most advanced and popular models still have errors owing to imperfect parameter sets and often cannot access longer timescales of many crucial biological processes. Experimental approaches offer more certainty but can struggle to detect and measure lightly populated conformations of target proteins and subtle allostery. An emerging solution is to integrate available experimental data into advanced molecular simulations. In the future, molecular simulation in combination with experimental data may be able to offer detailed models of important drug targets such that improved functional mechanisms or selectivity can be accessed.

Hybrid All-Atom/Coarse-Grained Simulations of Proteins by Direct Coupling of CHARMM and PRIMO Force Fields

Hybrid all-atom/coarse-grained (AA/CG) simulations of proteins offer a computationally efficient compromise where atomistic details are only applied to biologically relevant regions while benefiting from the speedup of treating the remaining parts of a given system at the CG level. The recently developed CG model, PRIMO, allows a direct coupling with an atomistic force field with no additional modifications or coupling terms and the ability to carry out dynamic simulations without any restraints on secondary or tertiary structures. A hybrid AA/CG scheme based on combining all-atom CHARMM and coarse-grained PRIMO representations was validated via molecular dynamics and replica exchange simulations of soluble and membrane proteins. The AA/CG scheme was also tested in the calculation of the free energy profile for the transition from the closed to the open state of adenylate kinase via umbrella sampling molecular dynamics method. The overall finding is that the AA/CG scheme generates dynamics and energetics that are qualitatively and quantitatively comparable to AA simulations while offering the computational advantages of coarse-graining. This model opens the door to challenging applications where high accuracy is required only in parts of large biomolecular complexes.

Characterization of Zebrafish Cardiac and Slow Skeletal Troponin C Paralogs by MD Simulation and ITC

Zebrafish, as a model for teleost fish, have two paralogous troponin C (TnC) genes that are expressed in the heart differentially in response to temperature acclimation. Upon Ca(2+) binding, TnC changes conformation and exposes a hydrophobic patch that interacts with troponin I and initiates cardiac muscle contraction. Teleost-specific TnC paralogs have not yet been functionally characterized. In this study we have modeled the structures of the paralogs using molecular dynamics simulations at 18°C and 28°C and calculated the different Ca(2+)-binding properties between the teleost cardiac (cTnC or TnC1a) and slow-skeletal (ssTnC or TnC1b) paralogs through potential-of-mean-force calculations. These values are compared with thermodynamic binding properties obtained through isothermal titration calorimetry (ITC). The modeled structures of each of the paralogs are similar at each temperature, with the exception of helix C, which flanks the Ca(2+) binding site; this region is also home to paralog-specific sequence substitutions that we predict have an influence on protein function. The short timescale of the potential-of-mean-force calculation precludes the inclusion of the conformational change on the ΔG of Ca(2+) interaction, whereas the ITC analysis includes the Ca(2+) binding and conformational change of the TnC molecule. ITC analysis has revealed that ssTnC has higher Ca(2+) affinity than cTnC for Ca(2+) overall, whereas each of the paralogs has increased affinity at 28°C compared to 18°C. Microsecond-timescale simulations have calculated that the cTnC paralog transitions from the closed to the open state more readily than the ssTnC paralog, an unfavorable transition that would decrease the ITC-derived Ca(2+) affinity while simultaneously increasing the Ca(2+) sensitivity of the myofilament. We propose that the preferential expression of cTnC at lower temperatures increases myofilament Ca(2+) sensitivity by this mechanism, despite the lower Ca(2+) affinity that we have measured by ITC.

Protein interaction evolution from promiscuity to specificity with reduced flexibility in an increasingly complex network

A key question regarding protein evolution is how proteins adapt to the dynamic environment in which they function and how in turn their evolution shapes the protein interaction network. We used extant and resurrected ancestral plant MADS-domain transcription factors to understand how SEPALLATA3, a protein with hub and glue properties, evolved and takes part in network organization. Although the density of dimeric interactions was saturated in the network, many new interactions became mediated by SEPALLATA3 after a whole genome triplication event. By swapping SEPALLATA3 and its ancestors between dimeric networks of different ages, we found that the protein lost the capacity of promiscuous interaction and acquired specificity in evolution. This was accompanied with constraints on conformations through proline residue accumulation, which made the protein less flexible. SHORT VEGETATIVE PHASE on the other hand (non-hub) was able to gain protein-protein interactions due to a C-terminal domain insertion, allowing for a larger interaction interface. These findings illustrate that protein interaction evolution occurs at the level of conformational dynamics, when the binding mechanism concerns an induced fit or conformational selection. Proteins can evolve towards increased specificity with reduced flexibility when the complexity of the protein interaction network requires specificity.

Energetics and structural characterization of the “DFG-flip” conformational transition of B-RAF kinase: a SITS molecular dynamics study

A combination of a homology modeling technique and an enhanced sampling molecular dynamics simulation implemented using the SITS method is employed to compute a detailed map of the free-energy landscape and explore the conformational transition pathway of B-RAF kinase.

Mutation L1196M-induced conformational changes and the drug resistant mechanism of anaplastic lymphoma kinase studied by free energy perturbation and umbrella sampling

Anaplastic lymphoma kinase (ALK) has been regarded as a promising drug target in the treatment of tumors and the mutation L1196M induces different levels of drug resistance toward the existing inhibitors.

Molecular Mechanism and Energy Basis of Conformational Diversity of Antibody SPE7 Revealed by Molecular Dynamics Simulation and Principal Component Analysis

More and more researchers are interested in and focused on how a limited repertoire of antibodies can bind and correspondingly protect against an almost limitless diversity of invading antigens. In this work, a series of 200-ns molecular dynamics (MD) simulations followed by principal component (PC) analysis and free energy calculations were performed to probe potential mechanism of conformational diversity of antibody SPE7. The results show that the motion direction of loops H3 and L3 is different relative to each other, implying that a big structural difference exists between these two loops. The calculated energy landscapes suggest that the changes in the backbone angles ψ and φ of H-Y101 and H-Y105 provide significant contributions to the conformational diversity of SPE7. The dihedral angle analyses based on MD trajectories show that the side-chain conformational changes of several key residues H-W33, H-Y105, L-Y34 and L-W93 around binding site of SPE7 play a key role in the conformational diversity of SPE7, which gives a reasonable explanation for potential mechanism of cross-reactivity of single antibody toward multiple antigens.

Clarifying binding difference of ATP and ADP to extracellular signal-regulated kinase 2 by using molecular dynamics simulations

Extracellular signal-regulated kinase 2 is a promising target for designs and development of anticancer drugs. Molecular dynamics simulations and molecular mechanics Poisson-Boltzmann method were applied to study binding difference of ADP and ATP to extracellular signal-regulated kinase 2. The results prove that the binding ability of ATP to extracellular signal-regulated kinase 2 is stronger than that of ADP. Principal component analysis performed by using molecular dynamics trajectories suggests that binding of ADP and ATP to extracellular signal-regulated kinase 2 change motion directions of two helices α1 and α2. Residue-based free energy decomposition method was adopted to calculate contributions of separate residues to associations of ADP and ATP with extracellular signal-regulated kinase 2. The results show that ADP and ATP produce strong CH-π interactions with five residues Ile29, Val37, Ala50, Leu105, and Leu154. In addition, five hydrogen bonding interactions of ADP and ATP with residues Lys52, Gln103, Asp104, and Met106 also stabilize bindings of ADP and ATP to extracellular signal-regulated kinase 2. Overall, the CH-π interactions of ATP with five residues Ile29, Val37, Ala50, Leu105, and Leu154 are stronger than ADP. This study is expected to contribute a significant theoretical hint for designs of anticancer drugs targeting extracellular signal-regulated kinase 2.

Exploring Flexibility of Progesterone Receptor Ligand Binding Domain Using Molecular Dynamics

Progesterone receptor (PR), a member of nuclear receptor (NR) superfamily, plays a vital role for female reproductive tissue development, differentiation and maintenance. PR ligand, such as progesterone, induces conformation changes in PR ligand binding domain (LBD), thus mediates subsequent gene regulation cascades. PR LBD may adopt different conformations upon an agonist or an antagonist binding. These different conformations would trigger distinct transcription events. Therefore, the dynamics of PR LBD would be of general interest to biologists for a deep understanding of its structure-function relationship. However, no apo-form (non-ligand bound) of PR LBD model has been proposed either by experiments or computational methods so far. In this study, we explored the structural dynamics of PR LBD using molecular dynamics simulations and advanced sampling tools in both ligand-bound and the apo-forms. Resolved by the simulation study, helix 11, helix 12 and loop 895–908 (the loop between these two helices) are quite flexible in antagonistic conformation. Several residues, such as Arg899 and Glu723, could form salt-bridging interaction between helix 11 and helix 3, and are important for the PR LBD dynamics. And we also propose that helix 12 in apo-form PR LBD, not like other NR LBDs, such as human estrogen receptor α (ERα) LBD, may not adopt a totally extended conformation. With the aid of umbrella sampling and metadynamics simulations, several stable conformations of apo-form PR LBD have been sampled, which may work as critical structural models for further large scale virtual screening study to discover novel PR ligands for therapeutic application.

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.

Quarterly intrinsic disorder digest (January-February-March, 2014)

This is the 5^th issue of the Digested Disorder series that represents a reader's digest of the scientific literature on intrinsically disordered proteins. We continue to use only 2 criteria for inclusion of a paper to this digest: The publication date (a paper should be published within the covered time frame) and the topic (a paper should be dedicated to any aspect of protein intrinsic disorder). The current digest issue covers papers published during the first quarter of 2014; i.e., during the period of January, February, and March of 2014. Similar to previous issues, the papers are grouped hierarchically by topics they cover, and for each of the included papers a short description is given on its major findings.

Applying high-performance computing in drug discovery and molecular simulation

In recent decades, high-performance computing (HPC) technologies and supercomputers in China have significantly advanced, resulting in remarkable achievements. Computational drug discovery and design, which is based on HPC and combines pharmaceutical chemistry and computational biology, has become a critical approach in drug research and development and is financially supported by the Chinese government. This approach has yielded a series of new algorithms in drug design, as well as new software and databases. This review mainly focuses on the application of HPC to the fields of drug discovery and molecular simulation at the Chinese Academy of Sciences, including virtual drug screening, molecular dynamics simulation, and protein folding. In addition, the potential future application of HPC in precision medicine is briefly discussed.

Drug resistance mechanisms of three mutations V32I, I47V and V82I in HIV-1 protease toward inhibitors probed by molecular dynamics simulations and binding free energy predictions

Molecular dynamics simulation and binding free energy calculations were used to probe drug resistance of HIV-1 protease mutations toward inhibitors.

New generation of elastic network models

The intrinsic flexibility of proteins and nucleic acids can be grasped from remarkably simple mechanical models of particles connected by springs. In recent decades, Elastic Network Models (ENMs) combined with Normal Model Analysis widely confirmed their ability to predict biologically relevant motions of biomolecules and soon became a popular methodology to reveal large-scale dynamics in multiple structural biology scenarios. The simplicity, robustness, low computational cost, and relatively high accuracy are the reasons behind the success of ENMs. This review focuses on recent advances in the development and application of ENMs, paying particular attention to combinations with experimental data. Successful application scenarios include large macromolecular machines, structural refinement, docking, and evolutionary conservation.

Conformational Equilibrium of CDK/Cyclin Complexes by Molecular Dynamics with Excited Normal Modes

Cyclin-dependent kinases (CDKs) and their associated regulatory cyclins are central for timely regulation of cell-cycle progression. They constitute attractive pharmacological targets for development of anticancer therapeutics, since they are frequently deregulated in human cancers and contribute to sustained, uncontrolled tumor proliferation. Characterization of their structural/dynamic features is essential to gain in-depth insight into structure-activity relationships. In addition, the identification of druggable pockets or key intermediate conformations yields potential targets for the development of novel classes of inhibitors. Structural studies of CDK2/cyclin A have provided a wealth of information concerning monomeric/heterodimeric forms of this kinase. There is, however, much less structural information for other CDK/cyclin complexes, including CDK4/cyclin D1, which displays an alternative (open) position of the cyclin partner relative to CDK, contrasting with the closed CDK2/cyclin A conformation. In this study, we carried out normal-mode analysis and enhanced sampling simulations with our recently developed method, molecular dynamics with excited normal modes, to understand the conformational equilibrium on these complexes. Interestingly, the lowest-frequency normal mode computed for each complex described the transition between the open and closed conformations. Exploration of these motions with an explicit-solvent representation using molecular dynamics with excited normal modes confirmed that the closed conformation is the most stable for the CDK2/cyclin A complex, in agreement with their experimentally available structures. On the other hand, we clearly show that an open↔closed equilibrium may exist in CDK4/cyclin D1, with closed conformations resembling that captured for CDK2/cyclin A. Such conformational preferences may result from the distinct distributions of frustrated contacts in each complex. Using the same approach, the putative roles of the Thr(160) phosphoryl group and the T-loop conformation were investigated. These results provide a dynamic view of CDKs revealing intermediate conformations not yet characterized for CDK members other than CDK2, which will be useful for the design of inhibitors targeting critical conformational transitions.

Analysis of the Free Energy Landscapes for the Opening-Closing Dynamics of the Maltose Transporter ATPase MalK2 Using Enhanced-Sampling Molecular Dynamics Simulation

Protein dynamics are considered significant for many physiological processes, such as metabolism, biomolecular recognition, and the regulation of several vital cellular processes. Due to their flexibility, proteins may stay in different substates with or without the existence of the cognate substrates. To describe these phenomena, two models have been proposed: the "induced fit" and the "conformational selection" mechanisms. In this study, we used MalK2, the subunits that mainly include the nucleotide-binding domains (NBDs) of the maltose transporter from Escherichia coli, as a target to understand the NBD dimerization mechanism. Accelerated and conventional molecular dynamics have been performed. The results revealed that Mg-ATP binding to MalK2 led to a significant change in the free energy profile and thus stabilized the closed conformation. On the contrary, when Mg-ATP was removed, the open conformation would be favored. The fact that ligand binding induces a drastic free energy change leads to a significant inference: MalK2 dimerization would occur through the induced-fit mechanism rather than the conformational selection mechanism. This study sheds new light on the NBD dimerization mechanism and would be of wide applicability to other ABC transporters.

Increasing the sampling efficiency of protein conformational transition using velocity-scaling optimized hybrid explicit/implicit solvent REMD simulation

The application of temperature replica exchange molecular dynamics (REMD) simulation on protein motion is limited by its huge requirement of computational resource, particularly when explicit solvent model is implemented. In the previous study, we developed a velocity-scaling optimized hybrid explicit/implicit solvent REMD method with the hope to reduce the temperature (replica) number on the premise of maintaining high sampling efficiency. In this study, we utilized this method to characterize and energetically identify the conformational transition pathway of a protein model, the N-terminal domain of calmodulin. In comparison to the standard explicit solvent REMD simulation, the hybrid REMD is much less computationally expensive but, meanwhile, gives accurate evaluation of the structural and thermodynamic properties of the conformational transition which are in well agreement with the standard REMD simulation. Therefore, the hybrid REMD could highly increase the computational efficiency and thus expand the application of REMD simulation to larger-size protein systems.

Revealing the favorable dissociation pathway of type II kinase inhibitors via enhanced sampling simulations and two-end-state calculations

How does a type II inhibitor bind to/unbind from a kinase target is still a confusing question because the small molecule occupies both the ATP pocket and the allosteric pocket of the kinase binding site. Here, by using enhanced sampling simulations (umbrella sampling, US) and two-end-state free energy calculations (MM/GSBA), we systemically studied the dissociation processes of two distinct small molecules escaping from the binding pocket of p38 MAP kinase through the allosteric channel and the ATP channel. The results show that the unbinding pathways along the allosteric channel have much lower PMF depths than those along the ATP channel, suggesting that the allosteric channel is more favorable for the dissociations of the two inhibitors and thereby supporting the general understanding that the largest channel of a target is usually the entry/exit pathway for the binding/dissociation of small molecules. Interestingly, the MM/GBSA approach yielded similar PMF profiles compared with those based on US, a much time consuming approach, indicating that for a general study, such as detecting the important transition state of a ligand binding/unbinding process, MM/GBSA may be a feasible choice.

A comparative study of trypsin specificity based on QM/MM molecular dynamics simulation and QM/MM GBSA calculation

Hydrogen bonding and polar interactions play a key role in identification of protein-inhibitor binding specificity. Quantum mechanics/molecular mechanics molecular dynamics (QM/MM MD) simulations combined with DFT and semi-empirical Hamiltonian (AM1d, RM1, PM3, and PM6) methods were performed to study the hydrogen bonding and polar interactions of two inhibitors BEN and BEN1 with trypsin. The results show that the accuracy of treating the hydrogen bonding and polar interactions using QM/MM MD simulation of PM6 can reach the one obtained by the DFT QM/MM MD simulation. Quantum mechanics/molecular mechanics generalized Born surface area (QM/MM-GBSA) method was applied to calculate binding affinities of inhibitors to trypsin and the results suggest that the accuracy of binding affinity prediction can be significantly affected by the accurate treatment of the hydrogen bonding and polar interactions. In addition, the calculated results also reveal the binding specificity of trypsin: (1) the amidinium groups of two inhibitors generate favorable salt bridge interaction with Asp189 and form hydrogen bonding interactions with Ser190 and Gly214, (2) the phenyl of inhibitors can produce favorable van der Waals interactions with the residues His58, Cys191, Gln192, Trp211, Gly212, and Cys215. This systematic and comparative study can provide guidance for the choice of QM/MM MD methods and the designs of new potent inhibitors targeting trypsin.