Free energies of solvation in the context of protein folding: Implications for implicit and explicit solvent models

Implicit solvent models for biomolecular simulations have been developed to use in place of more expensive explicit models; however, these models make many assumptions and approximations that are likely to affect accuracy. Here, the changes in free energies of solvation upon folding of several fast folding proteins are calculated from previously run μs–ms simulations with a number of implicit solvent models and compared to the values needed to be consistent with the explicit solvent model used in the simulations. In the majority of cases, there is a significant and substantial difference between the values calculated from the two approaches that is robust to the details of the calculations. These differences could only be remedied by selecting values for the model parameters—the internal dielectric constant for the polar term and the surface tension coefficient for the nonpolar term—that were system‐specific or physically unrealistic. We discuss the potential implications of our findings for both implicit and explicit solvent simulations. © 2015 Wiley Periodicals, Inc.     

Efficient approximate all-atom solvent accessible surface area method parameterized for folded and denatured protein conformations

Continuing advances in computer hardware and software are permitting atomic-resolution molecular simulations for longer time scales and on larger systems. Despite these advances, routinely performing atomistic simulations with explicit water for even small proteins, which reach the folding time of such proteins, remains intractable for the foreseeable future. An implicit approximation of the solvent environment using a solvent accessible surface area (SASA) term in a molecular mechanics potential function allows exclusion of the explicit water molecules in protein simulations. This reduces the number of particles by approximately an order of magnitude. We present a fast and acceptably accurate approximate all-atom SASA method parameterized using a set of folded and heat-denatured conformations of globular proteins. The parameters are shown to be transferable to folded and heat-denatured conformations for another set of proteins. Calculation of the approximate SASA and the associated derivatives with respect to atomic positions for a 4644 atom protein requires only 1/11th the CPU time required for calculation of the nonbonded interactions for this system. On a per atom basis, this algorithm is three times faster than the fastest previously published approximate SASA method and achieves the same level of accuracy.     

Breaking the polar‐nonpolar division in solvation free energy prediction

Implicit solvent models divide solvation free energies into polar and nonpolar additive contributions, whereas polar and nonpolar interactions are inseparable and nonadditive. We present a feature functional theory (FFT) framework to break this ad hoc division. The essential ideas of FFT are as follows: (i) representability assumption: there exists a microscopic feature vector that can uniquely characterize and distinguish one molecule from another; (ii) feature‐function relationship assumption: the macroscopic features, including solvation free energy, of a molecule is a functional of microscopic feature vectors; and (iii) similarity assumption: molecules with similar microscopic features have similar macroscopic properties, such as solvation free energies. Based on these assumptions, solvation free energy prediction is carried out in the following protocol. First, we construct a molecular microscopic feature vector that is efficient in characterizing the solvation process using quantum mechanics and Poisson–Boltzmann theory. Microscopic feature vectors are combined with macroscopic features, that is, physical observable, to form extended feature vectors. Additionally, we partition a solvation dataset into queries according to molecular compositions. Moreover, for each target molecule, we adopt a machine learning algorithm for its nearest neighbor search, based on the selected microscopic feature vectors. Finally, from the extended feature vectors of obtained nearest neighbors, we construct a functional of solvation free energy, which is employed to predict the solvation free energy of the target molecule. The proposed FFT model has been extensively validated via a large dataset of 668 molecules. The leave‐one‐out test gives an optimal root‐mean‐square error (RMSE) of 1.05 kcal/mol. FFT predictions of SAMPL0, SAMPL1, SAMPL2, SAMPL3, and SAMPL4 challenge sets deliver the RMSEs of 0.61, 1.86, 1.64, 0.86, and 1.14 kcal/mol, respectively. Using a test set of 94 molecules and its associated training set, the present approach was carefully compared with a classic solvation model based on weighted solvent accessible surface area. © 2017 Wiley Periodicals, Inc.     

Application of the frozen atom approximation to the GB/SA continuum model for solvation free energy

The generalized Born/surface area (GB/SA) continuum model for solvation free energy is a fast and accurate alternative to using discrete water molecules in molecular simulations of solvated systems. However, computational studies of large solvated molecular systems such as enzyme-ligand complexes can still be computationally expensive even with continuum solvation methods simply because of the large number of atoms in the solute molecules. Because in such systems often only a relatively small portion of the system such as the ligand binding site is under study, it becomes less attractive to calculate energies and derivatives for all atoms in the system. To curtail computation while still maintaining high energetic accuracy, atoms distant from the site of interest are often frozen; that is, their coordinates are made invariant. Such frozen atoms do not require energetic and derivative updates during the course of a simulation. Herein we describe methodology and results for applying the frozen atom approach to both the generalized Born (GB) and the solvent accessible surface area (SASA) parts of the GB/SA continuum model for solvation free energy. For strictly pairwise energetic terms, such as the Coulombic and van-der-Waals energies, contributions from pairs of frozen atoms can be ignored. This leaves energetic differences unaffected for conformations that vary only in the positions of nonfrozen atoms. Due to the nonlocal nature of the GB analytical form, however, excluding such pairs from a GB calculation leads to unacceptable inaccuracies. To apply a frozen-atom scheme to GB calculations, a buffer region within the frozen-atom zone is generated based on a user-definable cutoff distance from the nonfrozen atoms. Certain pairwise interactions between frozen atoms in the buffer region are retained in the GB computation. This allows high accuracy in conformational GB comparisons to be maintained while achieving significant savings in computational time compared to the full (nonfrozen) calculation. A similar approach for using a buffer region of frozen atoms is taken for the SASA calculation. The SASA calculation is local in nature, and thus exact SASA energies are maintained. With a buffer region of 8 A for the frozen-atom cases, excellent agreement in differences in energies for three different conformations of cytochrome P450 with a bound camphor ligand are obtained with respect to the nonfrozen cases. For various minimization protocols, simulations run 2 to 10.5 times faster and memory usage is reduced by a factor of 1.5 to 5. Application of the frozen atom method for GB/SA calculations thus can render computationally tractable biologically and medically important simulations such as those used to study ligand-receptor binding conformations and energies in a solvated environment.     

Calculation of solvation free energy using RISM theory for peptide in salt solution

We developed a robust, highly efficient algorithm for solving the full reference interaction site model (RISM) equations for salt solutions near a solute molecule with many atomic sites. It was obtained as an extension of our previously reported algorithm for pure water near the solute molecule. The algorithm is a judicious hybrid of the Newton–Raphson and Picard methods. The most striking advantage is that the Jacobian matrix is just part of the input data and need not be recalculated at all. To illustrate the algorithm, we solved the full RISM equations for a dipeptide (NH₂(SINGLE BOND)CHCH₃(SINGLE BOND)CONH(SINGLE BOND)CHCH₃(SINGLE BOND)COOH) in a 1 M NaCl solution. The extended simple point charge (SPC/E) model was employed for water molecules. Two different conformations of the dipeptide were considered. It was assumed for each conformation that the dipeptide was present either as an un-ionized form or as a zwitterion. The structure of the salt solution near the dipeptide and salt effects on the solvation free energy were also discussed. © 1998 John Wiley & Sons, Inc. J Comput Chem 19: 1724–1735, 1998     

A component analysis of the free energies of folding of 35 proteins: A consensus view on the thermodynamics of folding at the molecular level

What factors favor protein folding? This is a textbook question. Parsing the experimental free energies of folding/unfolding into diverse enthalpic and entropic components of solute and solvent favoring or disfavoring folding is not an easy task. In this study, we present a computational protocol for estimating the free energy contributors to protein folding semi‐quantitatively using ensembles of unfolded and native states generated via molecular dynamics simulations. We tested the methodology on 35 proteins with diverse structural motifs and sizes and found that the calculated free energies correlate well with experiment (correlation coefficient ∼ 0.85), enabling us to develop a consensus view of the energetics of folding. As a more sensitive test of the methodology, we also investigated the free energies of folding of an additional 33 single point mutants and obtained a correlation coefficient of 0.8. A notable observation is that the folding free energy components appear to carry signatures of the fold (SCOP classification) of the protein. © 2017 Wiley Periodicals, Inc.     

The treatment of solvation by a generalized Born model and a self-consistent charge-density functional theory-based tight-binding method

We present a model to calculate the free energies of solvation of small organic compounds as well as large biomolecules. This model is based on a generalized Born (GB) model and a self-consistent charge-density functional theory-based tight-binding (SCC-DFTB) method with the nonelectrostatic contributions to the free energy of solvation modeled in terms of solvent-accessible surface areas (SA). The parametrization of the SCC-DFTB/GBSA model has been based on 60 neutral and six ionic molecules composed of H, C, N, O, and S, and spanning a wide range of chemical groups. Effective atomic radii as parameters have been obtained through Monte Carlo Simulated Annealing optimization in the parameter space to minimize the differences between the calculated and experimental free energies of solvation. The standard error in the free energies of solvation calculated by the final model is 1.11 kcal mol(-1). We also calculated the free energies of solvation for these molecules using a conductor-like screening model (COSMO) in combination with different levels of theory (AM1, SCC-DFTB, and B3LYP/6-31G*) and compared the results with SCC-DFTB/GBSA. To assess the efficiency of our model for large biomolecules, we calculated the free energy of solvation for a HIV protease-inhibitor complex containing 3,204 atoms using the SCC-DFTB/GBSA and the SCC-DFTB/COSMO models, separately. The computed relative free energies of solvation are comparable, while the SCC-DFTB/GBSA model is three to four times more efficient, in terms of computational cost.     

The SGB/NP hydration free energy model based on the surface generalized born solvent reaction field and novel nonpolar hydration free energy estimators

The development and parameterization of a solvent potential of mean force designed to reproduce the hydration thermodynamics of small molecules and macromolecules aimed toward applications in conformation prediction and ligand binding free energy prediction is presented. The model, named SGB/NP, is based on a parameterization of the Surface Generalized Born continuum dielectric electrostatic model using explicit solvent free energy perturbation calculations and a newly developed nonpolar hydration free energy estimator motivated by the results of explicit solvent simulations of the thermodynamics of hydration of hydrocarbons. The nonpolar model contains, in addition to the more commonly used solvent accessible surface area term, a component corresponding to the attractive solute-solvent interactions. This term is found to be important to improve the accuracy of the model, particularly for cyclic and hydrogen bonding compounds. The model is parameterized against the experimental hydration free energies of a set of small organic molecules. The model reproduces the experimental hydration free energies of small organic molecules with an accuracy comparable or superior to similar models employing more computationally demanding estimators and/or a more extensive set of parameters.     

Prediction of internal standards in reversed-phase liquid chromatography: IV. Correlation and prediction of retention in reversed-phase ion-pair chromatography based on linear solvation energy relationships

This paper describes the results of the evaluation of a new solvation parameter model for reversed-phase ion-pair chromatography by linear gradient elution. This model is described as . The first six terms are the usual solvation parameter equation for neutral solutes, and the seventh term represents the contribution to retention from solute’s ionization. The last term describes the retention increase due to ion-pair effect. Retention times obtained for 60 solutes (neutral, acidic and basic) in acetonitrile/aqueous mobile phases with different ion-pair reagents (phosphoric acid, trifluoroacetic acid, heptafluorobutyric acid, perchloric acid, and hexafluorophosphoric acid) are used to evaluate the capability of the function. It is concluded that the model describes the retention of ionizable/ionized compounds under ion-pair conditions very well. Accordingly, the function extends the application of linear solvation energy relationships (LSERs) to ionizable compounds in ion-pair chromatography, and allows us to easily predict their retention for chromatographic optimization, including selectivity optimization and internal standard selection. Finally, the conclusion can be extended to ioscratic elution.     

Force field evaluation for biomolecular simulation: free enthalpies of solvation of polar and apolar compounds in various solvents.

Recently, the GROMOS biomolecular force field parameter set 53A6--which has been parametrized to reproduce experimentally determined free enthalpies of hydration and solvation in cyclohexane of amino acid side-chain analogs--was presented. To investigate the transferability of the new parameter set, we calculated free enthalpies of solvation of a range of polar and apolar compounds in different solvents (methanol, dimethyl sulfoxide (DMSO), acetonitrile, and acetone) from molecular dynamics simulations using the GROMOS 53A6 force field. For methanol and DMSO, parameters were used that are available in the 53A6 parameter set. For acetonitrile, a recently developed model was taken and for acetone, two models available in literature were used. We found that trends in and values for the solvation free enthalpies are in satisfactory agreement with experiment, except for the solvation in acetone for which deviations from experiment can be explained in terms of the properties of the models used.     

A biomolecular force field based on the free enthalpy of hydration and solvation: the GROMOS force-field parameter sets 53A5 and 53A6

Successive parameterizations of the GROMOS force field have been used successfully to simulate biomolecular systems over a long period of time. The continuing expansion of computational power with time makes it possible to compute ever more properties for an increasing variety of molecular systems with greater precision. This has led to recurrent parameterizations of the GROMOS force field all aimed at achieving better agreement with experimental data. Here we report the results of the latest, extensive reparameterization of the GROMOS force field. In contrast to the parameterization of other biomolecular force fields, this parameterization of the GROMOS force field is based primarily on reproducing the free enthalpies of hydration and apolar solvation for a range of compounds. This approach was chosen because the relative free enthalpy of solvation between polar and apolar environments is a key property in many biomolecular processes of interest, such as protein folding, biomolecular association, membrane formation, and transport over membranes. The newest parameter sets, 53A5 and 53A6, were optimized by first fitting to reproduce the thermodynamic properties of pure liquids of a range of small polar molecules and the solvation free enthalpies of amino acid analogs in cyclohexane (53A5). The partial charges were then adjusted to reproduce the hydration free enthalpies in water (53A6). Both parameter sets are fully documented, and the differences between these and previous parameter sets are discussed.     

Application of reference‐modified density functional theory: Temperature and pressure dependences of solvation free energy

Recently, we proposed a reference‐modified density functional theory (RMDFT) to calculate solvation free energy (SFE), in which a hard‐sphere fluid was introduced as the reference system instead of an ideal molecular gas. Through the RMDFT, using an optimal diameter for the hard‐sphere reference system, the values of the SFE calculated at room temperature and normal pressure were in good agreement with those for more than 500 small organic molecules in water as determined by experiments. In this study, we present an application of the RMDFT for calculating the temperature and pressure dependences of the SFE for solute molecules in water. We demonstrate that the RMDFT has high predictive ability for the temperature and pressure dependences of the SFE for small solute molecules in water when the optimal reference hard‐sphere diameter determined for each thermodynamic condition is used. We also apply the RMDFT to investigate the temperature and pressure dependences of the thermodynamic stability of an artificial small protein, chignolin, and discuss the mechanism of high‐temperature and high‐pressure unfolding of the protein. © 2017 Wiley Periodicals, Inc.     

Solvation free energies calculated using the GB/SA model: Sensitivity of results on charge sets,protocols, and force fields

The sensitivity of aqueous solvation free energies (SFEs), estimated using the GB/SA continuum solvent model, on charge sets, protocols, and force fields, was studied. Simple energy calculations using the GB/SA solvent model were performed on 11 monofunctional organic compounds. Results indicate that calculated SFEs are strongly dependent on the charge sets. Charges derived from electrostatic potential fitting to high level ab initio wave functions using the CHELPG procedure and “class IV” charges from AM1/CM1a or PM3/CM1p calculations yielded better results than the corresponding Mulliken charges. Calculated SFEs were similar to MC/FEP energies obtained in the presence of explicit TIP4P water. Further improvements were obtained by using GVB/6-31G** and MP2/6-31+G** (CHELPG) charge sets that included correlation effects. SFEs calculated using charge sets assigned by the OPLSA* force field gave the best results of all standard force fields (MM2*, MM3*, MMFF, AMBER*, and OPLSA*) implemented in MacroModel. Comparison of relative and absolute SFEs computed using either the GB/SA continuum model or MC/FEP calculations in the presence of explicit TIP4P water showed that, in general, relative SFEs can be estimated with greater accuracy. A second set of 20 mono- and difunctional molecules was also studied and relative SFEs estimated using energy minimization and thermodynamic cycle perturbation (TCP) protocols. SFEs calculated from TCP calculations using the GB/SA model were sensitive to bond lengths of dummy bonds (i.e., bonds involving dummy atoms). In such cases, keeping the bond lengths of dummy bonds close to the corresponding bond lengths of the starting structures improved the agreement of TCP-calculated SFEs with energy minimization results. Overall, these results indicate that GB/SA solvation free energy estimates from simple energy minimization calculations are of similar accuracy and value to those obtained using more elaborate TCP protocols. © 1998 John Wiley & Sons, Inc. J Comput Chem 19: 769–780, 1998     

Using the local elevation method to construct optimized umbrella sampling potentials: Calculation of the relative free energies and interconversion barriers of glucopyranose ring conformers in water

A method is proposed to combine the local elevation (LE) conformational searching and the umbrella sampling (US) conformational sampling approaches into a single local elevation umbrella sampling (LEUS) scheme for (explicit‐solvent) molecular dynamics (MD) simulations. In this approach, an initial (relatively short) LE build‐up (searching) phase is used to construct an optimized biasing potential within a subspace of conformationally relevant degrees of freedom, that is then used in a (comparatively longer) US sampling phase. This scheme dramatically enhances (in comparison with plain MD) the sampling power of MD simulations, taking advantage of the fact that the preoptimized biasing potential represents a reasonable approximation to the negative of the free energy surface in the considered conformational subspace. The method is applied to the calculation of the relative free energies of β‐D ‐glucopyranose ring conformers in water (within the GROMOS 45A4 force field). Different schemes to assign sampled conformational regions to distinct states are also compared. This approach, which bears some analogies with adaptive umbrella sampling and metadynamics (but within a very distinct implementation), is shown to be: (i) efficient (nearly all the computational effort is invested in the actual sampling phase rather than in searching and equilibration); (ii) robust (the method is only weakly sensitive to the details of the build‐up protocol, even for relatively short build‐up times); (iii) versatile (a LEUS biasing potential database could easily be preoptimized for small molecules and assembled on a fragment basis for larger ones). © 2009 Wiley Periodicals, Inc. J Comput Chem 2010