General Trends of Dihedral Conformational Transitions in a Globular ProteinProteins, Volume 84, Issue 4, pp 501-514 (2016)
Dihedral conformational transitions are analyzed systematically in a model globular protein, cytochrome P450cam, to examine their structural and chemical dependences through combined conventional molecular dynamics (cMD), accelerated molecular dynamics (aMD) and Adaptive Biasing Force (ABF) simulations. The aMD simulations are performed at two acceleration levels, using dihedral and dual boost, respectively. In comparison with cMD, aMD samples protein dihedral transitions ∼2 times faster on average using dihedral boost, and ∼3.5 times faster using dual boost. In the protein backbone, significantly higher dihedral transition rates are observed in the Bend, Coil and Turn flexible regions, followed by the β bridge and β sheet, and then the helices. Moreover, protein sidechains of greater length exhibit higher transition rates on average in the aMD-enhanced sampling. Sidechains of the same length (particularly Nχ = 2) exhibit decreasing transition rates with residues when going from hydrophobic to polar, then charged and aromatic chemical types. The reduction of dihedral transition rates is found to be correlated with increasing energy barriers as identified through ABF free energy calculations. These general trends of dihedral conformational transitions provide important insights into the hierarchical dynamics and complex free energy landscapes of functional proteins. This article is protected by copyright. All rights reserved
Computation of pH-dependent Binding Free EnergiesBiopolymers, Volume 105, Issue 1, pp 43-49 (2016).
Protein–ligand binding accompanies changes in the surrounding electrostatic environments of the two binding partners and may lead to changes in protonation upon binding. In cases where the complex formation results in a net transfer of protons, the binding process is pH-dependent. However, conventional free energy computations or molecular docking protocols typically employ fixed protonation states for the titratable groups in both binding partners set a priori, which are identical for the free and bound states. In this review, we draw attention to these important yet largely ignored binding-induced protonation changes in protein–ligand association by outlining physical origins and prevalence of the protonation changes upon binding. Following a summary of various theoretical methods for pKa prediction, we discuss the theoretical framework to examine the pH dependence of protein–ligand binding processes.
Substrate Channeling Between the Human Dihydrofolate Reductase and Thymidylate SynthaseProt. Sci 25, 79–86 (2016). Special issue honoring Ronald Levy.
In vivo, as an advanced catalytic strategy, transient non-covalently bound multi-enzyme complexes can be formed to facilitate the relay of substrates, i. e. substrate channeling, between sequential enzymatic reactions and to enhance the throughput of multi-step enzymatic pathways. The human thymidylate synthase and dihydrofolate reductase catalyze two consecutive reactions in the folate metabolism pathway, and experiments have shown that they are very likely to bind in the same multi-enzyme complex in vivo. While reports on the protozoa thymidylate synthase-dihydrofolate reductase bifunctional enzyme give substantial evidences of substrate channeling along a surface “electrostatic highway,” attention has not been paid to whether the human thymidylate synthase and dihydrofolate reductase, if they are in contact with each other in the multi-enzyme complex, are capable of substrate channeling employing surface electrostatics. This work utilizes protein–protein docking, electrostatics calculations, and Brownian dynamics to explore the existence and mechanism of the substrate channeling between the human thymidylate synthase and dihydrofolate reductase. The results show that the bound human thymidylate synthase and dihydrofolate reductase are capable of substrate channeling and the formation of the surface “electrostatic highway.” The substrate channeling efficiency between the two can be reasonably high and comparable to that of the protozoa.
Unconstrained Enhanced Sampling for Free Energy Calculations of Biomolecules: A ReviewMolec. Simul, Volume 42, Issue 13, pp1046-1055 (2016)
Free energy calculations are central to understanding the structure, dynamics and function of biomolecules. Yet insufficient sampling of biomolecular configurations is often regarded as one of the main sources of error. Many enhanced sampling techniques have been developed to address this issue. Notably, enhanced sampling methods based on biasing collective variables (CVs), including the widely used umbrella sampling, adaptive biasing force and metadynamics, have been discussed in a recent excellent review (Abrams and Bussi, Entropy, 2014). Here, we aim to review enhanced sampling methods that do not require predefined system-dependent CVs for biomolecular simulations and as such do not suffer from the hidden energy barrier problem as encountered in the CV-biasing methods. These methods include, but are not limited to, replica exchange/parallel tempering, self-guided molecular/Langevin dynamics, essential energy space random walk and accelerated molecular dynamics. While it is overwhelming to describe all details of each method, we provide a summary of the methods along with the applications and offer our perspectives. We conclude with challenges and prospects of the unconstrained enhanced sampling methods for accurate biomolecular free energy calculations.
Hybrid Finite Element and Brownian Dynamics Method for Charged Particles.J. Chem. Phys., 144, 164107 (2016).
Diffusion is often the rate-determining step in many biological processes. Currently, the two main computational methods for studying diffusion are stochastic methods, such as Brownian dynamics, and continuum methods, such as the finite element method. A previous study introduced a new hybrid diffusion method that couples the strengths of each of these two methods, but was limited by the lack of interactions among the particles; the force on each particle had to be from an external field. This study further develops the method to allow charged particles. The method is derived for a general multidimensional system and is presented using a basic test case for a one-dimensional linear system with one charged species and a radially symmetric system with three charged species.
Development of Potent and Selective Inhibitors for Group VIA Calcium-Independent Phospholipase A2 Guided by Molecular Dynamics and Structure-Activity Relationships.J. Med. Chem. 59 (9), pp 4403-4414 (2016)
The development of inhibitors for phospholipases A2 (PLA2s) is important in elucidating their implication in various biological pathways. PLA2 enzymes are an important pharmacological target implicated in various inflammatory diseases. Computational chemistry, organic synthesis and in vitro assays were employed to develop potent and selective inhibitors for Group VIA calcium-independent PLA2. A set of fluoroketone inhibitors, were studied for their binding mode with two human cytosolic PLA enzymes: Group IVA cPLA2 and Group VIA iPLA2. New compounds were synthesized and tested against three major PLA2 enzyme. This study led to the development of four potent and selective thioether fluoroketone inhibitors as well as a thioether keto-1,2,4-oxadiazole inhibitor for GVIA iPLA2, which will serve as lead compounds for future development and studies. The keto-1,2,4-oxadiazole functionality with a thioether is a novel structure and it will be used as a lead to develop inhibitors with higher potency and selectivity towards GVIA iPLA2.
Molecular dynamic study of MlaC protein in Gram-negative bacteria: conformational flexibility, solvent effect and protein-phospholipid bindingProt. Sci. Volume 25, Issue 8, pp 1430-1437 (2016)
The development of inhibitors for phospholipase A2 (PLA2) is important in elucidating the enzymes implication in various biological pathways. PLA2 enzymes are an important pharmacological target implicated in various inflammatory diseases. Computational chemistry, organic synthesis, and in vitro assays were employed to develop potent and selective inhibitors for group VIA calcium-independent PLA2. A set of fluoroketone inhibitors was studied for their binding mode with two human cytosolic PLA2 enzymes: group IVA cPLA2 and group VIA iPLA2. New compounds were synthesized and assayed toward three major PLA2s. This study led to the development of four potent and selective thioether fluoroketone inhibitors as well as a thioether keto-1,2,4-oxadiazole inhibitor for GVIA iPLA2, which will serve as lead compounds for future development and studies. The keto-1,2,4-oxadiazole functionality with a thioether is a novel structure, and it will be used as a lead to develop inhibitors with higher potency and selectivity toward GVIA iPLA2.
Computer-Aided Drug Design Guided by Hydrogen/Deuterium Exchange Mass Spectrometry: A Powerful Combination for the Development of Potent and Selective Inhibitors of Group VIA Calcium-Independent Phospholipase A2Bioorg. Med. Chem. (Jorgensen special issue) 24 (20), 4801-4811 (2016)
Potent and selective inhibitors for phospholipases A2 (PLA2) are useful for studying their intracellular functions. PLA2 enzymes liberate arachidonic acid from phospholipids activating eicosanoid pathways that involve cyclooxygenase (COX) and lipoxygenase (LOX) leading to inflammation. Anti-inflammatory drugs target COX and LOX; thus, PLA2 can also be targeted to diminish inflammation at an earlier stage in the process. This paper describes the employment of enzymatic assays, hydrogen/deuterium exchange mass spectrometry (DXMS) and computational chemistry to develop PLA2 inhibitors. Beta-thioether trifluoromethylketones (TFKs) were screened against human GVIA calcium-independent, GIVA cytosolic and GV secreted PLA2s. These compounds exhibited inhibition toward Group VIA calcium-independent PLA2 (GVIA iPLA2), with the most potent and selective inhibitor 3 (OTFP) obtaining an XI(50) of 0.0002 mole fraction (IC50 of 110 nM). DXMS binding experiments in the presence of OTFP revealed the peptide regions of GVIA iPLA2 that interact with the inhibitor. Molecular docking and dynamics simulations in the presence of a membrane were guided by the DXMS data in order to identify the binding mode of OTFP. Clustering analysis showed the binding mode of OTFP that occupied 70% of the binding modes occurring during the simulation. The resulted 3D complex was used for docking studies and a structure–activity relationship (SAR) was established. This paper describes a novel multidisciplinary approach in which a 3D complex of GVIA iPLA2 with an inhibitor is reported and validated by experimental data. The SAR showed that the sulfur atom is vital for the potency of beta-thioether analogues, while the hydrophobic chain is important for selectivity. This work constitutes the foundation for further design, synthesis and inhibition studies in order to develop new beta-thioether analogues that are potent and selective for GVIA iPLA2 exclusively.
G-Protein Coupled Receptors: Advances in Simulation and Drug DiscoveryCurr. Opin. Struct. Biol. 41: 83-89 (2016).
G-protein coupled receptors (GPCRs), the largest family of human membrane proteins, mediate cellular signaling and represent primary targets of about one third of currently marketed drugs. GPCRs undergo highly dynamic structural transitions during signal transduction, from binding of extracellular ligands to coupling with intracellular effector proteins. Molecular dynamics (MD) simulations have been utilized to investigate GPCR signaling mechanisms (such as pathways of ligand binding and receptor activation/deactivation) and to design novel small-molecule drug candidates. Future research directions point towards modeling cooperative binding of multiple orthosteric and allosteric ligands to GPCRs, GPCR oligomerization and interactions of GPCRs with different intracellular signaling proteins. Through methodological and supercomputing advances, MD simulations will continue to provide important insights into GPCR signaling mechanisms and further facilitate structure-based drug design.
Stochastic Level-Set Variational Implicit-Solvent Approach to Solute-Solvent Interfacial FluctuationsJ. Chem. Phys, 145 (5) (2016).
Recent years have seen the initial success of a variational implicit-solvent model (VISM), implemented with a robust level-set method, in capturing efficiently different hydration states and providing quantitatively good estimation of solvation free energies of biomolecules. The level-set minimization of the VISM solvation free-energy functional of all possible solute-solvent interfaces or dielectric boundaries predicts an equilibrium biomolecular conformation that is often close to an initial guess. In this work, we develop a theory in the form of Langevin geometrical flow to incorporate solute-solvent interfacial fluctuations into the VISM. Such fluctuations are crucial to biomolecular conformational changes and binding process. We also develop a stochastic level-set method to numerically implement such a theory. We describe the interfacial fluctuation through the "normal velocity" that is the solute-solvent interfacial force, derive the corresponding stochastic level-set equation in the sense of Stratonovich so that the surface representation is independent of the choice of implicit function, and develop numerical techniques for solving such an equation and processing the numerical data. We apply our computational method to study the dewetting transition in the system of two hydrophobic plates and a hydrophobic cavity of a synthetic host molecule cucurbituril. Numerical simulations demonstrate that our approach can describe an underlying system jumping out of a local minimum of the free-energy functional and can capture dewetting transitions of hydrophobic systems. In the case of two hydrophobic plates, we find that the wavelength of interfacial fluctuations has a strong influence to the dewetting transition. In addition, we find that the estimated energy barrier of the dewetting transition scales quadratically with the inter-plate distance, agreeing well with existing studies of molecular dynamics simulations. Our work is a first step toward the inclusion of fluctuations into the VISM and understanding the impact of interfacial fluctuations on biomolecular solvation with an implicit-solvent approach.
Accelerated structure-based design of chemically diverse allosteric modulators of a muscarinic G protein-coupled receptorProc. Natl. Acad. Sci. USA in press (2016).
Design of ligands that provide receptor selectivity has emerged as a new paradigm for drug discovery of G protein-coupled receptors, and may, for certain families of receptors, only be achieved via identification of chemically diverse allosteric modulators. Here, the extracellular vestibule of the M2 muscarinic acetylcholine receptor (mAChR) is targeted for structure-based design of allosteric modulators. Accelerated molecular dynamics (aMD) simulations were performed to construct structural ensembles that account for the receptor flexibility. Compounds obtained from the National Cancer Institute (NCI) were docked to the receptor ensembles. Retrospective docking of known ligands showed that combining aMD simulations with Glide induced fit docking (IFD) provided much improved enrichment factors compared with the Glide virtual screening workflow. Glide IFD was thus applied in receptor ensemble docking and 38 top-ranked NCI compounds were selected for experimental testing. In [3H]NMS radioligand dissociation assays, about half of the 38 lead compounds altered the radioligand dissociation rate, a hallmark of allosteric behavior. In further competition binding experiments, we identified 12 compounds with affinity ≤30µM. With final functional experiments on six selected compounds, we confirmed four of them as new negative allosteric modulators (NAMs) and one as positive allosteric modulator (PAM) of agonist-mediated response at the M2 mAChR. Two of the NAMs showed subtype selectivity without significant effect at the M1 and M3 mAChRs. To our knowledge, this study demonstrates for the first time an unprecedented successful structure-based approach to identify chemically diverse and selective GPCR allosteric modulators with outstanding potential for further structure-activity relationship studies.
Dynamic structure and inhibition of a malaria drug target: Geranylgeranyl diphosphate synthaseBiochemistry, 55 (36), pp 5180-5190 (2016)
We report a molecular dynamics investigation of the structure, function, and inhibition of geranylgeranyl diphosphate synthase (GGPPS), a potential drug target, from the malaria parasite Plasmodium vivax. We discovered several GGPPS inhibitors, benzoic acids, and determined their structures crystallographically. We then used molecular dynamics simulations to investigate the dynamics of three such inhibitors and two bisphosphonate inhibitors, zoledronate and a lipophilic analogue of zoledronate, as well as the enzyme’s product, GGPP. We were able to identify the main motions that govern substrate binding and product release as well as the molecular features required for GGPPS inhibition by both classes of inhibitor. The results are of broad general interest because they represent the first detailed investigation of the mechanism of action, and inhibition, of an important antimalarial drug target, geranylgeranyl diphosphate synthase, and may help guide the development of other, novel inhibitors as new drug leads.
Striking plasticity of CRISPR-Cas9 and key role of non-target DNA, as revealed by molecular simulationsACS Cent. Sci., (Web) September 9, 2016
The CRISPR (Clustered Regularly Interspaced Short Palindromic Repeats)-Cas9 system recently emerged as a transformative genome-editing technology that is innovating basic bioscience and applied medicine and biotechnology. The endonuclease Cas9 associates to a guide RNA to match and cleave complementary sequences in double stranded DNA, forming an RNA:DNA hybrid and a displaced non-target DNA strand. Although extensive structural studies are ongoing, the conformational dynamics of Cas9 and its interplay with the nucleic acids during association and DNA cleavage are largely unclear. Here, by employing multi-microsecond time scale molecular dynamics, we reveal the conformational plasticity of Cas9 and identify key determinants that allow its large-scale conformational changes during nucleic acid binding and processing. We show how the “closure” of the protein, which accompanies nucleic acid binding, fundamentally relies on highly coupled and specific motions of the protein domains, collectively initiating the prominent conformational changes needed for nucleic acid association. We further reveal a key role of the non-target DNA during the process of activation of the nuclease HNH domain, showing how the non-target DNA positioning triggers local conformational changes that favor the formation of a catalytically competent Cas9. Finally, a remarkable conformational plasticity is identified as an intrinsic property of the HNH domain, constituting a necessary element that allows for the HNH repositioning. These novel findings constitute a reference for future experimental studies aimed at a full characterization of the dynamic features of the CRISPR-Cas9 system, and – more importantly – call for novel structure engineering efforts that are of fundamental importance for the rational design of new genome-engineering applications.
Graded Activation and Free Energy Landscapes of a Muscarinic G protein-coupled ReceptorProc. Natl. Acad. Sci. USA Vol. 113 (43), pp. 12162-12167 (2016)
G protein-coupled receptors (GPCRs) recognize ligands of widely different efficacies, from inverse to partial and full agonists, which transduce cellular signals at differentiated levels. However, the mechanism of such graded activation remains unclear. Using the Gaussian accelerated molecular dynamics (GaMD) method that enables both unconstrained enhanced sampling and free energy calculation, we have performed extensive GaMD simulations (~19 μs in total) to investigate structural dynamics of the M2 muscarinic GPCR that is bound by the full agonist iperoxo (IXO), partial agonist arecoline (ARC) and inverse agonist 3-quinuclidinyl-benzilate (QNB), in the presence or absence of the G protein mimetic nanobody. In the receptor-nanobody complex, IXO binding leads to higher fluctuations in the protein coupling interface than ARC, especially in the receptor TM5, TM6 and TM7 intracellular domains that are essential elements for GPCR activation, but less flexibility in the receptor extracellular region due to stronger binding compared with ARC. Two different binding poses are revealed for ARC in the orthosteric pocket. Removal of the nanobody leads to GPCR deactivation that is characterized by inward movement of the TM6 intracellular end. Distinct low-energy intermediate conformational states are identified for the IXO- and ARC-bound M2 receptor. Both dissociation and binding of an orthosteric ligand are unprecedentedly observed in a single all-atom GPCR simulation in the case of partial agonist ARC binding to the M2 receptor. This study demonstrates the applicability of GaMD for exploring free energy landscapes of large biomolecules and the simulations provide important insights into GPCR functional mechanism.