Stability scale and atomic solvation parameters extracted from 1023 mutation experiments

The stability scale of 20 amino acid residues is derived from a database of 1023 mutation experiments on 35 proteins. The resulting scale of hydrophobic residues has an excellent correlation with the octanol-to-water transfer free energy corrected with an additional Flory-Huggins molar-volume term (correlation coefficient r = 0.95, slope = 1.05, and a near zero intercept). Thus, hydrophobic contribution to folding stability is characterized remarkably well by transfer experiments. However, no corresponding correlation is found for hydrophilic residues. Both the hydrophilic portion and the entire scale, however, correlate strongly with average burial accessible surface (r = 0.76 and 0.97, respectively). Such a strong correlation leads to a near uniform value of the atomic solvation parameters for atoms C, S, O/N, O(-0.5), and N(+0.5,1). All are in the range of 12-28 cal x mol(-1) A(-2), close to the original estimate of hydrophobic contribution of 25-30 cal x mol(-1) A(-2) to folding stability. Without any adjustable parameters, the new stability scale and new atomic solvation parameters yielded an accurate prediction of protein-protein binding free energy for a separate database of 21 protein-protein complexes (r = 0.80 and slope = 1.06, and r = 0.83 and slope = 0.93, respectively).     

Calculation of the free energy of association for protein complexes.

We have developed a method for calculating the association energy of quaternary complexes starting from their atomic coordinates. The association energy is described as the sum of two solvation terms and an energy term to account for the loss of translational and rotational entropy. The calculated solvation energy, using atomic solvation parameters and the solvent accessible surface areas, has a correlation of 96% with experimentally determined values. We have applied this methodology to examine intermediates in viral assembly and to assess the contribution isomerization makes to the association energy of molecular complexes. In addition, we have shown that the calculated association can be used as a predictive tool for analyzing modeled molecular complexes.     

A solvent model for simulations of peptides in bilayers. I. Membrane-promoting alpha-helix formation

We describe an efficient solvation model for proteins. In this model atomic solvation parameters imitating the hydrocarbon core of a membrane, water, and weak polar solvent (octanol) were developed. An optimal number of solvation parameters was chosen based on analysis of atomic hydrophobicities and fitting experimental free energies of gas-cyclohexane, gas-water, and octanol-water transfer for amino acids. The solvation energy term incorporated into the ECEPP/2 potential energy function was tested in Monte Carlo simulations of a number of small peptides with known energies of bilayer-water and octanol-water transfer. The calculated properties were shown to agree reasonably well with the experimental data. Furthermore, the solvation model was used to assess membrane-promoting alpha-helix formation. To accomplish this, all-atom models of 20-residue homopolypeptides-poly-Leu, poly-Val, poly-Ile, and poly-Gly in initial random coil conformation-were subjected to nonrestrained Monte Carlo conformational search in vacuo and with the solvation terms mimicking the water and hydrophobic parts of the bilayer. All the peptides demonstrated their largest helix-forming tendencies in a nonpolar environment, where the lowest-energy conformers of poly-Leu, Val, Ile revealed 100, 95, and 80% of alpha-helical content, respectively. Energetic and conformational properties of Gly in all environments were shown to be different from those observed for residues with hydrophobic side chains. Applications of the solvation model to simulations of peptides and proteins in the presence of membrane, along with limitations of the approach, are discussed.     

A procedure for detection and quantitation of cavity volumes proteins. Application to measure the strength of the hydrophobic driving force in protein folding

Accurate identification of cavities is important in the study of protein structure, stability, design, and ligand binding. Identification and quantitation of cavities is a nontrivial problem because most cavities are connected to the protein exterior. We describe a computational procedure for quantitating cavity volumes and apply this to derive an estimate of the hydrophobic driving force in protein folding. A grid-based Monte Carlo procedure is used to position water molecules on the surface of a protein. A Voronoi procedure is used to identify and quantitate empty space within the solvated protein. Additional cavities not detected by other existing procedures can be identified. Most of these are close to surface concavities. Residue volumes for both the interior and the surface residues as well as cavity volumes are in good agreement with volumes calculated from fully hydrated protein structures obtained from molecular dynamic simulations. We show that the loss of stability because of cavity-creating mutations correlates better with cavity volumes determined by this procedure than with cavity volumes determined by other methods. Available structural and thermodynamic data for a number of cavity-containing mutants were analyzed to obtain estimates of 26.1 cal x mol(-1) x A(-3) and 18.5 cal x mol(-1) x A(-2) for the relative contributions of cavity formation and the hydrophobic effect to the observed stability changes. The present estimate for the hydrophobic driving force is at the lower end of estimates derived from model compound studies and considerably lower than previous estimates of approximately 50 cal x mol(-1) x A(-2) derived from protein mutational data. In the absence of structural rearrangement, on average, deletion of a single methylene group is expected to result in losses in stability of 0.41 and 0.70 kcal x mol(-1) resulting from decrease in hydrophobicity and packing, respectively.     

A CFF 91-based Continuum Solvation Model: Solvation Free Energies of Small Organic Molecules and Conformations of the Alanine Dipeptide in Solution

Abstract

For molecular mechanics simulations of solvated molecules, it is important to use a consistent approach for calculating both the force field energy and the solvation free energy. A continuum solvation model based upon the atomic charges provided with the CFF91 force field is derived. The electrostatic component of the solvation free energy is described by the Poisson-Bolzmann equation while the nonpolar comonent of the solvation energy is assumed to be proportional to the solvent accessible surface area of the solute. Solute atomic radii used to describe the interface between the solute and solvent are fitted to reproduce the energies of small organic molecules. Data for 140 compounds are presented and compared to experiment and to the results from the well-characterized quantum mechanical solvation model AM1-SM2. In particular, accurate results are obtained for amino acid neutral analogues (mean unsigned error of 0.3 kcal/mol). The conformational energetics of the solvated alanine dipeptide is discussed.     

H-bonding in protein hydration revisited

H-bonding between protein surface polar/charged groups and water is one of the key factors of protein hydration. Here, we introduce an Accessible Surface Area (ASA) model for computationally efficient estimation of a free energy of water-protein H-bonding at any given protein conformation. The free energy of water-protein H-bonds is estimated using empirical formulas describing probabilities of hydrogen bond formation that were derived from molecular dynamics simulations of water molecules at the surface of a small protein, Crambin, from the Abyssinian cabbage (Crambe abyssinica) seed. The results suggest that atomic solvation parameters (ASP) widely used in continuum hydration models might be dependent on ASA for polar/charged atoms under consideration. The predictions of the model are found to be in qualitative agreement with the available experimental data on model compounds. This model combines the computational speed of ASA potential, with the high resolution of more sophisticated solvation methods.     

Extracting hydrophobic free energies from experimental data: relationship to protein folding and theoretical models

Solubility and vapor pressure measurements of hydrocarbons in water are generally thought to provide estimates of the strength of the hydrophobic effect in the range 20-30 cal/(mol.A2). Our reassessment of the solubility data on the basis of new developments in solution thermodynamics suggests that the hydrophobic surface free energy for hydrocarbon solutes is 46-47 cal/(mol.A2), although the actual value depends strongly on curvature effects [Nicholls et al. (1991) Proteins (in press); Sharp et al. (1991) Science 252, 106-109]. The arguments to support such a significant increase in the estimate of the hydrophobic effect stem partly from theoretical considerations and partly from the experimental results of De Young and Dill [(1990) J. Phys. Chem. 94, 801-809] on benzene partition between water and alkane solvents. Previous estimates of the hydrophobic effect derive from an analysis of solute partition data, which does not fully account for changes in volume entropy. We show here how the ideal gas equations, combined with experimental molar volumes, can account for such changes. Revised solubility scales for the 20 amino acids, based on cyclohexane to water and octanol to water transfer energies, are derived. The agreement between these scales, particularly the octanol scale, and mutant protein stability measurements from Kellis et al. [(1989) Biochemistry 28, 4914-4922] and Shortle et al. [(1990) Biochemistry 29, 8033-8041] is good. The increased strength of the hydrophobic interaction has implications for the energetics of protein folding, substrate binding, and nucleic acid base stacking and the interpretation of computer simulations.     

Energetics of protein folding

The energetics of protein folding determine the 3D structure of a folded protein. Knowledge of the energetics is needed to predict the 3D structure from the amino acid sequence or to modify the structure by protein engineering. Recent developments are discussed: major factors are reviewed and auxiliary factors are discussed briefly. Major factors include the hydrophobic factor (burial of non-polar surface area) and van der Waals interactions together with peptide hydrogen bonds and peptide solvation. The long-standing model for the hydrophobic factor (free energy change proportional to buried non-polar surface area) is contrasted with the packing-desolvation model and the approximate nature of the proportionality between free energy and apolar surface area is discussed. Recent energetic studies of forming peptide hydrogen bonds (gas phase) are reviewed together with studies of peptide solvation in solution. Closer agreement is achieved between the 1995 values for protein unfolding enthalpies in vacuum given by Lazaridis-Archontis-Karplus and Makhatadze-Privalov when the solvation enthalpy of the peptide group is taken from electrostatic calculations. Auxiliary factors in folding energetics include salt bridges and side-chain hydrogen bonds, disulfide bridges, and propensities to form alpha-helices and beta-structure. Backbone conformational entropy is a major energetic factor which is discussed only briefly for lack of knowledge.     

Empirical solvation models can be used to differentiate native from near-native conformations of bovine pancreatic trypsin inhibitor.

Several hydration models for peptides and proteins based on solvent accessible surface area have been proposed previously. We have evaluated some of these models as well as four new ones in the context of near-native conformations of a protein. In addition, we propose an empirical site-site distance-dependent correction that can be used in conjunction with any of these models. The set of near-native structures consisted of 39 conformations of bovine pancreatic trypsin inhibitor (BPTI) each of which was a local minimum of an empirical energy function (ECEPP) in the absence of solvent. Root-mean-square (rms) deviations from the crystallographically determined structure were in the following ranges: 1.06-1.94 A for all heavy atoms, 0.77-1.36 A for all backbone heavy atoms, 0.68-1.33 A for all alpha-carbon atoms, and 1.41-2.72 A for all side-chain heavy atoms. We have found that there is considerable variation among the solvent models when evaluated in terms of concordance between the solvation free energy and the rms deviations from the crystallographically determined conformation. The solvation model for which the best concordance (0.939) with the rms deviations of the C alpha atoms was found was derived from NMR coupling constants of peptides in water combined with an exponential site-site distance dependence of the potential of mean force. Our results indicate that solvation free energy parameters derived from nonpeptide free energies of hydration may not be transferrable to peptides. Parameters derived from peptide and protein data may be more applicable to conformational analysis of proteins. A general approach to derive parameters for free energy of hydration from ensemble-averaged properties of peptides in solution is described.     

Fast accurate evaluation of protein solvent exposure

Protein solvation energies are often taken to be proportional to solvent-accessible surface areas. Computation of these areas is numerically demanding and may become a bottleneck for folding and design applications. Fast graph-based methods, such as dead-end elimination (DEE), become possible if all energies, including solvation energies, are expressed as single-residue and pair-residue terms. To this end, Street and Mayo originated a pair-residue approximation for solvent-accessible surface areas (Street AG, Mayo SL. Pairwise calculation of protein solvent accessible surface areas. Fold Des 1998;3:253-258). The dominant source of error in this method is the overlapping burial of side-chain surfaces in the protein core. Here we report a new pair-residue approximation, which greatly reduces this overlap error by the use of optimized generic side-chains. We have tested the generic-side-chain method for the ten proteins studied by Street and Mayo and for 377 single-domain proteins from the CATH database (Orengo CA, Michie AD, Jones S, Jones DT, Swindells MB, Thornton JM. CATH-A hierarchic classification of protein domain structures. Structure 1997;5:1093-1108). With little additional cost in computation, the new method consistently reduces error for total areas and residue-by-residue areas by more than a factor of two. For example, the residue-by-residue error (for buried area) is reduced from 7.42 A(2) to 3.70 A(2). This difference translates into a solvation energy difference of approximately 0.2 kcal/mol per residue, amounting to a reduction in root-mean-square energy error of 2 kcal/mol for a 100 residue chain, a potentially critical difference for both protein folding and design applications.     

Evaluation of a fast implicit solvent model for molecular dynamics simulations.

A solvation term based on the solvent accessible surface area (SASA) is combined with the CHARMM polar hydrogen force field for the efficient simulation of peptides and small proteins in aqueous solution. Only two atomic solvation parameters are used: one is negative for favoring the direct solvation of polar groups and the other positive for taking into account the hydrophobic effect on apolar groups. To approximate the water screening effects on the intrasolute electrostatic interactions, a distance-dependent dielectric function is used and ionic side chains are neutralized. The use of an analytical approximation of the SASA renders the model extremely efficient (i.e., only about 50% slower than in vacuo simulations). The limitations and range of applicability of the SASA model are assessed by simulations of proteins and structured peptides. For the latter, the present study and results reported elsewhere show that with the SASA model it is possible to sample a significant amount of folding/unfolding transitions, which permit the study of the thermodynamics and kinetics of folding at an atomic level of detail.     

A microscopic view of peptide and protein solvation

The structure and dynamics of the water hydrating peptides and proteins are examined here at atomic resolution via molecular dynamics simulations. Detailed solvation density and residence time data for all 20 L-amino acids in an end-capped AXA tripeptide motif are presented. In addition, the solvation of the protein chymotrypsin inhibitor 2 is investigated as a point of comparison. Residues on the surface of proteins are not isolated; they interact both locally and non-locally in sequence space, and comparison of the solvation properties of each amino acid in both the peptide and protein allow us to distinguish inherent solvation properties from context-dependent perturbations due to neighboring residues. This work moves beyond traditional radial distribution functions and presents graphical representations of preferential solvation and orientation of water by side chains and the main chain. The combination of 0.3 micros of simulation data improves the statistical sampling over previous studies and reveals the significance of bridging water molecules that stabilize and mediate side chain-side chain, side chain-main chain and main chain-main chain interactions at the solvation interface.     

Comparison of atomic solvation parametric sets: applicability and limitations in protein folding and binding.

Atomic solvation parameters (ASP) are widely used to estimate the solvation contribution to the thermodynamic stability of proteins as well as the free energy of association for protein-ligand complexes. They are also included in several molecular mechanics computer programs. In this work, a total of eight atomic solvation parametric sets has been employed to calculate the solvation contribution to the free energy of folding delta Gs for 17 proteins. A linear correlation between delta Gs and the number of residues in each protein was found for each ASP set. The calculations also revealed a great variety in the absolute value and in the sign of delta Gs values such that certain ASP sets predicted the unfolded state to be more stable than the folded, whereas others yield precisely the opposite. Further, the solvation contribution to the free energy of association of helix pairs and to the disassociation of loops (connection between secondary structural elements in proteins) from the protein tertiary structures were computed for each of the eight ASP sets and discrepancies were evident among them.     

Modeling of peptides and proteins in membrane environment. I. A solvation model mimicking a lipid bilayer

A theoretical solvation model of peptides and proteins that mimics the heterogeneous membrane-water system was proposed. Our approach is based on the combined use of atomic parameters of solvation for water and hydrocarbons, which approximates the hydrated polar groups and acyl chains of lipids, respectively. This model was tested in simulations of several peptides: a nonpolar 20-mer polyleucine, a hydrophobic peptide with terminal polar groups, and a strongly amphiphilic peptide. The conformational space of the peptides in the presence of the membrane was studied by the Monte Carlo method. Unlike a polar solvent and vacuum, the membrane-like environment was shown to stabilize the alpha-helical conformation: low-energy structures have a helicity index of 100% in all cases. At the same time, the energetically most favorable orientations of the peptides relative to the membrane depend on their hydrophobic properties: nonpolar polyleucine is entirely immersed in the bilayer and the hydrophobic peptide with polar groups at the termini adopts a transbilayer orientation, whereas the amphiphilic peptide lies at the interface parallel to the membrane plane. The results of the simulations agree well with the available experimental data for these systems. In the following communications of this series, we plan to describe applications of the solvation model to membrane-bound proteins and peptides with biologically important functional activities.     

Solvation parameters for predicting the structure of surface loops in proteins: transferability and entropic effects

A new procedure for optimizing parameters of implicit solvation models introduced by us has been applied successfully first to cyclic peptides and more recently to three surface loops of ribonuclease A (Das and Meirovitch, Proteins 2001;43:303-314) using the simplified model E(tot) = E(FF)(epsilon = nr) + Sigma(i) sigma(i)A(i), where sigma(i) are atomic solvation parameters (ASPs) to be optimized, A(i) is the solvent accessible surface area of atom i, E(FF)(epsilon = nr) is the AMBER force-field energy of the loop-loop and loop-template interactions with a distance-dependent dielectric constant, epsilon = nr, where n is a parameter. The loop is free to move while the protein template is held fixed in its X-ray structure; an extensive conformational search for energy minimized loop structures is carried out with our local torsional deformation method. The optimal ASPs and n are those for which the structure with the lowest minimized energy [E(tot)(n,sigma(i))] becomes the experimental X-ray structure, or less strictly, the energy gap between these structures is within 2-3 kcal/mol. To check if a set of ASPs can be defined, which is transferable to a large number of loops, we optimize individual sets of ASPs (based on n = 2) for 12 surface loops from which an "averaged" best-fit set is defined. This set is then applied to the 12 loops and an independent "test" group of 8 loops leading in most cases to very small RMSD values; thus, this set can be useful for structure prediction of loops in homology modeling. For three loops we also calculate the free energy gaps to find that they are only slightly smaller than their energy counterparts, indicating that only larger n will enable reducing too large gaps. Because of its simplicity, this model allowed carrying out an extensive application of our methodology, providing thereby a large number of benchmark results for comparison with future calculations based on n > 2 as well as on more sophisticated solvation models with as yet unknown performance for loops.     

Thermodynamics of hydrogen bond and hydrophobic interactions in cyclodextrin complexes.

Values of K, delta G(o), delta H(o), delta S(o) and delta C(po) for the binding reaction of small organic ligands forming 1:1 complexes with either alpha- or beta-cyclodextrin were obtained by titration calorimetry from 15 degrees C to 45 degrees C. A hydrogen bond or hydrophobic interaction was introduced by adding a single functional group to the ligand. The thermodynamics of binding with and without the added group are compared to estimate the contribution of the hydrogen bond or hydrophobic interaction. A change in the environment of a functional group is required to influence the binding thermodynamics, but molecular size-dependent solute-solvent interactions have no effect. For phenolic O-H-O hydrogen bond formation, delta H(o) varies from -2 to -1.4 kcal mol(-1) from 15 degrees C to 45 degrees C, and delta C(p) is increased by 18 cal K(-1) mol(-1). The hydrophobic interaction has an opposite effect: in alpha-cyclodextrin, delta C(po) = -13.3 cal K(-1) mol(-1) per ligand -CH(2)-, identical to values found for the transfer of a -CH(2)-group from water to a nonpolar environment. At room temperature, the hydrogen bond and the -CH(2)-interaction each contribute about -600 cal mol(-1) to the stability (delta G(o)) of the complex. With increased temperature, the hydrogen bond stability decreases (i.e., hydrogen bonds "melt"), but the stability of the hydrophobic interaction remains essentially constant.     

Discrimination of the native from misfolded protein models with an energy function including implicit solvation.

An essential requirement for theoretical protein structure prediction is an energy function that can discriminate the native from non-native protein conformations. To date most of the energy functions used for this purpose have been extracted from a statistical analysis of the protein structure database, without explicit reference to the physical interactions responsible for protein stability. The use of the statistical functions has been supported by the widespread belief that they are superior for such discrimination to physics-based energy functions. An effective energy function which combined the CHARMM vacuum potential with a Gaussian model for the solvation free energy is tested for its ability to discriminate the native structure of a protein from misfolded conformations; the results are compared with those obtained with the vacuum CHARMM potential. The test is performed on several sets of misfolded structures prepared by others, including sets of about 650 good decoys for six proteins, as well as on misfolded structures of chymotrypsin inhibitor 2. The vacuum CHARMM potential is successful in most cases when energy minimized conformations are considered, but fails when applied to structures relaxed by molecular dynamics. With the effective energy function the native state is always more stable than grossly misfolded conformations both in energy minimized and molecular dynamics-relaxed structures. The present results suggest that molecular mechanics (physics-based) energy functions, complemented by a simple model for the solvation free energy, should be tested for use in the inverse folding problem, and supports their use in studies of the effective energy surface of proteins in solution. Moreover, the study suggests that the belief in the superiority of statistical functions for these purposes may be ill founded.     

Calorimetric determination of the energetics of the molten globule intermediate in protein folding: apo-alpha-lactalbumin

High-sensitivity differential scanning calorimetry has been used to characterize the energetics of the molten globule state of apo-alpha-lactalbumin. This characterization has been possible by performing temperature scans at different guanidine hydrochloride (GuHCl) concentrations in order to experimentally define the temperature-GuHCl stability surface of the protein. Multidimensional analysis of the heat capacity surface has allowed simultaneous resolution of the energetics of the unfolded and molten globule states. These experiments indicate that the intrinsic enthalpy difference (i.e., excluding additional contributions such as those arising from differential GuHCl binding) between the unfolded and native states is 31.8 kcal/mol at 25 degrees C whereas that of the molten globule and native states is only 7.7 kcal/mol. At the same temperature, the entropy changes are 99.2 and 23.7 cal/K.mol and the heat capacity changes are 1821 and 326 cal/K.mol, respectively. Analysis of the thermodynamic data indicates that in passing from the native to the molten globule state only approximately 19% of the hydrogen bonds are broken. In addition, the magnitude of delta Cp for the molten globule suggests that water does not largely penetrate into the interior of the molten globule, implying that significant hydrophobic interactions are still present in this state. These parameters provide precise energetic constraints to the allowed structural conformations of the molten globule.     

Calculations of electrostatic energies in proteins. The energetics of ionized groups in bovine pancreatic trypsin inhibitor

The importance of including different energy contributions in calculations of electrostatic energies in proteins is examined by calculating the intrinsic pKa values of the acidic groups of bovine pancreatic trypsin inhibitor. It appears that such calculations provide a powerful and revealing test; the relevant solvation energies of the ionized acids are of the order of -70 kcal/mol (1 cal = 4.184 J), and microscopic calculations that do not attempt to simulate the complete protein dielectric effect (including the surrounding solvent) can underestimate the solvation energy by as much as 50 kcal/mol. Reproducing correctly, by the same set of parameters, the solvation energies of ionized acids in different sites of a protein cannot be accomplished by including only part of the key energy contributions. The problems associated with macroscopic calculations are also considered and illustrated by the specific case of bovine pancreatic trypsin inhibitor. A promising approach is shown to be provided by a refinement of the previously developed Protein Dipoles Langevin Dipoles model. This model seems to represent consistently the microscopic dielectric of the protein and the surrounding water molecules. The model overcomes the problems associated with the macroscopic models (by treating explicitly the solvent molecules) and avoids the convergence problems associated with all-atom solvent models (by treating the average solvent polarization rather than averaging the actual polarization energy). This paper describes in detail the actual implementation of the model and examines its performance in evaluating intrinsic pKa values. Preliminary microscopic considerations of charge-charge interactions are presented.