To request any journal reprints, please contact us.
A nonlinear elasticity model of macromolecular conformational change induced by electrostatic forces
Y.C. Zhou, Michael Holst and J. Andrew McCammon
In this paper we propose a nonlinear elasticity model of macromolecular conformational change (deformation) induced by electrostatic forces generated by an implicit solvation model. The Poisson-Boltzmann equation for the electrostatic potential is analyzed in a domain varying with the elastic deformation of molecules, and a new continuous model of the electrostatic forces is developed to ensure solvability of the nonlinear elasticity equations. We derive the estimates of electrostatic forces corresponding to four types of perturbations to an electrostatic potential field, and establish the existance of an equilibrium configuration using a fixed-point argument, under the assumption that the change in the ionic strength and charges due to the additional molecules causing the deformation are sufficiently small.The results are valid for elastic models with arbitrarily complex dielectric interfaces and cavities, and can be generalized to large elastic deformation caused by high ionic strength, large charges, and strong external fields by using continuation methods.
Springs and Speeds in Free Energy Reconstruction from Irreversible Single-Molecule Pulling Experiments
David D.L. Minh and J. Andrew McCammon
The nonequilibrium work relation allows for the calculation of equilibrium free energy differences between states based on the exponential average of accumulated work from irreversible transitions. Here, we compare two distinct approaches of calculating free energy surfaces from unidirectional singlemolecule pulling experiments: the stiff spring approximation and the Hummer-Szabo method. First, we perform steered molecular dynamics simulations to mechanically stretch the model peptide deca-alanine using harmonic potentials with different spring stiffnesses and at various constant pulling velocities. We then calculate free energy surfaces based on the two methods and their variants, including the first and second cumulant expansion of the exponentially weighted work and the Gaussian position approximation for the delta function in Hummer and Szabo's expression. We find that with large harmonic force constants, the second cumulant expansion in conjunction with either the stiff spring approximation or the Hummer-Szabo method perform well. When interpreting dynamic force spectroscopy (pullings at different speeds), the second cumulant expansion of the stiff spring approximation performs the best when pulling velocities are similar, but variants of the Hummer-Szabo perform the best when they are spread over a large spectrum. While these conclusion are not definitive for all systems, the insights should prove useful for scientists interpreting nonequilibrium pulling experiments.
Continuum Simulations of Acetylcholine Consumption by Acetylcholinesterase: A Poisson-Nernst-Planck Approach
Y.C. Zhou, Benzhuo Lu, Gary A. Huber, Michael J. Holst and J. Andrew McCammon
The Poisson-Nernst-Planck (PNP) equation provides a continuum description of electrostatic-driven diffusion and is used here to model the diffusion and reaction of acetylcholine (ACh) with acetylcholinesterase (AChE) enzymes. This study focuses on the effects of ion and substrate concentrations on the reaction rate and rate coefficient. To this end, the PNP equations are numerically solved with a hybrid finite element and boundary element method at a wide range of ion and substrate concentrations, and the results are compared with the partially coupled Smoluchowski-Poisson-Boltzmann model. The reaction rate is found to depend strongly on the concentrations of both the substrate and ions; this is explained by the competition between the intersubstrate repulsion and the ionic screening effects. The reaction rate coefficient is independent of the substrate concentration only at very high ion concentrations, whereas at low ion concentrations the behavior of the rate depends strongly on the substrate concentration. Moreover, at physiological ion concentrations, variations in substrate concentration significantly affect the transient behavior of the reaction. Our results offer a reliable estimate of reaction rates at various conditions and imply that the concentrations of charged substrates must be coupled with the electrostatic computation to provide a more realistic description of neurotransmission and other electrodiffusion and reaction processes.
Dynamics of the Acetylcholinesterase Tetramer
Alemayehu A. Gorfe, Chia-en A. Chang, Ivaylo Ivanov and J. Andrew McCammon
Acetylcholinesterase rapidly hydrolyzes the neurotransmitter acetylcholine in cholinergic synapses, including the neuromuscular junction. The tetramer is the most important functional form of the enzyme. Two low-resolution crystal structures have been solved. One is compact with two of its four peripheral anionic sites (PAS) sterically blocked by complementary subunits. The other is a loose tetramer with all four subunits accessible to solvent. These structures lacked the C-terminal amphipathic t-peptide (WAT domain) that interacts with the proline-rich attachment domain (PRAD). A complete tetramer model (AChEt) was built based on the structure of the PRAD/WAT complex and the compact tetramer. Normal mode analysis suggested that AChEt could exist in multiple conformations with subunits fluctuating relative to one another. Here, a multiscale simulation involving all-atom molecular dynamics and Cbeta-based coarse-grained Brownian dynamics simulations was carried out to investigate the large scale inter-subunit dynamics in AChEt. We sampled the ns-micros time scale motions and found that the tetramer indeed constitutes a dynamic assembly of monomers. The inter-subunit fluctuation is correlated with the occlusion of the PAS. Such motions of the subunits "gate" ligand-protein association. The gates are open more than 80% of the time on average, which suggests a small reduction of ligand-protein binding. Despite the limitations in the starting model and approximations inherent in coarse graining, these results are consistent with experiments which suggest that binding of a substrate to the PAS is only somewhat hindered by the association of the subunits.
Computing accurate potentials of mean force in electrolyte solutions with the generalized gradient-augmented harmonic Fourier beads method
Ilja V. Khavrutskii, Joachim Dzubiella and J. Andrew McCammon
We establish the accuracy of the novel generalized gradient-augmented Harmonic Fourier Beads (ggaHFB) method in computing free-energy profiles or potentials of mean force (PMFs) through comparison with two independent conventional techniques. In particular, we employ umbrella sampling with 1D weighted histogram analysis method (WHAM) and free molecular dynamics simulation of radial distribution functions to compute the PMF for the Na+-Cl- ion pair separation to 16 Å in 1.0 M NaCl solution in water. The corresponding ggaHFB free-energy profile in 6D Cartesian space is in excellent agreement with the conventional benchmarks. We then explore changes in the PMF in response to lowering the NaCl concentration to physiological 0.3 and 0.1 M, and dilute 0.0 M concentrations. Finally, to expand the scope of the ggaHFB method, we formally develop the free-energy gradient approximation in arbitrary nonlinear coordinates. This formal development underscores the importance of the logarithmic Jacobian correction to reconstruct true PMFs from umbrella sampling simulations with either WHAM or ggaHFB techniques when nonlinear coordinate restraints are used with Cartesian propagators. The ability to employ nonlinear coordinates and high accuracy of the computed free-energy profiles further advocate the use of the ggaHFB method in studies of rare events in complex systems.
Molecular surface-free continuum model for electrodiffusion processes
Benzhuo Lu and J. Andrew McCammon
Incorporation of van der Waals interactions enables the continuum model of electrodiffusion in biomolecular system to avoid the artifacts of introducing a molecular surface and the painful task of the surface mesh generation. Calculation examples show that the electrostatics, diffusion-reaction kinetics, and molecular surface defined as an isosurface of a certain density distribution can be extracted from the solution of the Poisson-Nernst-Planck equations using this model. The molecular surface-free model enables a wider usage of some modern numerical methodologies such as finite element methods for biomolecular modeling, and sheds light on a new paradigm of continuum modeling for biomolecular systems.
An Improved Relaxed Complex Scheme for Receptor Flexibility in Computer-Aided Drug Design
Rommie E. Amaro, Riccardo Baron and J. Andrew McCammon
The interactions among associating (macro) molecules are dynamic, which adds to the complexity of molecular recognition. While ligand flexibility is well accounted for in computational drug design, the effective inclusion of receptor flexibility remains an important challenge. The relaxed complex scheme (RCS) is a promising computational methodology that combines the advantages of docking algorithms with dynamic structural information provided by molecular dynamics (MD) simulations, therefore explicitly accounting for the flexibility of both the receptor and the docked ligands. Here, we briefly review the RCS and discuss new extensions and improvements of this methodology in the context of ligand binding to two example targets: kinetoplastid RNA editing ligase 1 and the W191G cavity mutant of cytochrome c peroxidase. The RCS improvements include its extension to virtual screening, more rigorous characterization of local and global binding effects, and methods to improve its computational efficiency by reducing the receptor ensemble to a representative set of configurations. The choice of receptor ensemble, its influence on the predictive power of RCS, and the current limitations for an accurate treatment of the solvent contributions are also briefly discussed. Finally, we outline potential methodological improvements that we anticipate will assist future development.
Novel Druggable Hot Spots in Avian Influenza Neuraminidase H5N1 Revealed by Computational Solvent Mapping of a Reduced and Representative Receptor Ensemble
Melissa R. Landon, Rommie E. Amaro, Riccardo Baron, Chi Ho Ngan, David Ozonoff, J. Andrew McCammon and Sandor Vajda
The influenza virus subtype H5N1 has raised concerns of a possible human pandemic threat because of its high virulence and mutation rate. Although several approved anti-influenza drugs effectively target the neuraminidase, some strains have already acquired resistance to the currently available anti-influenza drugs. In this study, we present the synergistic application of extended explicit solvent molecular dynamics (MD) and computational solvent mapping (CS-Map) to identify putative “hot spots” within flexible binding regions of N1 neuraminidase. Using representative conformations of the N1 binding region extracted from a clustering analysis of four concatenated 40-ns MD simulations, CS-Map was utilized to assess the ability of small, solvent-sized molecules to bind within close proximity to the sialic acid binding region. Mapping analyses of the dominant MD conformations reveal the presence of additional hot spot regions in the 150- and 430-loop regions. Our hot spot analysis provides further support for the feasibility of developing high-affinity inhibitors capable of binding these regions, which appear to be unique to the N1 strain.
Electrostatic Free Energy and its Variations in Implicit Solvent Models
Jianwei Che, Joachim Dzubiella, Bo Li and J. Andrew McCammon
A mean-field approach to the electrostatics for solutes in electrolyte solution is revisited and rigorously justified. In this approach, an electrostatic free energy functional is constructed that depends solely on the local ionic concentrations. The unique set of such concentrations that minimize this free energy are given by the usual Boltzmann distributions through the electrostatic potential which is determined by the Poisson-Boltzmann equation. This approach is then applied to the variational implicit solvent description of the solvation of molecules [Dzubiella, Swanson, McCammon, Phys. Rev. Lett. 2006, 96, 087802; J. Chem. Phys. 2006, 124, 084905]. Care is taken for the singularities of the potential generated by the solute point charges. The variation of the electrostatic free energy with respect to the location change of solute-solvent interfaces, that is, dielectric boundaries, is derived. Such a variation gives rise to the normal component of the effective surface force per unit surface area that is shown to be attractive to the fixed point charges in the solutes. Two examples of applications are given to validate the analytical results. The first one is a one-dimensional model system resembling, for example, a charged solute or cavity in a one-dimensional channel. The second one, which is of its own interest, is the electrostatic free energy of a charged sphercal solute immersed in an ionic solution. An analytical formula is derived for the Debye-Hückel approximation of the free energy, extending the classical Born's formula to one that includes ionic concentrations. Variations of the nonlinear Poisson-Boltzmann free energy are also obtained.
Control of cation permeation through the nicotinic receptor channel
Hai-Long Wang, Xiaolin Cheng, Palmer Taylor, J. Andrew McCammon and Steven M. Sine
We used molecular dynamics (MD) simulations to explore the transport of single cations through the channel of the muscle nicotinic acetylcholine receptor (nAChR). Four MD simulations of 16 ns were performed at physiological and hyperpolarized membrane potentials, with and without restraints of the structure, but all without bound agonist. With the structure unrestrained and a potential of -100 mV, one cation traversed the channel during a transient period of channel hydration; at -200 mV, two cations traversed the channel while the channel was continuously hydrated. With the structure restrained, however, no cations traverse at either membrane potential, even though the channel was continuously hydrated. The overall results show that cation selective transport through nAChR channel is governed by electrostatics interactions to achieve charge selectivity, but relies on trans-membrane potential, channel hydration and protein dynamics to enable ion passage.
Recent Progress in Numerical Methods for the Poisson-Boltzmann Equation in Biophysical Applications
B.Z. Lu, Y.C. Zhou, M.J. Holst and J.A. McCammon
Efficiency and accuracy are two major concerns in numerical solutions of the Poisson-Boltzmann equation for applications in chemistry and biophysics. Recent developments in boundary element methods, interface methods, adaptive methods, finite element methods, and other approaches for the Poisson-Boltzmann equation as well as related mesh generation techniques are reviewed. We also discussed the challenging problems and possible future work, in particular, for the aim of biophysical applications.
A novel switch region regulates H-ras membrane orientation and signal output
Daniel Abankwa, Michael Hanzal-Bayer, Nicolas Ariotti, Sarah J. Plowman, Alemayehu A. Gorfe, Robert G. Parton, J. Andrew McCammon and John F. Hancock
The plasma membrane nanoscale distribution of H-ras is regulated by guanine nucleotide binding. To explore the structural basis of H-ras membrane organization, we combined molecular dynamic simulations and medium-throughput FRET measurements on live cells. We extracted a set of FRET values, termed a FRET vector, to describe the lateral segregation and orientation of H-ras with respect to a large set of nanodomain markers. We show that mutation of basic residues in helix α4 or the hypervariable region (HVR) selectively alter the FRET vectors of GTP- or GDP-loaded H-ras, demonstrating a critical role for these residues in stabilizing GTP- or GDP-H-ras interactions with the plasma membrane. By a similar analysis, we find that the β2-β3 loop and helix α5 are involved in a novel conformational switch that operates through helix α4 and the HVR to reorient the H-ras G-domain with respect to the plasma membrane. Perturbation of these switch elements enhances MAPK activation by stabilizing GTP-H-ras in a more productive signalling conformation. The results illustrate how the plasma membrane spatially constrains signalling conformations by acting as a semi-neutral interaction partner.
Inhibition of Cathepsin B by Au(I) Complexes: A Kinetic and Computational Study
Shamila S. Gunatilleke, Cesar Augusto F. de Oliveira, J. Andrew McCammon and Amy M. Barrios
Gold(I) compounds have been used in the treatment of rheumatoid arthritis for over 80 years, but the biological targets and the structure-activity relationships of these drugs are not well understood. Of particular interest is the molecular mechanism behind the antiarthritic activity of the orally available drug triethylphosphine(2,3,4,6-tetra-O-acetyl-β-1-d-thiopyranosato-S) gold(I) (auranofin, Ridaura). The cathepsin family of lysosomal, cysteine-dependent enzymes is an attractive biological target of Au(I) and is inhibited by auranofin and auranofin analogs with reasonable potency. Here we employ a combination of experimental and computational investigations into the effect of changes in the phosphine ligand of auranofin on its in vitro inhibition of cathepsin B. Sequential replacement of the ethyl substituents of triethylphosphine by phenyl groups leads to increasing potency in the resultant Au(I) complexes, due in large part to favorable interactions of the more sterically bulky Au(I)-PR3 fragments with the enzyme active site.
Feature-Preserving Adaptive Mesh Generation for Molecular Shape Modeling and Simulation
Zeyun Yu, Michael J. Holst, Yuhui Cheng and J. Andrew McCammon
We describe a chain of algorithms for molecular surface and volumetric mesh generation. We take as inputs the centers and radii of all atoms of a molecule and the toolchain outputs both triangular and tetrahedral meshes that can be used for molecular shape modeling and simulation. Experiments on a number of molecules are demonstrated, showing that our methods possess several desirable properties: feature-preservation, local adaptivity, high quality, and smoothness (for surface meshes). We also demonstrate an example of molecular simulation using the finite element method and the meshes generated by our method. The approaches presented and their implementations are also applicable to other types of inputs such as 3D scalar volumes and triangular surface meshes with low quality, and hence can be used for generation/improvment of meshes in a broad range of applications.
Catalytically Requisite Conformational Dynamics in the mRNA-Capping Enzyme Probed by Targeted Molecular Dynamics
Robert V. Swift and J. Andrew McCammon
The addition of a N7-methyl guanosine cap to the 5' end of nascent mRNA is carried out by the mRNA-capping enzyme, a two-domain protein that is a member of the nucleotidyltransferase superfamily. The mRNA-capping enzyme is composed of a catalytic nucleotidyltransferase domain and a noncatalytic oligonucleotide/oligosaccharide binding (OB) domain. Large-scale domain motion triggered by substrate binding mediates catalytically requisite conformational rearrangement of the GTP substrate prior to the chemical step. In this study, we employ targeted molecular dynamics (TMD) on the PBCV-1 capping enzyme to probe the global domain dynamics and internal dynamics of conserved residues during the conformational transformation from the open to the closed state. Analysis of the resulting trajectories along with structural and sequence homology to other members of the superfamily allows us to suggest a conserved mechanism of conformational rearrangements spanning all mRNA-capping enzymes and all ATP-dependent DNA ligases. Our results suggest that the OB domain moves quasi-statically toward the nucleotidyltransferase domain, pivoting about a short linker region. The approach of the OB domain brings a conserved RxDK sequence, an element of conserved motif VI, within proximity of the triphosphate of GTP, destabilizing the unreactive conformation and thereby allowing thermal fluctuations to partition the substrate toward the catalytically competent state.
(Thermo)dynamic role of receptor flexibility, entropy, and motional correlation in protein-ligand binding
Riccardo Baron and J. Andrew McCammon
The binding of 2-amino-5-methylthiazole to the W191G cavity mutant of cytochrome c peroxidase is an ideal test case to investigate the entropic contribution to the binding free energy due to changes of receptor flexibility. The dynamic and thermodynamic role of receptor flexibility were studied by 50-ns long explicit-solvent molecular dynamics simulations of three separate receptor ensembles: W191G binding a K+ ion, W191G-2a5mt complex with a closed 190-195 gating loop, and apo with an open loop. We employ a method recently proposed to estimate accurate absolute single-molecule configurational entropies and their differences for systems undergoing conformational transitions. We find that receptor flexibility plays a generally-underestimated role in protein-ligand binding (thermo)dynamics and that changes of receptor motional correlation determine such large entropy contributions.
Computer-aided Drug Discovery: Physics-based Simulations from the Molecular to the Cellular Level
J. Andrew McCammon
Mapping the nucleotide and isoform dependent structural and dynamical features of Ras proteins
Alemayehu A. Gorfe, Barry J. Grant and J. Andrew McCammon
Ras GTPases are conformational switches controlling cell proliferation, differentiation and development. Despite their prominent role in many forms of cancer, the mechanism of conformational transition between inactive GDP- and active GTP-bound states remains unclear. Here we describe a detailed analysis of available experimental structures and molecular dynamics simulations to quantitatively assess the structural and dynamical features of active and inactive states and their interconversion. We demonstrate that GTP-bound and nucleotide-free G12V H-ras sample a wide region of conformational space, and show that the inactive to active transition is a multiphase process defined by the relative rearrangement of the two switches and the orientation of Tyr32. We also modeled and simulated N- and K-ras proteins and found that K-ras is more flexible than N- and H-ras. We identified a number of isoform-specific long-range side chain interactions that define unique pathways of communication between the nucleotide binding site and the C-terminus.
Intrinsic Conformational Flexibility of Acetylcholinesterase
Jennifer M. Bui and J. Andrew McCammon
Proteins have been metaphorically described – due to the introduction and extraordinary advances in biomolecular dynamics and computational biophysics over the past decades – as "kicking and screaming" molecules [G. Weber, Adv. Protein Chem. 29 (1975) 1-83]. In fact, dynamic fluctuations in protein structural conformation have been known to play an important role in protein function. However, fundamental mechanisms by which protein fluctuations couple with catalytic function of particular enzymes remain poorly understood. To understand the dynamical properties of acetylcholinesterase (AChE) in rapid termination of cationic neurotransmitter, acetylcholine at neurosynaptic junctions, multiple molecular dynamics (MD) trajectories of AChE in the presence and absence of its inhibitors [J.M. Bui, J.A. McCammon, Proc. Natl. Acad. Sci. U.S.A. 103 (2006) 15451-15456; J.M. Bui, Z. Radić, P. Taylor, J.A. McCammon, Biophys. J. 90 (2006) 3280-3287; J.M. Bui, K. Tai, J.A. McCammon, J. Am. Chem. Soc. 126 (2004) 7198-7205; J.M. Bui, R.H. Henchman, J.A. McCammon, Biophys. J. 85 (2003) 2267-2272] have been conducted and correlated with its inhibitory mechanisms. The intrinsic flexibilities of AChE, particularly of the long omega loop, are important in facilitating the ligand's inhibition of the enzyme.
Entropic contributions and the influence of the hydrophobic environment in promiscuous protein-protein association
Chia-en A. Chang, William A. McLaughlin, Riccardo Baron, Wei Wang and J. Andrew McCammon
The mechanisms by which a promiscuous protein can strongly interact with several different proteins using the same binding interface are not completely understood. An example is protein kinase A (PKA), which uses a single face on its docking/dimerization domain to interact with multiple A-kinase anchoring proteins (AKAP) that localize it to different parts of the cell. In the current study, the configurational entropy contributions to the binding between the AKAP protein HT31 with the D/D domain of RII α-regulatory subunit of PKA were examined. The results show that the majority of configurational entropy loss for the interaction was due to decreased fluctuations within rotamer states of the side chains. The result is in contrast to the widely held approximation that the decrease in the number of rotamer states available to the side chains forms the major component. Further analysis showed that there was a direct linear relationship between total configurational entropy and the number of favorable, alternative contacts available within hydrophobic environments. The hydrophobic binding pocket of the D/D domain provides alternative contact points for the side chains of AKAP peptides that allow them to adopt different binding conformations. The increase in binding conformations provides an increase in binding entropy and hence binding affinity. We infer that a general strategy for a promiscuous protein is to provide alternative contact points at its interface to increase binding affinity while the plasticity required for binding to multiple partners is retained. Implications are discussed for understanding and treating diseases in which promiscuous protein interactions are used.
High-fidelity geometric modeling for biomedical applications
Zeyun Yu, Michael J. Holst and J. Andrew McCammon
We describe a combination of algorithms for high-fidelity geometric modeling and mesh generation. Although our methods and implementations are application-neutral, our primary target application is multiscale biomedical models that range in scales across the molecular, cellular, and organ levels. Our software toolchain implementing these algorithms is general in the sense that it can take as input a molecule in PDB/PQR forms, a 3D scalar volume, or a user-defined triangular surface mesh that may have very low quality. The main goal of our work presented is to generate high quality and smooth surface triangulations from the aforementioned inputs, and to reduce the mesh sizes by mesh coarsening. Tetrahedral meshes are also generated for finite element analysis in biomedical applications. Experiments on a number of bio-structures are demonstrated, showing that our approach possesses several desirable properties: feature-preservation, local adaptivity, high quality, and smoothness (for surface meshes). The availability of this software toolchain will give researchers in computational biomedicine and other modeling areas access to higher-fidelity geometric models.
Ensemble-based Virtual Screening Reveals Potential Novel Antiviral Compounds for Avian Influenza Neuraminidase
Lily S. Cheng, Rommie E. Amaro, Dong Xu, Wilfred W. Li, Peter W. Arzberger and J. Andrew McCammon
Avian influenza virus subtype H5N1 is a potential pandemic threat with human-adapted strains resistant to antiviral drugs. Although virtual screening (VS) against a crystal or relaxed receptor structure is an established method to identify potential inhibitors, the more dynamic changes within binding sites are neglected. To accommodate full receptor flexibility, we use AutoDock4 to screen the NCI diversity set against representative receptor ensembles extracted from explicitly solvated molecular dynamics simulations of the neuraminidase system. The top hits are redocked to the entire nonredundant receptor ensemble and rescored using the relaxed complex scheme (RCS). Of the 27 top hits reported, half ranked very poorly if only crystal structures are used. These compounds target the catalytic cavity as well as the newly identified 150- and 430-cavities, which exhibit dynamic properties in electrostatic surface and geometric shape. This ensemble-based VS and RCS approach may offer improvement over existing strategies for structure-based drug discovery.
One-Bead Coarse-Grained Models for Proteins
Valentina Tozzini and J. Andrew McCammon
Most biomolecular processes involve macromolecular aggregates (on the size scale of 10-100 nm or more, including the cell membranes) and occur on a time scale of microseconds to milliseconds (or even hours to days, including folding and amyloid aggregation, for instance). Computer simulations based on atomic force fields are not yet able to reach these scales, since they are currently restricted in most cases to systems comprising fewer than a million atoms for times of less than 1 μs. In consideration of these facts, the idea of simplifying the description of a macromolecular system by including groups of atoms in a single interaction center (coarse-graining, CG) in order to reduce the number of interanal degrees of freedom and, with them, the computational cost, is quite simple and somewhat natural, even considering the hierarchical structural organization of the proteins and nucleic acids.
Acetylcholinesterase: mechanisms of covalent inhibition of H447I mutant determined by computational analyses
Y.H. Cheng, X.L. Cheng, Z. Radić and J.A. McCammon
The reaction mechanisms of two inhibitor TFK+ and TFK0 binding to H447I mutant mouse acetylcholinesterase (mAChE) have been investigated by using a combined ab initio quantum mechanical/molecular mechanical (QM/MM) approach and classical molecular dynamics (MD) simulations. TFK+ binding to the H447I mutant may proceed with a different reaction mechanism from the wild type. A water molecule takes over the role of His447 and participates in the bond breaking and forming as a "charge relayer". Unlike in the wild-type mAChE case, Glu334, a conserved residue from the catalytic triad, acts as a catalytic base in the reaction. The calculated energy barrier for this reaction is about 8 kcal/mol. These predictions await experimental verification. In the case of the neutral ligand TFK0, however, multiple MD simulations on the TFK0/H447I complex reveal that none of the water molecules can be retained in the active site as a "catalytic" water. Taken together our computational studies confirm that TFK0 is almost inactive in the H447I mutant, and also provide detailed mechanistic insights into the experimental observations.
Hot-spot residues at the E9/Im9 interface help binding via different mechanisms
Sergio E. Wong, Riccardo Baron and J. Andrew McCammon
Protein-protein association involves many interface interactions, but they do not contribute equally. Ala scanning experiments reveal that only a few mutations significantly lower binding affinity. These key residues, which appear to drive protein-protein association, are called hot-spot residues. Molecular dynamics simulations of the Colicin E9/Im9 complex show Im9 Glu41 and Im9 Ser50, both hot-spots, bind via different mechanisms. The results suggest Im9 Ser50 restricts Glu41 in a conformation auspicious for salt-bridge formation across the interface. This type of model may be helpful in engineering hot-spot clusters at protein-protein interfaces and, consequently, the design of specificity.
Water-membrane partition thermodynamics of an amphiphilic lipopeptide: An enthalpy-driven hydrophobic effect
Alemayehu A. Gorfe, Riccardo Baron and J. Andrew McCammon
To shed light on the driving force for the hydrophobic effect that partitions amphiphilic lipoproteins between water and membrane, we carried out an atomically-detailed thermodynamic analysis of a triply lipid modified H-ras heptapeptide (ANCH) in water and in a DMPC bilayer. Combining molecular mechanical and continuum solvent approaches with an improved technique for solute entropy calculation, we obtained an overall transfer free energy of ~ -13 kcal/mol. This value is in qualitative agreement with free energy changes derived from a potential of mean force calculation and indirect experimental observations. Changes in free energies of solvation and ANCH conformational reorganization are unfavorable while ANCH-DMPC interactions – especially van der Waals – favor insertion. These results are consistent with an enthalpy-driven hydrophobic effect, in accord with earlier calorimetric data on the membrane partition of other amphiphiles. Furthermore, structural and entropic analysis of molecular dynamics (MD)-generated ensembles suggests that conformational selection may play a hitherto unappreciated role in membrane insertion of lipid-modified peptides and proteins.
Intrinsic free energy of the conformational transition of the KcsA signature peptide from conducting to non-conducting state
Ilja V. Khavrutskii, Mikolai Fajer and J. Andrew McCammon
We explore a conformational transition of the TATTVGYG signature peptide of the KcsA ion selectivity filter and its GYG to AYA mutant from the conducting α-strand state into the non-conducting pII-like state using a novel technique for multidimensional optimization of transition path ensembles and free energy calculations. We find that the wild type peptide, unlike the mutant, intrinsically favors the conducting state due to G77 backbone propensities and additional hydrophobic interaction between the V76 and Y78 sidechains in water. The molecular mechanical free energy profiles in explicit water are in very good agreement with the corresponding adiabatic energies from the Generalized Born Molecular Volume (GBMV) implicit solvent model. However comparisons of the energies to higher level B3LYP/6-31G(d) Density Functional Theory calculations with Polarizable Continuum Model (PCM) suggest that the non-conducting state might be more favorable than predicted by molecular mechanics simulations. By extrapolating the single peptide results to the tetrameric channel, we propose a novel hypothesis for the ion selectivity mechanism.
Thermodynamics of Peptide Insertion and Aggregation in a Lipid Bilayer
Arneh Babakhani, Alemayehu A. Gorfe, Judy E. Kim and J. Andrew McCammon
A variety of biomolecules mediate physiological processes by inserting and reorganizing in cell membranes, and the thermodynamic forces responsible for their partitioning are of great interest. Recent experiments provided valuable data on the free energy changes associated with the transfer of individual amino acids from water to membrane. However, a complete picture of the pathways and the associated changes in energy of peptide insertion into a membrane remains elusive. To this end, computational techniques supplement the experimental data with atomic-level details and shed light on the energetics of insertion. Here, we employed the technique of umbrella sampling in a total 850 ns of all-atom molecular dynamics simulations to study the free energy profile and the pathway of insertion of a model hexapeptide consisting of a tryptophan and five leucines (WL5). The computed free energy profile of the peptide as it travels from bulk solvent through the membrane core exhibits two minima: a local minimum at the water-membrane interface or the head group region; and a global minimum at the hydrophobic-hydrophilic interface close to the lipid glycerol region. A rather small barrier of roughly 1 kcal/mol exists at the membrane core, which is explained by the enhanced flexibility of the peptide when deeply-inserted. Combining our results with those in the literature, we present a thermodynamic model for peptide insertion and aggregation which involves peptide aggregation upon contact with the membrane at the solvent-lipid head group interface.
Coupling Accelerated Molecular Dynamics Methods with Thermodynamic Integration Simulations
César Augusto F. de Oliveira, Donald Hamelberg and J. Andrew McCammon
In this work we propose a straightforward and efficient approach to improve accuracy and convergence of free energy simulations in condensed-phase systems. We also introduce a new accelerated Molecular Dynamics (MD) approach in which molecular conformational transitions are accelerated by lowering the energy barriers while the potential surfaces near the minima are left unchanged. All free energy calculations were performed on the propane-to-propane model system. The accuracy of free energy simulations was significantly improved when sampling of internal degrees of freedom of solute was enhanced. However, accurate and converged results were only achieved when the solvent interactions were taken into account in the accelerated MD approaches. The analysis of the distribution of boost potential along the free energy simulations showed that the new accelerated MD approach samples efficiently both low- and high-energy regions of the potential surface. Since this approach also maintains substantial populations in regions near the minima, the statistics are not compromised in the thermodynamic integration calculations, and, as a result, the ensemble average can be recovered.
Molecular dynamics of a κB DNA element: base flipping via cross-strand intercalative stacking in a microsecond-scale simulation
Cameron Mura and J. Andrew McCammon
The sequence-dependent structural variability and conformational dynamics of DNA play pivotal roles in many biological milieus, such as in the site-specific binding of transcription factors to target regulatory elements. To better understand DNA structure, function, and dynamics in general, and protein-DNA recognition in the κB family of genetic regulatory elements in particular, we performed molecular dynamics simulations of a 20-bp DNA encompassing a cognate κB site recognized by the proto-oncogenic 'c-Rel' subfamily of NF-κB transcription factors. Simulations of the κB DNA in explicit water were extended to microsecond duration, providing a broad, atomically detailed glimpse into the structural and dynamical behavior of double helical DNA over many timescales. Of particular note, novel (and structurally plausible) conformations of DNA developed only at the long times sampled in this simulation -- including a peculiar state arising at ~0.7 μs and characterized by cross-strand intercalative stacking of nucleotides within a longitudinally sheared base pair, followed (at ~1 μs) by spontaneous base flipping of a neighboring thymine within the A-rich duplex. Results and predictions from the microsecond-scale simulation include implications for a dynamical NF-κB recognition motif, and are amenable to testing and further exploration via specific experimental approaches that are suggested herein.
Replica Exchange Accelerated Molecular Dynamics (REXAMD) Applied to Thermodynamic Integration
Mikolai Fajer, Donald Hamelberg and J. Andrew McCammon
Accelerated molecular dynamics (AMD) is an efficient strategy for accelerating the sampling of molecular dynamics simulations, and observable quantities such as free energies derived on the biased AMD potential can be reweighted to yield results consistent with the original, unmodified potential. In conventional AMD the reweighting procedure has an inherent statistical problem in systems with large acceleration, where the points with the largest biases will dominate the reweighted result and reduce the effective number of data points. We propose a replica exchange of various degrees of acceleration (REXAMD) to retain good statistics while achieving enhanced sampling. The REXAMD method is validated and benchmarked on two simple gas-phase model systems, and two different strategies for computing reweighted averages over a simulation are compared.
Similar Membrane Affinity of Mono- and Di-S-acylated Ras Membrane Anchors: A New Twist in the Role of Protein Lipidation
Alemayehu A. Gorfe and J. Andrew McCammon
The functionally required membrane attachment of Ras is achieved through an invariant isoprenylation of a C-terminal Cys, supplemented by further lipid modification of adjacent Cys residues by one (N-ras) or two (H-ras) palmitoyls. However, whether the triply lipidated membrane anchor of H-ras has a higher membrane affinity than its doubly lipidated counterpart, or whether the affinity contribution of the two palmitates and the farnesyl is additive, was not known. To address this issue, we carried out potential of mean force (PMF or free energy profile) calculations on a hexadecylated but nonpalmitoylated anchor (Cys186-HD), hexadecylated and monopalmitoylated anchors (Cys181-monopalmitate and Cys184-monopalmitate), and a nonlipid-modified anchor. We found that the overall insertion free energy follows the trend Cys181/Cys184-bipalmitate (wild type) ≈ Cys181-monopalmitate > Cys184-monopalmitate >> nonpalmitoylated anchor. Consistent with suggestions from recent cell biological experiments, the computed PMFs, coupled with structural analysis, demonstrate that membrane affinity of the Ras anchor depends on both the hydrophobicity of the palmitate and the prenyl groups and the spacing between them. The data further suggest that while Cys181-palmitate and Cys186-farnesyl together provide sufficient hydrophobic force for tight membrane binding, the palmitoyl at Cys184 is likely designed to serve another, presumably functional, role.
Diffusional Channeling in the Sulfate Activating Complex: Combined Continuum Modeling and Coarse-grained Brownian Dynamics Studies
Yuhui Cheng, Chia-en A. Chang, Zeyun Yu, Yongjie Zhang, Meihao Sun, Thomas S. Leyh, Michael J. Holst and J. Andrew McCammon
Enzymes required for sulfur metabolism have been suggested to gain effciency by restricted diffusion ("channeling") of an intermediate APS2- between active sites. This article describes modeling of the whole channeling process by numerical solution of the Smoluchowski diffusion equation, as well as by coarse-grained Brownian dynamics. The results suggest that electrostatics plays an essential role in the APS2- channeling. Furthermore, with coarse-grained Brownian dynamics, the substrate channeling process has been studied with reactions in multiple active sites. Our simulations provide a bridge for numerical modeling with Brownian dynamics to simulate the complicated reaction and diffusion and raise important questions relating to the electrostatically mediated substrate channeling in vitro, in situ and in vivo.
E9-Im9 Colicin DNase-Immunity Protein Biomolecular Association in Water: A Multiple-Copy and Accelerated Molecular Dynamics Simulation Study
Riccardo Baron, Sergio E. Wong, César A.F. de Oliveira and J. Andrew McCammon
Protein-protein transient and dynamic interactions underlie all biological processes. The molecular dynamics (MD) of the E9 colicin DNase protein, its Im9 inhibitor protein, and their E9-Im9 recognition complex are investigated by combining multiple-copy (MC) MD and accelerated MD (aMD) explicit-solvent simulation approaches, after validation with crystalline-phase and solution experiments. Im9 shows higher flexibility than its E9 counterpart. Im9 displays a significant reduction of backbone flexibility and increase in motional correlation upon E9 association. Im9 loops 23-31 and 54-64 open with respect to the E9-Im9 X-ray structure and show high conformational diversity. Upon association a large fraction (~20 nm2) of E9 and Im9 protein surfaces become inaccessible to water. Numerous salt bridges transiently occurring throughout our six 50-ns long MC-MD simulations are not present in the X-ray model. Among these Im9 Glu31–E9 Arg96 and Im9 Glu41–Lys89 involve interface interactions. Using 10-ns of Im9 aMD simulation we reconcile the largest thermodynamic impact measured for Asp51Ala mutation with Im9 structure and dynamics. Lys57 acts as an essential molecular switch to shift Im9 surface loop towards an ideal configuration for E9 inhibition. This is achieved by switching Asp60–Lys57 and Asp62–Lys57 hydrogen bonds to Asp51–Lys57 salt bridge. E9-Im9 recognition involves shifts of conformational distributions, re-organization of intra-molecular hydrogen bond patterns, and formation of new inter- and intra-molecular interactions. The description of key transient biological interactions can be significantly enriched by the dynamic and atomic-level information provided by computer simulations.
Three-Dimensional Geometric Modeling of Membrane-bound Organelles in Ventricular Myocytes: Bridging the Gap between Microscopic Imaging and Mathematical Simulation
Zeyun Yu, Michael J. Holst, Takeharu Hayashi, Chandrajit L. Bajaj, Mark H. Ellisman, J. Andrew McCammon and Masahiko Hoshijima
A general framework of image-based geometric processing is presented to bridge the gap between three-dimensional (3D) imaging that provides structural details of a biological system and mathematical simulation where high-quality surface or volumetric meshes are required. A 3D density map is processed in the order of image pre-processing (contrast enhancement and anisotropic filtering), feature extraction (boundary segmentation and skeletonization), and high-quality and realistic surface (triangular) and volumetric (tetrahedral) mesh generation. While the tool-chain described is applicable to general types of 3D imaging data, the performance is demonstrated specifically on membrane-bound organelles in ventricular myocytes that are imaged and reconstructed with electron microscopic (EM) tomography and two-photon microscopy (T-PM). Of particular interest in this study are two types of membrane-bound Ca2+-handling organelles, namely, transverse tubules (T-tubules) and junctional sarcoplasmic reticulum (jSR), both of which play an important role in regulating the excitation-contraction (E-C) coupling through dynamic Ca2+ mobilization in cardiomyocytes.
Discovery of drug-like inhibitors of an essential RNA-editing ligase in Trypanosoma brucei
Rommie E. Amaro, Achim Schnaufer, Heidrun Interthal, Wim Hol, Kenneth D. Stuart and J. Andrew McCammon
Trypanosomatid RNA editing is a unique process and essential for these organisms. It therefore represents a drug target for a group of protozoa that includes the causative agents for African sleeping sickness and other devastating tropical and subtropical diseases. Here, we present drug-like inhibitors of a key enzyme in the editing machinery, RNA-editing ligase 1 (REL1). These inhibitors were identified through a strategy employing molecular dynamics to account for protein flexibility. A virtual screen of the REL1 crystal structure against the National Cancer Institute Diversity Set was performed by using AutoDock4. The top 30 compounds, predicted to interact with REL1's ATP-binding pocket, were further refined by using the relaxed complex scheme (RCS), which redocks the compounds to receptor structures extracted from an explicitly solvated molecular dynamics trajectory. The resulting reordering of the ligands and filtering based on drug-like properties resulted in an initial recommended set of 8 ligands, 2 of which exhibited micromolar activity against REL1. A subsequent hierarchical similarity search with the most active compound over the full National Cancer Institute database and RCS rescoring resulted in an additional set of 6 ligands, 2 of which were confirmed as REL1 inhibitors with IC50 values of ≈1 μM. Tests of the 3 most promising compounds against the most closely related bacteriophage T4 RNA ligase 2, as well as against human DNA ligase IIIβ, indicated a considerable degree of selectivity for RNA ligases. These compounds are promising scaffolds for future drug design and discovery efforts against these important pathogens.