To request any journal reprints, please contact us.
The Low Dielectric Interior of Proteins is Sufficient to Cause Major Structural Changes in DNA on Association
Adrian H. Elcock and J. Andrew McCammon
Structural changes in DNA resulting from protein association were studied by computer simulation. A model protein defined simply as a region of low dielectric was found to cause a dramatic opening of the minor groove upon close approach to the DNA. Since in these simulations no direct forces acted between the protein and the DNA, the structural changes observed can result only from changes in the solvent and ionic environment of the DNA. From these results it seems likely that a significant driving force for DNA structural changes accompanying complex formation may be provided by the replacement of regions of high dielectric solvent by the low dielectric protein.
Structural Fluctuations of a Cryptophane Host: A Molecular Dynamics Simulation
P.D. Kirchhoff, M.B. Bass, B.A. Hanks, J.M. Briggs, Andre Collet and J.A. McCammon
Cryptophanes are aromatic hosts which bind a variety of guests. Here, we describe a 20 nanosecond molecular dynamics simulation of a particular cryptophane in water. This cryptophane features three pores which open onto a cavity where the guests bind. The molecular dynamics simulation in combination with a surfacing algorithm provides information on the frequency with which these pores open wide enough to admit guest molecules of any given size. We discuss these fluctuations and their possible consequences for binding kinetics.
Surface Titration: A Continuum Electrostatics Model
Jason L. Smart and J. Andrew McCammon
This communication presents a computational study of a highly simplified model of a titrating docosyl amine monolayer at an air/water interface. The specific goal of the study is to reproduce experimentally measured pKa shifts associated with the formation of the monolayer. The calculations also demonstrate the strong dependence of the pKa shift on distance of the amine groups from the interface and the weak dependence of the pKa shift on the size and type of hydrophobic tailgroup. By comparing computation with experiment, an estimate of the amine-interface separation can be made. The dependence of the pKa shift on the ionic strength is calculated at this separation.
Time-Dependent Rate Coefficients from Brownian Dynamics Simulations
Michael J. Potter, Brock Luty, Huan-Xiang Zhou and J. Andrew McCammon
The kinetics of diffusion-influenced reaction is described by a time-dependent rate coefficient k(t). In this paper, we implement an algorithm for calculating k(t), introduced by one of us, in the UHBD software package. In the implementation, the electrostatic force exerted on a diffusing substrate by an enzyme is obtained from a finite-difference solution of the Poisson equation. The rate coefficient is obtained by starting the substrate in the active site and monitoring its survival probability using Brownian dynamics simulations. The technique is applied to the binding of superoxide to Cu/Zn superoxide dismutase. The long-time limit of k(t) is found to be in agreement with both the experimental value and that calculated using an algorithm designed specifically for finding k(t=infinity).
Free Energy Simulations: Correcting for Electrostatic Cutoffs by Use of the Poisson Equation
Haluk Resat and J. Andrew McCammon
The use of electrostatic cutoffs in calculations of free energy changes by molecular dynamics or Monte Carlo simulation is known to introduce errors, which can be quite large when the net charge of the system is changed. The Born equation has often been used to correct for such errors, but this and other analytical methods cannot be used for many systems with complicated structures. Here, we show that numerical methods for solving the Poisson equation, which have been extensively developed recently for studies of solvation thermodynamics, provide a more generally applicable alternative to the traditional Born-type corrections.
The Determinants of pKas in Proteins
Jan Antosiewicz, J. Andrew McCammon and Michael K. Gilson
The formulation of realistic and predictive models for protonation equilibria in proteins represents a fundamental challenge in physical chemistry.
A number of laboratories have used detailed solutions of the linearized Poisson- Boltzmann (PB) equation for proteins in solution to estimate the energetics of ionization processes. Although certain components of the PB model, such as the ionic strength and the dielctric constant of the solvent, are fixed by the experimental conditions to be simulated, other components of the model are less certain. These components include the atomic charges and radii used in the electrostatics calculations; the conformation or conformational distribution of the protein; and the dielectric constant of the protein. The present study makes use of the substantial body of structural and pKa data for the proteins hen egg white lysozyme (HEWL), bovine pancreatic ribonuclease A (RNAse A), bovine pancreatic trypsin inhibitor (PTI), and turkey ovomucoid third domain (OMTKY3), to examine the influence of the dielectric constant of the protein, of atomic parameters, and of protein conformation, upon the accuracy of the calculations. This part of the paper is based upon cumulative comparisons of computed pKas with measurements, and aims to establish models that are increasingly accurate and thus, it is presumed, increasingly realistic.
However, closer examination of some of these systems, and of the experimental data themselves, is also of considerable interest. Accordingly, an additional study of the influence of bound phosphate upon the pKas of histidine side chains in RNAse A is provided, and the basis for the pKa shifts or carboxylic acids in OMTKY3 is analyzed. In addition, observations resulting from a cumulative analysis of 63 measured and calculated pKas are presented.
Quantum-Classical Molecular Dynamics Simulations of Proton Transfer Processes in Molecular Complexes and in Enzymes
P. Bała, P. Grochowski, B. Lesyng and J.A. McCammon
A quantum-classical molecular dynamics model (QCMD) designed for simulations of proton or electron transfer processes in molecular systems is described and applied to several model problems. The primary goal of this work is the elucidation of enzymatic reactions. For example, using the QCMD model, the dynamics of key protons in an enzyme's active site might be described by the time-dependent Schroedinger equation while the dynamics of the remaining atoms are described using MD. The coupling between the quantum proton(s) and the classical atoms is accomplished via extended Hellman-Feynman forces as well as the time dependence of the potential energy function in the Schroedinger equation. The potential energy function is either parameterized prior to the simulations or can be computed using a parameterized valence bond (VB) method (QCMD/VB model). The QCMD method was used to simulate proton transfer in a proton bound ammonia-ammonia dimer as well as to simulate dissociation of a Xe-HI complex in its electronic excited state. The simulation results are compared with data obtained using a quantum-classical time-dependent self-consistent field method (Q/C TDSCF) and with results of fully quantum-dynamical simulations. Finally QCMD/VB simulations of a hydrolytic process catalyzed by phospholipase A(2), including quantum-dynamical dissociation of a water molecule in the active site, are reported. To the best of our knowledge, these are the first simulations that explicitly use the time-dependent Schroedinger equation to describe enzyme catalytic activity.
Orientational Steering in Enzyme-Substrate Association: Ionic Strength Dependence of Hydrodynamic Torque Effects
Jan Antosiewicz, James M. Briggs and J. Andrew McCammon
The effect of hydrodynamic torques on the association rate constants for enzyme-ligand complexation is investigated by Brownian dynamics simulations. Our hydrodynamic models of the enzyme and ligand are composed of spherical elements with friction forces acting at their centers. A quantitative measure of hydrodynamic torque orientational effects is introduced by choosing, as a reference system, an enzyme-ligand model with the same average hydrodynamic interactions but without orientational dependence. Our simple models show a 15% increase in the rate constant caused by hydrodynamic torques at physiological ionic strength. For more realistic hydrodynamic models, which are not computationally feasible at present, this effect is probably higher. The most important finding of this work is that hydrodynamic complementarity in shape (i.e. like the fitting together of pieces of a puzzle) is most effective for interactions between molecules at physiological ionic strength.
Comparison of Continuum and Explicit Models of Solvation: Potentials of Mean Force for Alanine Dipeptide
Tami J. Marrone, Michael K. Gilson and J. Andrew McCammon
We compute the potential of mean force (PMF) around the phi and psi torsions of alanine dipeptide with a Poisson-Boltzmann (PB) method and compare these results to simulations in explicit water. The PB methods, which includes an apolar solvation term, qualitatively reproduces the PMF profiles generated in the explicit solvent simulation at a markedly lower computational cost. These results motivate more extensive testing of continuum methods for the study of conformational and binding equilibria in solution.
Binding of Tacrine and 6-Chlorotacrine by Acetylcholinesterase
S.T. Wlodek, J. Antosiewicz, J.A. McCammon, T.P. Straatsma, M.K. Gilson, J.M. Briggs, C. Humblet and J.L. Sussman
Multiconfiguration thermodynamic integration was used to determine the relative binding strength of tacrine and 6-chlorotacrine by Torpedo californica acetylcholinesterase. 6-Chlorotacrine appears to be bound stronger by 0.7+/-0.4 kcal/mol than unsubstituted tacrine when the active site triad residue His-440 is deprotonated. This result is in excellent agreement with experimental inhibition data on electric eel acetylcholinesterase. Electrostatic Poisson-Boltzmann calculations confirm that order of binding strength, resulting in ΔG of binding of -2.9 and -3.3 kcal/mol for tacrine and chlorotacrine, respectively, and suggest inhibitor binding does not occur when His-440 is charged. Our results suggest that electron density redistribution upon tacrine chlorination is mainly responsible for the increased attraction potential between protonated inhibitor molecule and adjacent aromatic groups of Phe-330 and Trp-84.
Use of the Grand Canonical Ensemble in Potential of Mean Force Calculations
H. Resat, M. Mezei and J.A. McCammon
Understanding and predicting the thermodynamics of association reactions at the microscopic level requires that it be possible to sample respresentative configurations of the reactants and solvent as a function of the reaction pathways. Because of geometric effects, certain methodological improvements in molecular simulation techniques are necessary before the reaction thermodynamics of complicated systems such as biopolymers with interlocking shapes can be investigated. Here, we propose the use of the grand canonical ensemble in molecular simulations when the traditional canonical ensemble based methods cannot appropriately account for the confined space effects. The success of the grand canonical ensemble molecular simulations in studying the association reaction profile is shown by testing it on simpler systems. Implications for future work and various possible application areas of the grand canonical ensemble simulations are discussed.
Study of Global Motions in Proteins by Weighted Masses Molecular Dynamics: Adenylate Kinase as a Test Case
Samir Elamrani, Michael B. Berry, George N. Phillips, Jr. and J. Andrew McCammon
The weighted masses molecular dynamics (WMMD) technique is applied to the protein adenylate kinase. A novel set of restraints has been developed to allow the use of this technique with proteins. The WMMD simulation is successful in predicting the flexibility of the two mobile domains of the protein. The end product of the simulation is similar to the known open and AMP bound forms of the enzyme. The biological relevance of the restraints used and potential methods of improving the technique are discussed.
Theoretical Study of Inhibition of Adenosine Deaminase by 8R-Coformycin and 8R-Deoxycoformycin
Tami J. Marrone, T.P. Straatsma, James M. Briggs, David K. Wilson, Florante A. Quiocho and J. Andrew McCammon
Molecular dynamics and free energy simulations were performed to examine the binding of (8R)-deoxycoformycin and (8R)-coformycin to adenosine deaminase. The two inhibitors differ only at the 2' position of the sugar ring; the sugar moiety of conformycin is ribose, while it is deoxyribose for deoxycoformycin. The 100 ps molecular dynamics trajectories reveal that Asp 19 and His 17 interact strongly with the 5' hydroxyl group of the sugar moiety of both inhibitors and appear to play an important role in binding the sugar. The 2' and 3' groups of the sugars are near the protein-water interface and can be stabilized by either protein residues or water. The flexibility of the residues at the opening of the active site helps to explain the modest difference in binding of the two inhibitors and how substrates/inhibitors can enter an otherwise inaccessible binding site.
Density Functional Based Parametrization of a Valence Bond Method and its Applications in Quantum-Classical Molecular Dynamics Simulations of Enzymatic Reactions
P. Grochowski, B. Lesyng, P. Bała and J.A. McCammon
An approximate valence bond (AVB) method was parametrized at a microscopic level for proton transfer and hydroxyanion nucleophilic reactions in enzyme catalytic processes. The method was applied to describe hydrolytic activity of phospholipase A2. The AVB parametrization is based on density functional and conventional ab initio calculations calibrated with respect to experimental data in the gas phase. The method was used as a fast generator of the potential energy function in a quantum-classical molecular dynamics (QCMD) simulations describing atomic motions as well as propagation of the proton wave function in the enzyme active site. The protein environment surrounding the active site and solvent effects are included in the model via electrostatic interactions perturbing the original AVB Hamiltonian.
Electrostatic Channeling in the Bifunctional Enzyme Dihydrofolate Reductase-Thymidylate Synthase
Adrian H. Elcock, Michael J. Potter, David A. Matthews, Daniel R. Knighton and J. Andrew McCammon
The bifunctional enzyme dihydrofolate reductase-thymidylate synthase (DHFR-TS) carries out two distinct reactions, with the dihydrofolate produced by the TS-catalysed reaction acting as the substrate for the DHFR-catalysed reaction. Brownian dynamics simulation techniques were used to investigate the possible role of electrostatics in determining efficient channeling of the substrate, by explicitly simulating substrate diffusion between the two active sites. With a substrate charge of -2, almost all (>95%) substrate molecules leaving the TS active site reached the DHFR active site at zero ionic strength. Under the same conditions, but in the absence of electrostatic effects, successful channeling was reduced to only around 6%: electrostatic effects therefore appear essential to explain the efficient channeling observed experimentally. The importance of substrate charge, the relative contributions of specific basic residues in the protein, the role played by the second monomer of the dimer in channeling and the effects of changing ionic strength were all investigated. Simulations performed for substrate transfer in the opposite direction suggest that channeling in DHFR-TS is not strongly directional and that the role of electrostatics is perhaps more one of restricting diffusion of the substrate than one of actively guiding it from the TS to the DHFR active site. The results demonstrate that electrostatic channeling can be a highly efficient means of transferring charged substrates between active sites in solvent-exposed environments.
Computational Science: New Horizons and Relevance to Pharmaceutical Design
James M. Briggs, Tami J. Marrone and J. Andrew McCammon
Computer methods are used extensively in the design and refinement of drug leads. A short summary is given for several computational methods followed by a description of how some of these methods have been applied to design drugs targeted to the renin-angiotensin system and to cholinergic synapses. These methods include quantitative structure activity relationship (QSAR) methods, comparative molecular field analyses (CoMFA), 3D database searching, de novo design of ligands, docking, and computational alchemy (free energy perturbation, FEP, thermodynamic integration, MCTI). Most of these methods can be used whether or not detailed structural information about the binding site is available, although without an X-ray structure, the analyses are more qualitative. All of these methods are used extensively in the commercial design of pharmaceuticals. The main problem with most of these mthods is in the scoring (ranking) of interactions or matches. Advances in this area and others (methods development and increase in capabilities of computers) will increase the predictive power of these mthods and help to speed the time to market of new pharmaceuticals.
A Speed Limit for Protein Folding (Commentary)
J. Andrew McCammon
Recent observations of fast steps in protein folding address the question: What is the shortest time that would be required for folding? The paper by Hagen et al. in the current issue of the Proceedings provides strong experimental support for a speed limit that is governed by diffusion. Diffusion is known to limit the rates of many of the activities of proteins, including enzymatic activity, electron transfer, and binding to other macromolecules. The new evidence suggests that folding can be added to this list.
Evidence for Electrostatic Channeling in a Fusion Protein of Malate Dehydrogenase and Citrate Synthase
Adrian H. Elcock and J. Andrew McCammon
Brownian dynamics simulations were performed to investigate a possible role for electrostatic channeling in transferring substrate between two of the enzymes of the citric acid cycle. The diffusion of oxaloacetate from one of the active sites of malate dehydrogenase (MDH) to the active sites of citrate synthase (CS) was simulated in the presence and absence of electrostatic forces using a modeled structure for a MDH-CS fusion protein. In the absence of electrostatic forces, fewer than 1% of substrate molecules leaving the MDH active site are transferred to CS. When electrostatic forces are present at zero ionic stength however, around 45% of substrate molecules are successfully channeled. As expected for an electrostatic mechanism of transfer, increasing the ionic strength in the simulations reduces the calculated transfer efficiency. Even at 150mM however, the inclusion of electrostatic forces results in an increase in transfer efficiency of more than one order of magnitude. The simulations therefore provide evidence for the involvement of electrostatic channeling in guiding substrate transfer between two of the enzymes of the citric acid cycle. Similar effects may opeate between other members of the citric acid metabolon.
Acetylcholinesterase: Role of the Enzyme's Charge Distribution in Steering Charged Ligands Toward the Active Site
Jan Antosiewicz, Stanislaw T. Wlodek and J. Andrew McCammon
The electrostatic steering of charged ligands toward the active site of Torpedo californica acetylcholinesterase is investigated by Brownian dynamics simulations of wild type enzyme and several mutated forms, in which some normally charged residues are neutralized.
The simulations reveal that the total ligand influx through a surface of 42A radius centered in the enzyme monomer and separated from the protein surface by 1-14A is not significantly influenced by electrostatic interactions. Electrostatic effects are visible for encounters with a surface of 32A radius, which is partially hidden inside the protein, but mostly within the solvent. A clear accumulation of encounter events for that sphere is observed in the area directly above the entrance to the active site gorge. In this area, the encounter events are increased by 40% compared to the case of a netural ligand, However, the differences among the encounter rates for the various mutants considered here are not pronouned, all rate constants being within plus or minus 10% of the average value.
The enzyme charge distribution becomes more important as the charged ligand moves toward the bottom of the gorge, where the active site is located. We show that neither the enzyme's total charge, nor its dipole moment, fully account for the electrostatic steering of ligand to the active site. Higher moments of the enzyme's charge distribution are also important. However, for a series of mutations for which the direction of the enzyme dipole moment is constant within a few degrees, one observes a gradual decrease in the diffusional encounter rate constant with the number of neutralized residues. On the other hand, for other mutants that change the direction of the dipole moment from that of the wild type, the calculated encounter rate constants can be very close to that of the wild type.
The present work yields two new insights to the kinetics of acetylcholinesterase. First, evolution appears to have built a redundant electrostatic steering capability into this important enzyme through the overall distribution of its thousands of partially charged atoms. And second, roughly half of the rate enhancement due to electrostatics arises from steering of the substrate outside the enzyme; the other half of the rate enhancement arises from improved trapping of the substrate after it had entered the gorge. The computational results reproduce qualitatively, and help to rationalize, many surprising experimental results obtained recently for human acetylcholinesterase.
Computing the Ionization States of Proteins with a Detailed Charge Model
Jan Antosiewicz, James M. Briggs, Adrian H. Elcock, Michael K. Gilson and J. Andrew McCammon
A convenient computational approach for the calculation of the pKas of ionizable groups in a protein is described. The method uses detailed models of the charges in both the neutral and ionized form of each ionizable group. A full derivation of the theoretical framework is presented, as are details of its implementation in the UHBD program. Application to four proteins whose crystal structures are known shows that the detailed charge model improves agreement with experimentally determined pKas when a low protein dielectric constant is assumed, relative to the results with a simpler single-site ionization model. It is also found that use of the detailed charge model increases the sensitivity of the computed pKas to the details of proton placement.
A 240-Fold Electrostatic Rate-Enhancement for Acetylcholinesterase - Substrate Binding can be Predicted by the Potential within the Active Site
Huan-Xiang Zhou, James M. Briggs and J. Andrew McCammon
An approximate relation that allows easy prediction of the enhancement of enzyme-substrate binding rates by an interaction potential was derived recently by one of us (Zhou, J. Chem. Phys., 1996, 105, 7235). This is given by the average of the Boltzmann factor over the region in which the substrate can effectively bind to the enzyme. We test this relation on the enzyme Torpedo californica acetylcholinesterase. The interaction potential was found by solving the Poisson-Boltzmann equation and treating the substrate as a test charge. At an ionic strength of 150 mM, the electrostatic rate-enhancement as calculated by Brownian dynamics simulation is enormous (~240 fold), yet it can be predicted by the average Boltzmann factor (~1000) to within a factor of 4. The simulated binding rate constants are in reasonable agreement with experimental results.