To request any journal reprints, please contact us.
Mechanistic insight into the role of transition state stabilization in Cyclophilin A
Donald Hamelberg and J. Andrew McCammon
Peptidyl prolyl cis–trans isomerases (PPIases) are ubiquitous enzymes in biology that catalyze the cis–trans isomerization of the proline imide peptide bond in many cell signaling pathways. The local change of the isomeric state of the prolyl peptide bond acts as a switching mechanism in altering the conformation of proteins. A complete understanding of the mechanism of PPIases is still lacking, and current experimental techniques have not been able to provide a detailed atomistic picture. Here we have carried out several accelerated molecular dynamics simulations with explicit solvent, and we have provided a detailed description of cis–trans isomerization of the free and cyclophilin A-catalyzed process. We show that the catalytic mechanism of cyclophilin is due mainly to the stabilization and preferential binding of the transition state that is achieved by a favorable hydrogen bond interaction with a backbone NH group. We also show that the substrate in the transition state interacts more favorably with the enzyme than the cis isomer, which in turn interacts more favorably than the trans isomer. The stability of the enzyme–substrate complex is directly correlated with the interaction the substrate makes with a highly conserved arginine residue. Finally, we show that catalysis is achieved through the rotation of the carbonyl oxygen on the N-terminal of the prolyl peptide bond in a predominately unidirectional fashion.
Multi-Scale Modeling of Ventricular Myocytes: Contributions of structural and functional heterogeneities to excitation-contraction coupling in the normal and failing rodent heart
Shaoying Lu, Anushka P. Michailova, Jeffrey J. Saucerman, Yuhui Cheng, Zeyun Yu, Timothy H. Kaiser, Wilfred W. Li, Randolph E. Bank, Michael J. Holst, J. Andrew McCammon, Takeharu Hayashi, Masahiko Hoshijima, Peter Arzberger, Andrew D. McCulloch
There is a growing body of experimental evidence suggesting that the Ca2+ signaling in ventricular myocytes is characterized by high gradient near the cell membrane and a more uniform Ca2+ distribution in the cell interior. An important reason for this phenomenon might be that in these cells the t-tubular system forms a network of extracellular space, extending deep into the cell interior. This allows the electrical signal, that propagates rapidly along the cell membrane, to reach the vicinity of the sarcoplasmic reticulum (SR) where intracellular Ca2+ required for myofilament activation is stored. Early studies of cardiac muscle showed that the t-tubules are found at intervals of ~2 μm along the longitudinal cell axis in close proximity to the Z-disks of the sarcomeres. Subsequent studies have demonstrated that the t-tubular system has also longitudinal extensions.
Nathan A. Baker and J. Andrew McCammon
Computational Studies of Protein Dynamics
J. Andrew McCammon
Theoretical and computational studies of protein function have reached the point at which they are making important contributions to drug discovery and other practical applications. At the same time, they are deepening our understanding of the principles of protein activity, including the dynamical features that give rise to NMR and other experimental measurements, and the time-dependent aspects of biological function.
Conformational Dynamics of the Flexible Catalytic Loop in Mycobacterium tuberculosis 1-Deoxy-D-xylulose 5-Phosphate Reductoisomerase
Sarah L. Williams and J. Andrew McCammon
In mycobacteria, the biosynthesis of the precursors to the essential isoprenoids, isopentenyl diphosphate and dimethylallyl pyrophosphate is carried out by the methylerythritol phosphate pathway. This route of synthesis is absent in humans, who utilize the alternative mevalonate acid route, thus making the enzymes of the methylerythritol phosphate pathway of chemotherapeutic interest. One such identified target is the second enzyme of the pathway, 1-deoxy-D-xylulose 5-phosphate reductoisomerase. Only limited information is currently available concerning the catalytic mechanism and structural dynamics of this enzyme, and only recently has a crystal structure of Mycobacterium tuberculosis species of this enzyme been resolved including all factors required for binding. Here, the dynamics of the enzyme is studied in complex with NADPH, Mn2+, in the presence and absence of the fosmidomycin inhibitor using conventional molecular dynamics and an enhanced sampling technique, reversible digitally filtered molecular dynamics. The simulations reveal significant differences in the conformational dynamics of the vital catalytic loop between the inhibitor-free and inhibitor-bound enzyme complexes and highlight the contributions of conserved residues in this region. The substantial fluctuations observed suggest that 1-deoxy-D-xylulose 5-phosphate reductoisomerase may be a promising target for computer-aided drug discovery through the relaxed complex method.
Coupling the Level-Set Method with Molecular Mechanics for Variational Implicit Solvation of Nonpolar Molecules
Li-Tien Cheng, Yang Xie, Joachim Dzubiella, J. Andrew McCammon, Jianwei Che, Bo Li
We construct a variational explicit-solute implicit-solvent model for the solvation of molecules. Central in this model is an effective solvation free-energy functional that depends solely on the position of solute-solvent interface and solute atoms. The total free energy couples altogether the volume and interface energies of solutes, the solute-solvent van der Waals interactions, and the solute-solute mechanical interactions. A curvature dependent surface tension is incorporated through the so-called Tolman length which serves as the only fitting parameter in the model. Our approach extends the original variational implicit-solvent model of Dzubiella, Swanson, and McCammon [Phys. Rev. Lett. 2006, 96, 087802 and J. Chem. Phys. 2006, 124, 084905] to include the solute molecular mechanics. We also develop a novel computational method that combines the level-set technique with optimization algorithms to determine numerically the equilibrium conformation of nonpolar molecules. Numerical results demonstrate that our new model and methods can capture essential properties of nonpolar molecules and their interactions with the solvent. In particular, with a suitable choice of the Tolman length for the curvature correction to the surface tension, we obtain the solvation free energy for a benzene molecule in a good agreement with experimental results.
Finite Element Analysis of Drug Electrostatic Diffusion: Inhibition Rate Studies in N1 Neuraminidase
Yuhui Cheng, Michael J. Holst, J.A. McCammon
This article describes the numerical solution of the steady-state Poisson-Boltzmann-Smoluchowski (PBS) and Poisson-Nernst-Planck (PNP) equations to study diffusion in biomolecular systems. Specifically, finite element methods have been developed to calculate electrostatic interactions and ligand binding rate constants for large biomolecules. The resulting software has been validated and applied to the wild-type and several mutated avian influenza neurominidase crystal structures. The calculated rates show very good agreement with recent experimental studies. Furthermore, these finite element methods require significantly fewer computational resources than existing particle-based Brownian dynamics methods and are robust for complicated geometries. The key nding of biological importance is that the electrostatic steering plays the important role in the drug binding process of the neurominidase.
Toward Understanding the Conformational Dynamics of RNA Ligation
Robert V. Swift, Jacob Durrant, Rommie E. Amaro, J. Andrew McCammon
Members of the genus Trypanosoma, which include the pathogenic species Trypanosoma brucei and Trypanosoma cruzi, edit their post-transcriptional mitochondrial RNA via a multiprotein complex called the editosome. In T. brucei, the RNA is nicked prior to uridylate insertion and deletion. Following editing, nicked RNA is religated by one of two RNA-editing ligases (TbREL). This study describes a recent 70 ns molecular dynamics simulation of TbREL1, an ATP-dependent RNA-editing ligase of the nucleotidyltransferase superfamily that is required for the survival of T. brucei insect and bloodstream forms. In this work, a model of TbREL1 in complex with its full double-stranded RNA (dsRNA) substrate is created on the basis of the homologous relation between TbREL1 and T4 Rnl2. The simulation captures TbREL1 dynamics in the state immediately preceding RNA ligation, providing insights into the functional dynamics and catalytic mechanism of the kinetoplastid ligation reaction. Important features of RNA binding and specificity are revealed for kinetoplastid ligases and the broader nucleotidyltransferase superfamily.
AutoGrow: A Novel Algorithm for Protein Inhibitor Design
Jacob D. Durrant, Rommie E. Amaro, J. Andrew McCammon
Due in part to the increasing availability of crystallographic protein structures as well as rapid improvements in computing power, the past few decades have seen an explosion in the field of computer-based rational drug design. Several algorithms have been developed to identify or generate potential ligands in silico by optimizing the ligand–receptor hydrogen bond, electrostatic, and hydrophobic interactions. We here present AutoGrow, a novel computer-aided drug design algorithm that combines the strengths of both fragment-based growing and docking algorithms. To validate AutoGrow, we recreate three crystallographically resolved ligands from their constituent fragments.
Characterizing Loop Dynamics and Ligand Recognition in Human- and Avian-type Influenza Neuraminidases via Generalized Born Molecular Dynamics and End-point Free Energy Calculations
Rommie E. Amaro, Xiaolin Cheng, Ivaylo Ivanov, Dong Xu, J. Andrew McCammon
The comparative dynamics and inhibitor binding free energies of group-1 and group-2 pathogenic influenza A subtype neuraminidase (NA) enzymes are of fundamental biological interest and relevant to structure-based drug design studies for antiviral compounds. In this work, we present seven generalized Born molecular dynamics simulations of avian (N1)- and human (N9)-type NAs in order to probe the comparative flexibility of the two subtypes, both with and without the inhibitor oseltamivir bound. The enhanced sampling obtained through the implicit solvent treatment suggests several provocative insights into the dynamics of the two subtypes, including that the group-2 enzymes may exhibit similar motion in the 430-binding site regions but different 150-loop motion. End-point free energy calculations elucidate the contributions to inhibitor binding free energies and suggest that entropic considerations cannot be neglected when comparing across the subtypes. We anticipate the findings presented here will have broad implications for the development of novel antiviral compounds against both seasonal and pandemic influenza strains.
MM-PBSA captures key role of intercalating water at a protein-protein interface
Sergio Wong, Rommie E. Amaro, J. Andrew McCammon
The calculation of protein interaction energetics is of fundamental interest, yet accurate quantities are difficult to obtain due to the complex and dynamic nature of protein interfaces. This is further complicated by the presence of water molecules, which can exhibit transient interactions of variable duration and strength with the protein surface. The T-cell receptor (TCR) and its staphylococcal enterotoxin 3 (SEC3) binding partner are well-characterized examples of a protein-protein interaction system exhibiting interfacial plasticity, cooperativity, and additivity among mutants. Specifically engineered mutants induce intercalating interfacial water molecules, which subsequently enhance protein-protein binding affinity. In this work, we perform a set of molecular mechanics (MM) Poisson-Boltzmann (PB) surface area (SA) calculations on the wild type and two mutant TCR-SEC3 systems and show that the method is able to discriminate between weak and strong binders only when key explicit water molecules are included in the analysis. The results presented here point to the promise of MM-PBSA toward rationalizing molecular recognition at protein-protein interfaces, while establishing a general approach to handle explicit interfacial water molecules in such calculations.
Free Energy for the Permeation of Na+ and Cl- Ions and their Ion-Pair through Zwitterionic Dimyristoyl Phosphatidylcholine Lipid Bilayer by Umbrella Integration with Harmonic Fourier Beads
Ilja V. Khavrutskii, Alemayehu A. Gorfe, Benzhuo Lu, J. Andrew McCammon
Understanding the mechanism of ion permeation across lipid bilayers is key to controlling osmotic pressure and developing new ways of delivering charged, drug-like molecules inside cells. Recent reports suggest ion-pairing as the mechanism to lower the free energy barrier for the ion permeation in disagreement with predictions from the simple electrostatic models. In this paper we quantify the effect of ion-pairing or charge quenching on the permeation of Na+ and Cl- ions across DMPC lipid bilayer by computing the corresponding potentials of mean force (PMFs) using fully atomistic molecular dynamics simulations. We find that the free energy barrier to permeation reduces in the order Na+–Cl- ion-pair (27.6 kcal/mol) > Cl- (23.6 kcal/mol) > Na+ (21.9 kcal/mol). Furthermore, with the help of these PMFs we derive the change in the binding free energy between the Na+ and Cl- with respect to that in water as a function of the bilayer permeation depth. Despite the fact that the bilayer boosts the Na+–Cl- ion binding free energy by as high as 17.9 kcal/mol near its center, ion-pairing between such hydrophilic ions as Na+ and Cl- does not assist their permeation. However, based on a simple thermodynamic cycle, we suggest that ion-pairing between ions of opposite charge and solvent philicity could enhance ion permeation. Comparison of the computed permeation barriers for Na+ and Cl- ions with available experimental data supports this notion. This work establishes general computational methodology to address ion-pairing in fluid anisotropic media and details the ion permeation mechanism on atomic level.
A virtual screening study of the acetylcholine binding protein using a relaxed-complex approach
Arneh Babakhani, Todd T. Talley, Palmer Taylor, J. Andrew McCammon
The nicotinic acetylcholine receptor (nAChR) is a member of the ligand-gated ion channel family and is implicated in many neurological events. Yet, the receptor is difficult to target without high-resolution structures. In contrast, the structure of the acetylcholine binding protein (AChBP) has been solved to high resolution, and it serves as a surrogate structure of the extra-cellular domain in nAChR. Here we conduct a virtual screening study of the AChBP using the relaxed-complex method, which involves a combination of molecular dynamics simulations (to achieve receptor structures) and ligand docking. The library screened through comes from the National Cancer Institute, and its ligands show great potential for binding AChBP in various manners. These ligands mimic the known binders of AChBP; a significant subset docks well against all species of the protein and some distinguish between the various structures. These novel ligands could serve as potential pharmaceuticals in the AChBP/nAChR systems.
Single-channel current through nicotinic receptor produced by closure of binding site C-loop
Hai-Long Wang, Reza Toghraee, David Papke, Xiao-Lin Cheng, J. Andrew McCammon, Umberto Ravaioli, Steven M. Sine
We investigated the initial coupling of agonist binding to channel gating of the nicotinic acetylcholine receptor using targeted molecular-dynamics (TMD) simulation. After TMD simulation to accelerate closure of the C-loops at the agonist binding sites, the region of the pore that passes through the cell membrane expands. To determine whether the structural changes in the pore result in ion conduction, we used a coarse-grained ion conduction simulator, Biology Boltzmann transport Monte Carlo, and applied it to two structural frames taken before and after TMD simulation. The structural model before TMD simulation represents the channel in the proposed "resting" state, whereas the model after TMD simulation represents the channel in the proposed "active" state. Under external voltage biases, the channel in the "active" state was permeable to cations. Our simulated ion conductance approaches that obtained experimentally and recapitulates several functional properties characteristic of the nicotinic acetylcholine receptor. Thus, closure of the C-loop triggers a structural change in the channel sufficient to account for the open channel current. This approach of applying Biology Boltzmann transport Monte Carlo simulation can be used to further investigate the binding to gating transduction mechanism and the structural bases for ion selection and translocation.
Thioamide Hydroxypyrothiones Supersede Amide Hydroxypyrothiones in Potency Against Anthrax Lethal Factor
Arpita Agrawal, César Augusto F. de Oliveira, Yuhui Cheng, Jennifer A. Jacobsen, J. Andrew McCammon, Seth M. Cohen
Anthrax lethal factor (LF) is a critical virulence factor in the pathogenesis of anthrax. A structure-activity relationship (SAR) of potential lethal factor inhibitors (LFi) is presented in which the zinc-binding group (ZBG), linker, and backbone moieties for a series of hydroxypyrone-based compounds were systematically varied. It was found that hydroxypyrothione ZBGs generate more potent inhibitors than hydroxypyrone ZBGs. Furthermore, coupling the hydroxypyrothione to a backbone group via a thioamide bond improves potency when compared to an amide linker. QM/MM studies show that the thioamide bond in these inhibitors allows for the formation of two additional hydrogen bonds with the protein active site. In both types of hydroxypyrothione compounds, ligand efficiencies of 0.29–0.54 kcal/mol per heavy atom were achieved. The results highlight the need for a better understanding to optimize the interplay between the ZBG, linker, and backbone to get improved LFi.
Distinct Glycan Topology for Avian and Human Sialo-Pentasaccharide Receptor Analogues upon Binding Different Hemagglutinins: A Molecular Dynamics Perspective
Dong Xu, E. Irene Newhouse, Rommie E. Amaro, Hsing C. Pao, Lily S. Cheng, Phineus R.L. Markwick, J. Andrew McCammon, Wilfred W. Li, Peter W. Arzberger
Hemagglutinin (HA) binds to sialylated glycans exposed on the host cell surface in the initial stage of avian influenza virus infection. It has been previously hypothesized that glycan topology plays a critical role in the human adaptation of avian flu viruses, such as the potentially pandemic H5N1. Comparative molecular dynamics studies are complementary to experimental techniques, including glycan microarray, to understand the mechanism of species-specificity switch better. The examined systems comprise explicitly solvated trimeric forms of avian H3, H5, and swine H9 in complex with avian and human glycan receptor analogues—LSTa (α-2,3-linked lactoseries tetrasaccharide a) and LSTc (α-2,6-linked lactoseries tetrasaccharide c), respectively. The glycans adopted distinct topological profiles with inducible torsional angles when bound to different HAs. The corresponding receptor binding domain amino acid contact profiles were also distinct. Avian H5 was able to accommodate LSTc in a tightly “folded umbrella”-like topology through interactions with all five sugar residues. After considering conformational entropy, the relative binding free-energy changes, calculated using the molecular mechanics-generalized Born surface area technique, were in agreement with previous experimental findings and provided insights on electrostatic, van der Waals, desolvation, and entropic contributions to HA–glycan interactions. The topology profile and the relative abundance of free glycan receptors may influence receptor binding kinetics. Glycan composition and topological changes upon binding different HAs may be important determinants in species-specificity switch.
Ras Conformational Switching: Simulating Nucleotide Dependent Conformational Transitions with Accelerated Molecular Dynamics
Barry J. Grant, Alemayehu A. Gorfe, J. Andrew McCammon
Ras mediates signaling pathways controlling cell proliferation and development by cycling between GTP- and GDP-bound active and inactive conformational states. Understanding the complete reaction path of this conformational change and its intermediary structures is critical to understanding Ras signaling. We characterize nucleotide-dependent conformational transition using multiple-barrier-crossing accelerated molecular dynamics (aMD) simulations. These transitions, achieved for the first time for wild-type Ras, are impossible to observe with classical molecular dynamics (cMD) simulations due to the large energetic barrier between end states. Mapping the reaction path onto a conformer plot describing the distribution of the crystallographic structures enabled identification of highly populated intermediate structures. These structures have unique switch orientations (residues 25–40 and 57–75) intermediate between GTP and GDP states, or distinct loop3 (46–49), loop7 (105–110), and α5 C-terminus (159–166) conformations distal from the nucleotide-binding site. In addition, these barrier-crossing trajectories predict novel nucleotide-dependent correlated motions, including correlations of α2 (residues 66–74) with α3-loop7 (93–110), loop2 (26–37) with loop10 (145–151), and loop3 (46–49) with α5 (152–167). The interconversion between newly identified Ras conformations revealed by this study advances our mechanistic understanding of Ras function. In addition, the pattern of correlated motions provides new evidence for a dynamic linkage between the nucleotide-binding site and the membrane interacting C-terminus critical for the signaling function of Ras. Furthermore, normal mode analysis indicates that the dominant collective motion that occurs during nucleotide-dependent conformational exchange, and captured in aMD (but absent in cMD) simulations, is a low-frequency motion intrinsic to the structure.
Independent-Trajectories Thermodynamic-Integration Free-Energy Changes for Biomolecular Systems: Determinants of H5N1 Avian Influenza Virus Neuraminidase Inhibition by Peramivir
Morgan Lawrenz, Riccardo Baron, J. Andrew McCammon
Free-energy changes are essential physicochemical quantities for understanding most biochemical processes. Yet, the application of accurate thermodynamic-integration (TI) computation to biological and macromolecular systems is limited by finite-sampling artifacts. In this paper, we employ independent-trajectories thermodynamic-integration (IT-TI) computation to estimate improved free-energy changes and their uncertainties for (bio)molecular systems. IT-TI aids sampling statistics of the thermodynamic macrostates for flexible associating partners by ensemble averaging of multiple, independent simulation trajectories. We study peramivir (PVR) inhibition of the H5N1 avian influenza virus neuraminidase flexible receptor (N1). Binding site loops 150 and 119 are highly mobile, as revealed by N1-PVR 20-ns molecular dynamics. Due to such heterogeneous sampling, standard TI binding free-energy estimates span a rather large free-energy range, from a 19% underestimation to a 29% overestimation of the experimental reference value (-62.2 ± 1.8 kJ/mol). Remarkably, our IT-TI binding free-energy estimate (-61.1 ± 5.4 kJ/mol) agrees with a 2% relative difference. In addition, IT-TI runs provide a statistics-based free-energy uncertainty for the process of interest. Using ~800 ns of overall sampling, we investigate N1-PVR binding determinants by IT-TI alchemical modifications of PVR moieties. These results emphasize the dominant electrostatic contribution, particularly through the N1 E277–PVR guanidinium interaction. Future drug development may be also guided by properly tuning ligand flexibility and hydrophobicity. IT-TI will allow estimation of relative free energies for systems of increasing size, with improved reliability by employing large-scale distributed computing.
Substrate Induced Population Shifts and Stochastic Gating in the PBCV-1 mRNA Capping Enzyme
Robert V. Swift and J. Andrew McCammon
The 317 residue PBCV-1 mRNA capping enzyme catalyzes the second enzymatic reaction in the formation of an N-7-methyl-GMP cap on the 5'-end of the nascent mRNA. It is composed of two globular domains bound by a short flexible peptide linker, which have been shown to undergo opening and closing events. The small size and experimentally demonstrated domain mobility make the PBCV-1 capping enzyme an ideally suited model system to explore domain mobility in context of substrate binding. Here, we specifically address the following four questions: (1) How does substrate binding affect relative domain mobility: is the system better described by an induced fit or population shift mechanism? (2) What are the gross characteristics of a conformation capable of binding substrate? (3) Does “domain gating” of the active site affect the rate of substrate binding? (4) Does the magnitude of receptor conformational fluctuations confer substrate specificity by sterically occluding molecules of a particular size or geometry? We answer these questions using a combination of theory, Brownian dynamics, and molecular dynamics. Our results show that binding efficiency is a function of conformation but that isomerization between efficient and inefficient binding conformations does not impact the substrate association rate. Additionally, we show that conformational flexibility alone is insufficient to explain single stranded mRNA specificity. While our results are specific to the PBCV-1 mRNA capping enzyme, they provide a useful context within which the substrate binding behavior of similarly structured enzymes or proteins may be considered.
Molecular Dynamics Simulations of ELIC - a Prokaryotic Homologue of the Nicotinic Acetylcholine Receptor
Xiaolin Cheng, Ivaylo Ivanov, Hailong Wang, Steven M. Sine, J. Andrew McCammon
The ligand-gated ion channel from Erwinia chrysanthemi (ELIC) is a prokaryotic homologue of the eukaryotic nicotinic acetylcholine receptor (nAChR) that responds to the binding of neurotransmitter acetylcholine and mediates fast signal transmission. ELIC is similar to the nAChR in its primary sequence and overall subunit organization, but despite their structural similarity, it is not clear whether these two ligand-gated ion channels (LGICs) operate in a similar manner. Further, it is not known to what extent mechanistic insights gleaned from the ELIC structure translate to eukaryotic counterparts such as the nAChR. Here we use molecular dynamics simulation to probe the conformational dynamics and hydration of the transmembrane pore of ELIC. The results are compared with those from our previous simulation of the human α7 nAChR. Overall, ELIC displays increased stability compared to the nAChR, whereas the two proteins exhibit remarkable similarity in their global motion and flexibility patterns. The majority of the increased stability of ELIC does not stem from the deficiency of the models used in the simulations, and but rather seems to have a structural basis. Slightly altered dynamical correlation features are also observed among several loops within the membrane region. In sharp contrast to the nAChR, ELIC is completely dehydrated from the pore center to the extracellular end throughout the simulation. Finally, the simulation of an ELIC mutant substantiates the important role of F246 on the stability, hydration and possibly function of the ELIC channel.
Using multistate free energy techniques to improve the efficiency of replica exchange accelerated molecular dynamics
Mikolai Fajer, Robert V. Swift, J. Andrew McCammon
Replica exchange accelerated molecular dynamics (REXAMD) is a method that enhances conformational sampling while retaining at least one replica on the original potential, thus avoiding the statistical problems of exponential reweighting. In this article, we study three methods that can combine the data from the accelerated replicas to enhance the estimate of properties on the original potential: weighted histogram analysis method (WHAM), pairwise multistate Bennett acceptance ratio (PBAR), and multistate Bennett acceptance ratio (MBAR). We show that the method that makes the most efficient use of equilibrium data from REXAMD simulations is the MBAR method. This observation holds for both alchemical free energy and structural observable prediction. The combination of REXAMD and MBAR should allow for more efficient scaling of the REXAMD method to larger biopolymer systems.
Elucidating the Inhibition Mechanism of HIV-1 Non-Nucleoside Reverse Transcriptase Inhibitors through Multi-Copy Molecular Dynamics Simulations
Anthony Ivetac and J. Andrew McCammon
Human immunodeficiency virus-1 reverse transcriptase (RT) inhibition is a major focus of current anti-AIDS drug discovery and development programs, comprising 17 of the 31 Food and Drug Administration-approved compounds. The emergence of the non-nucleoside RT inhibitor (NNRTI) class of compounds provides a highly specific and structurally diverse set of drugs, which act noncompetitively to perturb normal RT function. Despite a relatively rich set of crystallographic data of RT in various states, details of the allosteric modulation of RT dynamics by NNRTIs are lacking. Capturing this inhibitory mechanism could fuel the design of more effective inhibitors at the NNRTI site and also drive the identification of novel allosteric sites. To address this, we have performed multicopy molecular dynamics (MD) simulations of RT in the presence and absence of the NNRTI nevirapine (cumulative total simulation time, 360 ns). By comparing the collective motions of the MD and crystallographic structures, we demonstrate that the chief effect of NNRTIs is to constrain a key rigid-body motion between the "fingers" and "thumb" subdomains of the p66 subunit. We show that the NNRTI binding pocket (NNIBP) is proximal to the hinge points for this essential motion, and NNRTIs therefore act as "molecular wedges," sterically blocking the full range of motion. To explain how this impaired movement might result in the experimentally observed loss of polymerase activity, we show that the motion influences the geometry of key catalytic residues on opposite faces of the NNIBP. From a methodological point of view, our results suggest that the multicopy MD simulation approach is very useful when studying proteins that perform such large conformational changes.
Darwinian biophysics: Electrostatics and evolution in the kinetics of molecular binding (Commentary)
J. Andrew McCammon
At the molecular level of biology, the competition for favorable outcomes has been shaped by evolution, just as in more familiar examples from ecological biology. At both levels, this competition is often based on raw speed. There are differences, of course. Most notably, a race between molecules is more often determined by diffusional dynamics than by inertial dynamics. The driving forces on molecules typically comprise electrostatic nudges rather than thundering hooves digging into soil. Electrostatic interactions can be surprisingly effective, however. The rate of degradation of the neurotransmitter acetylcholine by the synaptic enzyme acetylcholinesterase is known to be increased by a factor of up to a few hundred as a result of "electrostatic steering" of the positively charged acetylcholine molecule toward the predominantly negative active-site region of the enzyme. This tends to optimize the clearing and resetting of neuromuscular junctions and other cholinergic synapses, which offered a clear competitive advantage to our successful ancestors, relative to more sluggish individuals of their species who faced the same predators. Such selective pressures are also recorded in proteins at the next level of a hierarchy, in some of the venom molecules of snakes such as the green mamba that prey on small mammals in sub-Saharan East Africa. The green mamba toxin fasciculin-2 is a small protein whose positively charged surface is attracted to, and clamps down on, the active-site entrance of acetylcholinesterase, causing muscular activity of the unfortunate rat or other prey to cease.
The role of conserved water molecules in the catalytic domain of protein kinases.
James D.R. Knight, Donald Hamelberg, J. Andrew McCammon, Rashmi Kothary
Protein kinases are essential signaling molecules with a characteristic bilobal shape that has been studied for over 15 years. Despite the number of crystal structures available, little study has been directed away from the prototypical functional elements of the kinase domain. We have performed a structural alignment of 13 active-conformation kinases and discovered the presence of six water molecules that occur in conserved locations across this group of diverse kinases. Molecular dynamics simulations demonstrated that these waters confer a great deal of stability to their local environment and to a key catalytic residue. Our results highlight the importance of novel elements within the greater kinase family and suggest that conserved water molecules are necessary for efficient kinase function.
An adaptive fast multipole boundary element method for Poisson-Boltzmann electrostatics.
Benzhuo Lu, Xiaolin Cheng, Jingfang Huang, J. Andrew McCammon
The numerical solution of the Poisson-Boltzmann (PB) equation is a useful but a computationally demanding tool for studying electrostatic solvation effects in chemical and biomolecular systems. Recently, we have described a boundary integral equation-based PB solver accelerated by a new version of the fast multipole method (FMM). The overall algorithm shows an order N complexity in both the computational cost and memory usage. Here, we present an updated version of the solver by using an adaptive FMM for accelerating the convolution type matrix-vector multiplications. The adaptive algorithm, when compared to our previous nonadaptive one, not only significantly improves the performance of the overall memory usage but also remarkably speeds the calculation because of an improved load balancing between the local- and far-field calculations. We have also implemented a node-patch discretization scheme that leads to a reduction of unknowns by a factor of 2 relative to the constant element method without sacrificing accuracy. As a result of these improvements, the new solver makes the PB calculation truly feasible for large-scale biomolecular systems such as a 30S ribosome molecule even on a typical 2008 desktop computer.
Location of Inhibitors Bound to Group IVA Phospholipase A2 Determined by Molecular Dynamics and Deuterium Exchange Mass Spectrometry
John E. Burke, Arneh Babakhani, Alemayehu A. Gorfe, George Kokotos, Sheng Li, Virgil L. Woods, Jr., J. Andrew McCammon, Edward A. Dennis
An analysis of group IVA (GIVA) phospholipase A2 (PLA2) inhibitor binding was conducted using a combination of deuterium exchange mass spectrometry (DXMS) and molecular dynamics (MD). Models of the GIVA PLA2 inhibitors pyrrophenone and the 2-oxoamide AX007 docked into the protein were designed on the basis of deuterium exchange results, and extensive molecular dynamics simulations were run to determine protein-inhibitor contacts. The models show that both inhibitors interact with key residues that also exhibit changes in deuterium exchange upon inhibitor binding. Pyrrophenone is bound to the protein through numerous hydrophobic residues located distal from the active site, while the oxoamide is bound mainly through contacts near the active site. We also show differences in protein dynamics around the active site between the two inhibitor-bound complexes. This combination of computational and experimental methods is useful in defining more accurate inhibitor binding sites and can be used in the generation of better inhibitors against GIVA PLA2.
Multiple pathways guide oxygen diffusion into flavoenzyme active sites
Riccardo Baron, Conor Riley, Pirom Chenprakhon, Kittisak Thotsaporn, Remko T. Winter, Andrea Alfieri, Federico Forneris, Willem J. H. van Berkel, Pimchai Chaiyen, Marco W. Fraaije, Andrea Mattevi, J. Andrew McCammon
Dioxygen (O2) and other gas molecules have a fundamental role in a variety of enzymatic reactions. However, it is only poorly understood which O2 uptake mechanism enzymes employ to promote efficient catalysis and how general this is. We investigated O2 diffusion pathways into monooxygenase and oxidase flavoenzymes, using an integrated computational and experimental approach. Enhanced-statistics molecular dynamics simulations reveal spontaneous protein-guided O2 diffusion from the bulk solvent to preorganized protein cavities. The predicted protein-guided diffusion paths and the importance of key cavity residues for oxygen diffusion were verified by combining site-directed mutagenesis, rapid kinetics experiments, and high-resolution X-ray structures. This study indicates that monooxygenase and oxidase flavoenzymes employ multiple funnel-shaped diffusion pathways to absorb O2 from the solvent and direct it to the reacting C4a atom of the flavin cofactor. The difference in O2 reactivity among dehydrogenases, monooxygenases, and oxidases ultimately resides in the fine modulation of the local environment embedding the reactive locus of the flavin.
Enzymatic Activity versus Structural Dynamics: The Case of Acetylcholinesterase Tetramer
Alemayehu A. Gorfe, Benzhuo Lu, Zeyun Yu, J. Andrew McCammon
The function of many proteins, such as enzymes, is modulated by structural fluctuations. This is especially the case in gated diffusion-controlled reactions (where the rates of the initial diffusional encounter and of structural fluctuations determine the overall rate of the reaction) and in oligomeric proteins (where function often requires a coordinated movement of individual subunits). A classic example of a diffusion-controlled biological reaction catalyzed by an oligomeric enzyme is the hydrolysis of synaptic acetylcholine (ACh) by tetrameric acetylcholinesterase (AChEt). Despite decades of efforts, the extent to which enzymatic efficiency of AChEt (or any other enzyme) is modulated by flexibility is not fully determined. This article attempts to determine the correlation between the dynamics of AChEt and the rate of reaction between AChEt and ACh. We employed equilibrium and nonequilibrium electro-diffusion models to compute rate coefficients for an ensemble of structures generated by molecular dynamics simulation. We found that, for the static initial model, the average reaction rate per active site is ~22–30% slower in the tetramer than in the monomer. However, this effect of tetramerization is modulated by the intersubunit motions in the tetramer such that a complex interplay of steric and electrostatic effects either guides or blocks the substrate into or from each of the four active sites. As a result, the rate per active site calculated for some of the tetramer structures is only ~15% smaller than the rate in the monomer. We conclude that structural dynamics minimizes the adverse effect of tetramerization, allowing the enzyme to maintain similar enzymatic efficiency in different oligomerization states.
Functional Dynamics of the Folded Ankyrin Repeats of IκBα Revealed by Nuclear Magnetic Resonance
Carla F. Cervantes, Phineus R. L. Markwick, Shih-Che Sue, J. Andrew McCammon, H. Jane Dyson, Elizabeth A. Komives
Inhibition of nuclear factor κB (NF-κB) is mainly accomplished by IκBα, which consists of a signal response sequence at the N-terminus, a six-ankyrin repeat domain (ARD) that binds NF-κB, and a C-terminal PEST sequence. Previous studies with the ARD revealed that the fifth and sixth repeats are only partially folded in the absence of NF-κB. Here we report NMR studies of a truncated version of IκBα, containing only the first four ankyrin repeats, IκBα(67–206). This four-repeat segment is well-structured in the free state, enabling full resonance assignments to be made. H–D exchange, backbone dynamics, and residual dipolar coupling (RDC) experiments reveal regions of flexibility. In addition, regions consistent with the presence of micro- to millisecond motions occur periodically throughout the repeat structure. Comparison of the RDCs with the crystal structure gave only moderate agreement, but an ensemble of structures generated by accelerated molecular dynamics gave much better agreement with the measured RDCs. The regions showing flexibility correspond to those implicated in entropic compensation for the loss of flexibility in ankyrin repeats 5 and 6 upon binding to NF-κB. The regions showing micro- to millisecond motions in the free protein are the ends of the β-hairpins that directly interact with NF-κB in the complex.
Interfaces and hydrophobic interactions in receptor-ligand systems: A level-set variational implicit solvent approach
Li-Tien Cheng, Zhongming Wang, Piotr Setny, Joachim Dzubiella, Bo Li, J. Andrew McCammon
A model nanometer-sized hydrophobic receptor-ligand system in aqueous solution is studied by the recently developed level-set variational implicit solvent model (VISM). This approach is compared to all-atom computer simulations. The simulations reveal complex hydration effects within the (concave) receptor pocket, sensitive to the distance of the (convex) approaching ligand. The ligand induces and controls an intermittent switching between dry and wet states of the hosting pocket, which determines the range and magnitude of the pocket-ligand attraction. In the level-set VISM, a geometric free-energy functional of all possible solute-solvent interfaces coupled to the local dispersion potential is minimized numerically. This approach captures the distinct metastable states that correspond to topologically different solute-solvent interfaces, and thereby reproduces the bimodal hydration behavior observed in the all-atom simulation. Geometrical singularities formed during the interface relaxation are found to contribute significantly to the energy barrier between different metastable states. While the hydration phenomena can thus be explained by capillary effects, the explicit inclusion of dispersion and curvature corrections seems to be essential for a quantitative description of hydrophobically confined systems on nanoscales. This study may shed more light onto the tight connection between geometric and energetic aspects of biomolecular hydration and may represent a valuable step toward the proper interpretation of experimental receptor-ligand binding rates.
Absolute single-molecule entropies from quasi-harmonic analysis of microsecond molecular dynamics: Correction terms and convergence properties
Riccardo Baron, P. H. Hünenberger, J. Andrew McCammon
The convergence properties of the absolute single-molecule configurational entropy and the correction terms used to estimate it are investigated using microsecond molecular dynamics simulation of a peptide test system and an improved methodology. The results are compared with previous applications for systems of diverse chemical nature. It is shown that (i) the effect of anharmonicity is small, (ii) the effect of pairwise correlation is typically large and, (iii) the latter affects to a larger extent the entropy estimate of thermodynamic states characterized by a higher motional correlation. The causes of such deviations from a quasi-harmonic behavior are explained. This improved approach provides entropies also for molecular systems undergoing conformational transitions and characterized by highly frustrated energy surfaces, thus not limitedly to systems sampling a single quasi-harmonic basin. Overall, this study emphasizes the need for extensive phase-space sampling in order to obtain a reliable estimation of entropic contributions.
Mechanism of Glycan Receptor Recognition and Specificity Switch for Avian, Swine and Human Adapted Influenza Virus Hemagglutinins: A Molecular Dynamics Perspective
E. Irene Newhouse, Dong Xu, Phineus R. L. Markwick, Rommie E. Amaro, Hsing C. Pao, Kevin J. Wu, Maqsudul Alam, J. Andrew McCammon, Wilfred W. Li
Hemagglutinins (HA's) from duck, swine and human influenza viruses have previously been shown to prefer avian and human glycan receptor analogues with distinct topological profiles, pentasaccharides LSTa (α-2,3 linkage) and LSTc (α-2,6 linkage), in comparative molecular dynamics studies. Based upon detailed analyses of the dynamic motions of the receptor binding domains (RBD's) and interaction energy profiles with individual glycan residues, we have identified approximately 30 residue positions in the RBD that present distinct profiles with the receptor analogues. Glycan binding constrained the conformational space sampling by the HA. Electrostatic steering appeared to play a key role in glycan binding specificity. The complex dynamic behaviors of the major SSE and trimeric interfaces with or without bound glycans suggested that networks of interactions might account for species specificity in these low affinity and high avidity (multivalent) interactions between different HA and glycans. Contact frequency, energetic decomposition and H-bond analyses revealed species-specific differences in HA-glycan interaction profiles, not readily discernable from crystal structures alone. Interaction energy profiles indicated that mutation events at the set of residues such as 145, 156, 158 and 222 would favor human or avian receptor analogues, often through interactions with distal asialo-residues. These results correlate well with existing experimental evidence, and suggest new opportunities for simulation-based vaccine and drug development.
The oxygen-binding vs. oxygen-consuming paradigm in biocatalysis: Structural biology and biomolecular simulation
Riccardo Baron, J. Andrew McCammon, Andrea Mattevi
Oxygen biocatalysis and regulation is crucial to a variety of biochemical processes in nature. Oxygen-binding proteins cover only a limited part of oxygen biocatalysis, which involves numerous examples of biocatalysts consuming oxygen coupled to low binding affinities. The integration of experiments with powerful biomolecular simulation opens appealing possibilities to investigate crucial questions on the fascinating relationship between enzyme dynamics and oxygen biocatalysis in new protein structures.
Dewetting-controlled binding of ligands to hydrophobic pockets
P. Setny, Z. Wang, L.T. Cheng, B. Li, J. A. McCammon, J. Dzubiella
We report on a combined atomistic molecular dynamics simulation and implicit solvent analysis of a generic hydrophobic pocket-ligand (host-guest) system. The approaching ligand induces complex wetting-dewetting transitions in the weakly solvated pocket. The transitions lead to bimodal solvent fluctuations which govern magnitude and range of the pocket-ligand attraction. A recently developed implicit water model, based on the minimization of a geometric functional, captures the sensitive aqueous interface response to the concave-convex pocket-ligand configuration semiquantitatively.
Toward a Unified Representation of Protein Structural Dynamics in Solution
Phineus R. L. Markwick, Guillaume Bouvignies, Loic Salmon, J. Andrew McCammon, Michael Nilges, Martin Blackledge
An atomic resolution description of protein flexibility is essential for understanding the role that structural dynamics play in biological processes. Despite the unique dependence of nuclear magnetic resonance (NMR) to motional averaging on different time scales, NMR-based protein structure determination often ignores the presence of dynamics, representing rapidly exchanging conformational equilibria in terms of a single static structure. In this study, we use the rich dynamic information encoded in experimental NMR parameters to develop a molecular and statistical mechanical characterization of the conformational behavior of proteins in solution. Critically, and in contrast to previously proposed techniques, we do not use empirical energy terms to restrain a conformational search, a procedure that can strongly perturb simulated dynamics in a non-predictable way. Rather, we use accelerated molecular dynamics simulation to gradually increase the level of conformational sampling and to identify the appropriate level of sampling via direct comparison of unrestrained simulation with experimental data. This constraint-free approach thereby provides an atomic resolution free-energy weighted Boltzmann description of protein dynamics occurring on time scales over many orders of magnitude in the protein ubiquitin.
A transition path ensemble study reveals a linchpin role for Mg2+ during rate-limiting ADP release from protein kinase A
Ilja V. Khavrutskii, Barry Grant, Susan S. Taylor, J. Andrew McCammon
Protein kinases are key regulators of diverse signaling networks critical for growth and development. Protein kinase A (PKA) is an important kinase prototype that phosphorylates protein targets at Ser and Thr residues by converting ATP to ADP. Mg2+ ions play a crucial role in regulating phosphoryl transfer and can limit overall enzyme turnover by affecting ADP release. However, the mechanism by which Mg2+ participates in ADP release is poorly understood. Here we use a novel transition path ensemble technique - the Harmonic Fourier Beads method - to explore the atomic and energetic details of the Mg2+-dependent ADP binding and release. Our studies demonstrate that adenine-driven ADP binding to PKA creates three ion-binding sites at the ADP/PKA interface that are absent otherwise. Two of these sites bind the previously characterized Mg2+ ions whereas the third site binds a monovalent cation with high affinity. This third site can bind the P-3 residue of substrate proteins and may serve as a reporter of the active site occupation. Binding of Mg2+ ions restricts mobility of the Gly-rich loop that closes over the active site. We find that simultaneous release of ADP with Mg2+ ions from the active site is unfeasible. Thus we conclude that Mg2+ ions act as a linchpin and that at least one ion must be removed prior to ADP release. Hence, the likely mechanism requires dissociation of at least one Mg2+ ion followed by the pyrophosphate-driven ADP release. The results of the present study enhance understanding of Mg2+-dependent association of nucleotides with protein kinases.