Open science grid study of the coupling between conformation and water content in the interior of a protein

Reformulation of the self-guided molecular simulation method

Self-guided molecular/Langevin dynamics (SGMD/SGLD) simulation methods were developed to enhance conformational sampling through promoting low frequency motion of molecular systems and have been successfully applied in many simulation studies. Quantitative understanding of conformational distribution in SGLD has been achieved by separating microscopic properties according to frequency. However, a missing link between the guiding factors and conformational distributions makes it highly empirical and system dependent when choosing the values of the guiding parameters. Based on the understanding that molecular interactions are the source of energy barriers and diffusion friction, this work reformulates the equation of the low frequency motion to resemble Langevin dynamics. This reformulation leads to new forms of guiding forces and establishes a relation between the guiding factors and conformational distributions. We call simulations with these new guiding forces the generalized self-guided molecular/Langevin dynamics (SGMDg/SGLDg). In addition, we present a new way to calculate low frequency properties and an efficient algorithm to implement SGMDg/SGLDg that minimizes memory usage and inter-processor communication. Through example simulations with a skewed double well system, an argon fluid, and a cryo-EM map flexible fitting case, we demonstrate the guiding effects on conformational distributions and conformational searching.

q-Canonical Monte Carlo Sampling for Modeling the Linkage between Charge Regulation and Conformational Equilibria of Peptides

The overall charge content and the patterning of charged residues have a profound impact on the conformational ensembles adopted by intrinsically disordered proteins. These parameters can be altered by charge regulation, which refers to the effects of post-translational modifications, pH-dependent changes to charge, and conformational fluctuations that modify the pKa values of ionizable residues. Although atomistic simulations have played a prominent role in uncovering the major sequence-ensemble relationships of IDPs, most simulations assume fixed charge states for ionizable residues. This may lead to erroneous estimates for conformational equilibria if they are linked to charge regulation. Here, we report the development of a new method we term q-canonical Monte Carlo sampling for modeling the linkage between charge regulation and conformational equilibria. The method, which is designed to be interoperable with the ABSINTH implicit solvation model, operates as follows: For a protein sequence with n ionizable residues, we start with all 2n charge microstates and use a criterion based on model compound pKa values to prune down to a subset of thermodynamically relevant charge microstates. This subset is then grouped into mesostates, where all microstates that belong to a mesostate have the same net charge. Conformational distributions, drawn from a canonical ensemble, are generated for each of the charge microstates that make up a mesostate using a method we designate as proton walk sampling. This method combines Metropolis Monte Carlo sampling in conformational space with an auxiliary Markov process that enables interconversions between charge microstates along a mesostate. Proton walk sampling helps identify the most likely charge microstate per mesostate. We then use thermodynamic integration aided by the multistate Bennett acceptance ratio method to estimate the free energies for converting between mesostates. These free energies are then combined with the per-microstate weights along each mesostate to estimate standard state free energies and pH-dependent free energies for all thermodynamically relevant charge microstates. The results provide quantitative estimates of the probabilities and preferred conformations associated with every thermodynamically accessible charge microstate. We showcase the application of q-canonical sampling using two model systems. The results establish the soundness of the method and the importance of charge regulation in systems characterized by conformational heterogeneity.

Self-guided Langevin dynamics via generalized Langevin equation

Self-guided Langevin dynamics (SGLD) is a molecular simulation method that enhances conformational search and sampling via acceleration of the low frequency motions of the system. This acceleration is produced via introduction of a guiding force which breaks down the detailed-balance property of the dynamics, implying that some reweighting is necessary to perform equilibrium sampling. Here, we eliminate the need of reweighing and show that the NVT and NPT ensembles are sampled exactly by a new version of self-guided motion involving a generalized Langevin equation (GLE) in which the random force is modified so as to restore detailed-balance. Through the examples of alanine dipeptide and argon liquid, we show that this SGLD-GLE method has enhanced conformational sampling capabilities compared with regular Langevin dynamics (LD) while being of comparable computational complexity. In particular, SGLD-GLE is fully size extensive and can be used in arbitrarily large systems, making it an appealing alternative to LD. © 2015 Wiley Periodicals, Inc.

Understanding the basis of a class of paradoxical mutations in AraC through simulations

Most mutations at position 15 in the N-terminal arm of the regulatory protein AraC leave the protein incapable of responding to arabinose and inducing the proteins required for arabinose catabolism. Mutations at other positions of the arm do not have this behavior. Simple energetic analysis of the interactions between the arm and bound arabinose do not explain the uninducibility of AraC with mutations at position 15. Extensive molecular dynamics (MD) simulations, carried out largely on the Open Science Grid, were done of the wild-type protein with and without bound arabinose and of all possible mutations at position 15, many of which were constructed and measured for this work. Good correlation was found for deviation of arm position during the simulations and inducibility as measured in vivo of the same mutant proteins. Analysis of the MD trajectories revealed that preservation of the shape of the arm is critical to inducibility. To maintain the correct shape of the arm, the strengths of three interactions observed to be strong in simulations of the wild-type AraC protein need to be preserved. These interactions are between arabinose and residue 15, arabinose and residues 8-9, and residue 13 and residue 15. The latter interaction is notable because residues L9, Y13, F15, W95, and Y97 form a hydrophobic cluster which needs to be preserved for retention of the correct shape.

Replica exchanging self-guided Langevin dynamics for efficient and accurate conformational sampling

This work presents a replica exchanging self-guided Langevin dynamics (RXSGLD) simulation method for efficient conformational searching and sampling. Unlike temperature-based replica exchanging simulations, which use high temperatures to accelerate conformational motion, this method uses self-guided Langevin dynamics (SGLD) to enhance conformational searching without the need to elevate temperatures. A RXSGLD simulation includes a series of SGLD simulations, with simulation conditions differing in the guiding effect and/or temperature. These simulation conditions are called stages and the base stage is one with no guiding effect. Replicas of a simulation system are simulated at the stages and are exchanged according to the replica exchanging probability derived from the SGLD partition function. Because SGLD causes less perturbation on conformational distribution than high temperatures, exchanges between SGLD stages have much higher probabilities than those between different temperatures. Therefore, RXSGLD simulations have higher conformational searching ability than temperature based replica exchange simulations. Through three example systems, we demonstrate that RXSGLD can generate target canonical ensemble distribution at the base stage and achieve accelerated conformational searching. Especially for large systems, RXSGLD has remarkable advantages in terms of replica exchange efficiency, conformational searching ability, and system size extensiveness.

Enhanced Sampling in Free Energy Calculations: Combining SGLD with the Bennett's Acceptance Ratio and Enveloping Distribution Sampling Methods

One of the key requirements for the accurate calculation of free energy differences is proper sampling of conformational space. Especially in biological applications, molecular dynamics simulations are often confronted with rugged energy surfaces and high energy barriers, leading to insufficient sampling and, in turn, poor convergence of the free energy results. In this work, we address this problem by employing enhanced sampling methods. We explore the possibility of using self-guided Langevin dynamics (SGLD) to speed up the exploration process in free energy simulations. To obtain improved free energy differences from such simulations, it is necessary to account for the effects of the bias due to the guiding forces. We demonstrate how this can be accomplished for the Bennett's acceptance ratio (BAR) and the enveloping distribution sampling (EDS) methods. While BAR is considered among the most efficient methods available for free energy calculations, the EDS method developed by Christ and van Gunsteren is a promising development that reduces the computational costs of free energy calculations by simulating a single reference state. To evaluate the accuracy of both approaches in connection with enhanced sampling, EDS was implemented in CHARMM. For testing, we employ benchmark systems with analytical reference results and the mutation of alanine to serine. We find that SGLD with reweighting can provide accurate results for BAR and EDS where conventional molecular dynamics simulations fail. In addition, we compare the performance of EDS with other free energy methods. We briefly discuss the implications of our results and provide practical guidelines for conducting free energy simulations with SGLD.

Efficient and Unbiased Sampling of Biomolecular Systems in the Canonical Ensemble: A Review of Self-Guided Langevin Dynamics

This review provides a comprehensive description of the self-guided Langevin dynamics (SGLD) and the self-guided molecular dynamics (SGMD) methods and their applications. Example systems are included to provide guidance on optimal application of these methods in simulation studies. SGMD/SGLD has enhanced ability to overcome energy barriers and accelerate rare events to affordable time scales. It has been demonstrated that with moderate parameters, SGLD can routinely cross energy barriers of 20 kT at a rate that molecular dynamics (MD) or Langevin dynamics (LD) crosses 10 kT barriers. The core of these methods is the use of local averages of forces and momenta in a direct manner that can preserve the canonical ensemble. The use of such local averages results in methods where low frequency motion "borrows" energy from high frequency degrees of freedom when a barrier is approached and then returns that excess energy after a barrier is crossed. This self-guiding effect also results in an accelerated diffusion to enhance conformational sampling efficiency. The resulting ensemble with SGLD deviates in a small way from the canonical ensemble, and that deviation can be corrected with either an on-the-fly or a post processing reweighting procedure that provides an excellent canonical ensemble for systems with a limited number of accelerated degrees of freedom. Since reweighting procedures are generally not size extensive, a newer method, SGLDfp, uses local averages of both momenta and forces to preserve the ensemble without reweighting. The SGLDfp approach is size extensive and can be used to accelerate low frequency motion in large systems, or in systems with explicit solvent where solvent diffusion is also to be enhanced. Since these methods are direct and straightforward, they can be used in conjunction with many other sampling methods or free energy methods by simply replacing the integration of degrees of freedom that are normally sampled by MD or LD.

The pKa Cooperative: a collaborative effort to advance structure-based calculations of pKa values and electrostatic effects in proteins

The pK(a) Cooperative (http://www.pkacoop.org) was organized to advance development of accurate and useful computational methods for structure-based calculation of pK(a) values and electrostatic energies in proteins. The Cooperative brings together laboratories with expertise and interest in theoretical, computational, and experimental studies of protein electrostatics. To improve structure-based energy calculations, it is necessary to better understand the physical character and molecular determinants of electrostatic effects. Thus, the Cooperative intends to foment experimental research into fundamental aspects of proteins that depend on electrostatic interactions. It will maintain a depository for experimental data useful for critical assessment of methods for structure-based electrostatics calculations. To help guide the development of computational methods, the Cooperative will organize blind prediction exercises. As a first step, computational laboratories were invited to reproduce an unpublished set of experimental pK(a) values of acidic and basic residues introduced in the interior of staphylococcal nuclease by site-directed mutagenesis. The pK(a) values of these groups are unique and challenging to simulate owing to the large magnitude of their shifts relative to normal pK(a) values in water. Many computational methods were tested in this first Blind Prediction Challenge and critical assessment exercise. A workshop was organized in the Telluride Science Research Center to objectively assess the performance of many computational methods tested on this one extensive data set. This volume of Proteins: Structure, Function, and Bioinformatics introduces the pK(a) Cooperative, presents reports submitted by participants in the Blind Prediction Challenge, and highlights some of the problems in structure-based calculations identified during this exercise.

pH replica-exchange method based on discrete protonation states

We propose a new algorithm for obtaining proton titration curves of ionizable residues. The algorithm is a pH replica-exchange method (PHREM), which is based on the constant pH algorithm of Mongan et al. (J Comput Chem 2004;25:2038-2048). In the original replica-exchange method, simulations of different replicas are performed at different temperatures, and the temperatures are exchanged between the replicas. In our PHREM, simulations of different replicas are performed at different pH values, and the pHs are exchanged between the replicas. The PHREM was applied to a blocked amino acid and to two protein systems (snake cardiotoxin and turkey ovomucoid third domain), in conjunction with a generalized Born implicit solvent. The performance and accuracy of this algorithm and the original constant pH method (PHMD) were compared. For a single set of simulations at different pHs, the use of PHREM yields more accurate Hill coefficients of titratable residues. By performing multiple sets of constant pH simulations started with different initial states, the accuracy of predicted pK(a) values and Hill coefficients obtained with PHREM and PHMD methods becomes comparable. However, the PHREM algorithm exhibits better samplings of the protonation states of titratable residues and less scatter of the titration points and thus better precision of measured pK(a) values and Hill coefficients. In addition, PHREM exhibits faster convergence of individual simulations than the original constant pH algorithm.

Toward canonical ensemble distribution from self-guided Langevin dynamics simulation

This work derives a quantitative description of the conformational distribution in self-guided Langevin dynamics (SGLD) simulations. SGLD simulations employ guiding forces calculated from local average momentums to enhance low-frequency motion. This enhancement in low-frequency motion dramatically accelerates conformational search efficiency, but also induces certain perturbations in conformational distribution. Through the local averaging, we separate properties of molecular systems into low-frequency and high-frequency portions. The guiding force effect on the conformational distribution is quantitatively described using these low-frequency and high-frequency properties. This quantitative relation provides a way to convert between a canonical ensemble and a self-guided ensemble. Using example systems, we demonstrated how to utilize the relation to obtain canonical ensemble properties and conformational distributions from SGLD simulations. This development makes SGLD not only an efficient approach for conformational searching, but also an accurate means for conformational sampling.

Combined computational and biochemical study reveals the importance of electrostatic interactions between the "pH sensor" and the cation binding site of the sodium/proton antiporter NhaA of Escherichia coli

Sodium proton antiporters are essential enzymes that catalyze the exchange of sodium ions for protons across biological membranes. The crystal structure of NhaA has provided a basis to explore the mechanism of ion exchange and its unique regulation by pH. Here, the mechanism of the pH activation of the antiporter is investigated through functional and computational studies of several variants with mutations in the ion-binding site (D163, D164). The most significant difference found computationally between the wild type antiporter and the active site variants, D163E and D164N, are low pK(a) values of Glu78 making them insensitive to pH. Although in the variant D163N the pK(a) of Glu78 is comparable to the physiological one, this variant cannot demonstrate the long-range electrostatic effect of Glu78 on the pH-dependent structural reorganization of trans-membrane helix X and, hence, is proposed to be inactive. In marked contrast, variant D164E remains sensitive to pH and can be activated by alkaline pH shift. Remarkably, as expected computationally and discovered here biochemically, D164E is viable and active in Na(+)/H(+) exchange albeit with increased apparent K(M). Our results unravel the unique electrostatic network of NhaA that connect the coupled clusters of the "pH sensor" with the binding site, which is crucial for pH activation of NhaA.