Quantitatively Accounting for Protein Reorganization in Computer-Aided Drug Design.

Structure-based drug design frequently operates under the assumption that a single holo structure is relevant. However, a large number of crystallographic examples clearly show that multiple conformations are possible. In those cases, the protein reorganization free energy must be known to accurately predict binding free energies for ligands. Only then can the energetic preference among these multiple protein conformations be utilized to design ligands with stronger binding potency and selectivity. Here, we present a computational method to quantify these protein reorganization free energies. We test it on two retrospective drug design cases, Abl kinase and HSP90, and illustrate how alternative holo conformations can be derisked and lead to large boosts in affinity. This method will allow computer-aided drug design to better support complex protein targets.

Structural and dynamic effects of paraoxon binding to human acetylcholinesterase by X-ray crystallography and inelastic neutron scattering.

Organophosphorus (OP) compounds, including nerve agents and some pesticides, covalently bind to the catalytic serine of human acetylcholinesterase (hAChE), thereby inhibiting acetylcholine hydrolysis necessary for efficient neurotransmission. Oxime antidotes can reactivate the OP-conjugated hAChE, but reactivation efficiency can be low for pesticides, such as paraoxon (POX). Understanding structural and dynamic determinants of OP inhibition and reactivation can provide insights to design improved reactivators. Here, X-ray structures of hAChE with unaged POX, with POX and oximes MMB4 and RS170B, and with MMB4 are reported. A significant conformational distortion of the acyl loop was observed upon POX binding, being partially restored to the native conformation by oximes. Neutron vibrational spectroscopy combined with molecular dynamics simulations showed that picosecond vibrational dynamics of the acyl loop soften in the ∼20-50 cm-1 frequency range. The acyl loop structural perturbations may be correlated with its picosecond vibrational dynamics to yield more comprehensive template for structure-based reactivator design.

Computational study of non-conductive selectivity filter conformations and C-type inactivation in a voltage-dependent potassium channel.

C-type inactivation is a time-dependent process of great physiological significance that is observed in a large class of K+ channels. Experimental and computational studies of the pH-activated KcsA channel show that the functional C-type inactivated state, for this channel, is associated with a structural constriction of the selectivity filter at the level of the central glycine residue in the signature sequence, TTV(G)YGD. The structural constriction is allosterically promoted by the wide opening of the intracellular activation gate. However, whether this is a universal mechanism for C-type inactivation has not been established with certainty because similar constricted structures have not been observed for other K+ channels. Seeking to ascertain the general plausibility of the constricted filter conformation, molecular dynamics simulations of a homology model of the pore domain of the voltage-gated potassium channel Shaker were performed. Simulations performed with an open intracellular gate spontaneously resulted in a stable constricted-like filter conformation, providing a plausible nonconductive state responsible for C-type inactivation in the Shaker channel. While there are broad similarities with the constricted structure of KcsA, the hypothetical constricted-like conformation of Shaker also displays some subtle differences. Interestingly, those are recapitulated by the Shaker-like E71V KcsA mutant, suggesting that the residue at this position along the pore helix plays a pivotal role in determining the C-type inactivation behavior. Free energy landscape calculations show that the conductive-to-constricted transition in Shaker is allosterically controlled by the degree of opening of the intracellular activation gate, as observed with the KcsA channel. The behavior of the classic inactivating W434F Shaker mutant is also characterized from a 10-µs MD simulation, revealing that the selectivity filter spontaneously adopts a nonconductive conformation that is constricted at the level of the second glycine in the signature sequence, TTVGY(G)D.

Covalent inhibition of hAChE by organophosphates causes homodimer dissociation through long-range allosteric effects.

Acetylcholinesterase (EC 3.1.1.7), a key acetylcholine-hydrolyzing enzyme in cholinergic neurotransmission, is present in a variety of states in situ, including monomers, C-terminally disulfide-linked homodimers, homotetramers, and up to three tetramers covalently attached to structural subunits. Could oligomerization that ensures high local concentrations of catalytic sites necessary for efficient neurotransmission be affected by environmental factors? Using small-angle X-ray scattering (SAXS) and cryo-EM, we demonstrate that homodimerization of recombinant monomeric human acetylcholinesterase (hAChE) in solution occurs through a C-terminal four-helix bundle at micromolar concentrations. We show that diethylphosphorylation of the active serine in the catalytic gorge or isopropylmethylphosphonylation by the RP enantiomer of sarin promotes a 10-fold increase in homodimer dissociation. We also demonstrate the dissociation of organophosphate (OP)-conjugated dimers is reversed by structurally diverse oximes 2PAM, HI6, or RS194B, as demonstrated by SAXS of diethylphosphoryl-hAChE. However, binding of oximes to the native ligand-free hAChE, binding of high-affinity reversible ligands, or formation of an SP-sarin-hAChE conjugate had no effect on homodimerization. Dissociation monitored by time-resolved SAXS occurs in milliseconds, consistent with rates of hAChE covalent inhibition. OP-induced dissociation was not observed in the SAXS profiles of the double-mutant Y337A/F338A, where the active center gorge volume is larger than in wildtype hAChE. These observations suggest a key role of the tightly packed acyl pocket in allosterically triggered OP-induced dimer dissociation, with the potential for local reduction of acetylcholine-hydrolytic power in situ. Computational models predict allosteric correlated motions extending from the acyl pocket toward the four-helix bundle dimerization interface 25 Å away.

The Activation of c-Src Tyrosine Kinase: Conformational Transition Pathway and Free Energy Landscape.

Tyrosine kinases are important cellular signaling allosteric enzymes that regulate cell growth, proliferation, metabolism, differentiation, and migration. Their activity must be tightly controlled, and malfunction can lead to a variety of diseases, particularly cancer. The nonreceptor tyrosine kinase c-Src, a prototypical model system and a representative member of the Src-family, functions as complex multidomain allosteric molecular switches comprising SH2 and SH3 domains modulating the activity of the catalytic domain. The broad picture of self-inhibition of c-Src via the SH2 and SH3 regulatory domains is well characterized from a structural point of view, but a detailed molecular mechanism understanding is nonetheless still lacking. Here, we use advanced computational methods based on all-atom molecular dynamics simulations with explicit solvent to advance our understanding of kinase activation. To elucidate the mechanism of regulation and self-inhibition, we have computed the pathway and the free energy landscapes for the "inactive-to-active" conformational transition of c-Src for different configurations of the SH2 and SH3 domains. Using the isolated c-Src catalytic domain as a baseline for comparison, it is observed that the SH2 and SH3 domains, depending upon their bound orientation, promote either the inactive or active state of the catalytic domain. The regulatory structural information from the SH2-SH3 tandem is allosterically transmitted via the N-terminal linker of the catalytic domain. Analysis of the conformational transition pathways also illustrates the importance of the conserved tryptophan 260 in activating c-Src, and reveals a series of concerted events during the activation process.

Full Atom Simulations of Spin Label Conformations.

Interpretation of EPR measurables from spin labels in terms of structure and dynamics requires knowledge of label behavior. General strategies were developed for simulation of labels used in EPR of proteins. The criteria for those simulations are (a) exhaustive sampling of rotamer space, (b) consensus of results independent of starting points, and (c) inclusion of entropy. These criteria are satisfied only when the number of transitions in any dihedral angle exceeds 100 and the simulation maintains thermodynamic equilibrium. Methods such as conventional MD do not efficiently cross energetic barriers, simulated anealing, Monte Carlo or popular Rotamer Library methodologies are potential energy based and ignore entropy (in addition to their specific shortcomings: environment fluctuations, fixed environment, or electrostatics). The Simulated scaling method avoids the above flaws by modulating the force fields between a reduced (allowing crossing energy barriers) and full potential (sampling minima). Spin label diffuses on this surface while remaining in thermodynamic equilibrium. Simulations show that (a) adopting a single conformation is rare, often there are two to four populated rotamers and (b) position of the NO varies up to 16 Å. These results illustrate necessity for caution when interpreting EPR signals in terms of molecular structure. For example, the 10-16 Å distance change in DEER should not be interpreted as a large conformational change, it can well be a flip about Cα-Cß bond. Rigorous exploration of possible rotamer structures of a spin label is paramount in signal interpretation. We advocate use of bifunctional labels, motion of which is restricted 10,000-fold and the NO position is restricted to 2-5 Å.

Generalized Scalable Multiple Copy Algorithms for Molecular Dynamics Simulations in NAMD.

Computational methodologies that couple the dynamical evolution of a set of replicated copies of a system of interest offer powerful and flexible approaches to characterize complex molecular processes. Such multiple copy algorithms (MCAs) can be used to enhance sampling, compute reversible work and free energies, as well as refine transition pathways. Widely used examples of MCAs include temperature and Hamiltonian-tempering replica-exchange molecular dynamics (T-REMD and H-REMD), alchemical free energy perturbation with lambda replica-exchange (FEP/λ-REMD), umbrella sampling with Hamiltonian replica exchange (US/H-REMD), and string method with swarms-of-trajectories conformational transition pathways. Here, we report a robust and general implementation of MCAs for molecular dynamics (MD) simulations in the highly scalable program NAMD built upon the parallel programming system Charm++. Multiple concurrent NAMD instances are launched with internal partitions of Charm++ and located continuously within a single communication world. Messages between NAMD instances are passed by low-level point-to-point communication functions, which are accessible through NAMD's Tcl scripting interface. The communication-enabled Tcl scripting provides a sustainable application interface for end users to realize generalized MCAs without modifying the source code. Illustrative applications of MCAs with fine-grained inter-copy communication structure, including global lambda exchange in FEP/λ-REMD, window swapping US/H-REMD in multidimensional order parameter space, and string method with swarms-of-trajectories were carried out on IBM Blue Gene/Q to demonstrate the versatility and massive scalability of the present implementation.

w-REXAMD: A Hamiltonian Replica Exchange Approach to Improve Free Energy Calculations for Systems with Kinetically Trapped Conformations.

Free energy governs the equilibrium extent of many biological processes. High barriers separating free energy minima often limit the sampling in molecular dynamics (MD) simulations, leading to inaccurate free energies. Here, we demonstrate enhanced sampling and improved free energy calculations, relative to conventional MD, using windowed accelerated MD within a Hamiltonian replica exchange framework (w-REXAMD). We show that for a case in which multiple conformations are separated by large free energy barriers, w-REXAMD is a useful enhanced sampling technique, without any necessary reweighting.

On the Role of Dewetting Transitions in Host-Guest Binding Free Energy Calculations.

We use thermodynamic integration (TI) and explicit solvent molecular dynamics (MD) simulation to estimate the absolute free energy of host-guest binding. In the unbound state, water molecules visit all of the internally accessible volume of the host, which is fully hydrated on all sides. Upon binding of an apolar guest, the toroidal host cavity is fully dehydrated; thus, during the intermediate λ stages along the integration, the hydration of the host fluctuates between hydrated and dehydrated states. Estimating free energies by TI can be especially challenging when there is a considerable difference in hydration between the two states of interest. We investigate these aspects using the popular TIP3P and TIP4P water models. TI free energy estimates through MD largely depend on water-related interactions, and water dynamics significantly affect the convergence of binding free energy calculations. Our results indicate that wetting/dewetting transitions play a major role in slowing the convergence of free energy estimation. We employ two alternative approaches-one analytical and the other empirically based on actual MD sampling-to correct for the standard state free energy. This correction is sizable (up to 4 kcal/mol), and the two approaches provide corrections that differ by about 1 kcal/mol. For the system considered here, the TIP4P water model combined with an analytical correction for the standard state free energy provides higher overall accuracy. This observation might be transferable to other systems in which water-related contributions dominate the binding process.

Using multistate free energy techniques to improve the efficiency of replica exchange accelerated molecular dynamics.

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.

Intrinsic Free Energy of the Conformational Transition of the KcsA Signature Peptide from Conducting to Nonconducting State.

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 nonconducting 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 side chains 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 nonconducting 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.

Replica-Exchange Accelerated Molecular Dynamics (REXAMD) Applied to Thermodynamic Integration.

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.

Mapping electron paramagnetic resonance spin label conformations by the simulated scaling method.

In order to efficiently simulate spin label behavior when attached to the protein backbone we developed a novel approach that enhances local conformational sampling. The simulated scaling (SS) approach (Li, H., et al. J. Chem. Phys. 2007, 126, 24106) couples the random walk of a potential scaling parameter and molecular dynamics in the framework of hybrid Monte Carlo. This approach allows efficient barrier crossings between conformations. The method retains the thermodynamic detailed balance allowing for determination of relative free energies between various conformations. The accuracy of our method was validated by comparison with the recently resolved X-ray crystal structure of a spin labeled T4 lysozyme in which the spin label was in the interior of the protein. Consistent potentials of mean force (PMF) are obtained for the spin label torsion angles to illustrate their behavior in various protein environments: surface, semiburied, and buried. These PMFs reflect the experimentally observed trends and provide the rationale for the spin label dynamics. We have used this method to compare an implicit and explicit solvent model in spin label modeling. The implicit model, which is computationally faster, was found to be in excellent agreement with the explicit solvent treatment. Based on this collection of results, we believe that the presented approach has great potential in the general strategy of describing the behavior of the spin label using molecular modeling and using this information in the interpretation of EPR measurements in terms of protein conformation and dynamics.

Simulated scaling method for localized enhanced sampling and simultaneous "alchemical" free energy simulations: a general method for molecular mechanical, quantum mechanical, and quantum mechanical/molecular mechanical simulations.

A potential scaling version of simulated tempering is presented to efficiently sample configuration space in a localized region. The present "simulated scaling" method is developed with a Wang-Landau type of updating scheme in order to quickly flatten the distributions in the scaling parameter lambdam space. This proposal is meaningful for a broad range of biophysical problems, in which localized sampling is required. Besides its superior capability and robustness in localized conformational sampling, this simulated scaling method can also naturally lead to efficient "alchemical" free energy predictions when dual-topology alchemical hybrid potential is applied; thereby simultaneously, both of the chemically and conformationally distinct portions of two end point chemical states can be efficiently sampled. As demonstrated in this work, the present method is also feasible for the quantum mechanical and quantum mechanical/molecular mechanical simulations.