An improved protein decoy set for testing energy functions for protein structure prediction

We have improved the original Rosetta centroid/backbone decoy set by increasing the number of proteins and frequency of near native models and by building on sidechains and minimizing clashes. The new set consists of 1,400 model structures for 78 different and diverse protein targets and provides a challenging set for the testing and evaluation of scoring functions. We evaluated the extent to which a variety of all-atom energy functions could identify the native and close-to-native structures in the new decoy sets. Of various implicit solvent models, we found that a solvent-accessible surface area-based solvation provided the best enrichment and discrimination of close-to-native decoys. The combination of this solvation treatment with Lennard Jones terms and the original Rosetta energy provided better enrichment and discrimination than any of the individual terms. The results also highlight the differences in accuracy of NMR and X-ray crystal structures: a large energy gap was observed between native and non-native conformations for X-ray structures but not for NMR structures.     

Recent advances in implicit solvent-based methods for biomolecular simulations

Implicit solvent-based methods play an increasingly important role in molecular modeling of biomolecular structure and dynamics. Recent methodological developments have mainly focused on the extension of the generalized Born (GB) formalism for variable dielectric environments and accurate treatment of nonpolar solvation. Extensive efforts in parameterization of GB models and implicit solvent force fields have enabled ab initio simulation of protein folding to native or near-native structures. Another exciting area that has benefited from the advances in implicit solvent models is the development of constant pH molecular dynamics methods, which have recently been applied to the calculations of protein pK(a) values and the studies of pH-dependent peptide and protein folding.     

Protein design simulations suggest that side-chain conformational entropy is not a strong determinant of amino acid environmental preferences

Loss of side-chain conformational entropy is an important force opposing protein folding and the relative preferences of the amino acids for being buried or solvent exposed may be partially determined by which amino acids lose more side-chain entropy when placed in the core of a protein. To investigate these preferences, we have incorporated explicit modeling of side-chain entropy into the protein design algorithm, RosettaDesign. In the standard version of the program, the energy of a particular sequence for a fixed backbone depends only on the lowest energy side-chain conformations that can be identified for that sequence. In the new model, the free energy of a single amino acid sequence is calculated by evaluating the average energy and entropy of an ensemble of structures generated by Monte Carlo sampling of amino acid side-chain conformations. To evaluate the impact of including explicit side-chain entropy, sequences were designed for 110 native protein backbones with and without the entropy model. In general, the differences between the two sets of sequences are modest, with the largest changes being observed for the longer amino acids: methionine and arginine. Overall, the identity between the designed sequences and the native sequences does not increase with the addition of entropy, unlike what is observed when other key terms are added to the model (hydrogen bonding, Lennard-Jones energies, and solvation energies). These results suggest that side-chain conformational entropy has a relatively small role in determining the preferred amino acid at each residue position in a protein.     

Free energy landscape of protein folding in water: explicit vs. implicit solvent

The Generalized Born (GB) continuum solvent model is arguably the most widely used implicit solvent model in protein folding and protein structure prediction simulations; however, it still remains an open question on how well the model behaves in these large-scale simulations. The current study uses the beta-hairpin from C-terminus of protein G as an example to explore the folding free energy landscape with various GB models, and the results are compared to the explicit solvent simulations and experiments. All free energy landscapes are obtained from extensive conformation space sampling with a highly parallel replica exchange method. Because solvation model parameters are strongly coupled with force fields, five different force field/solvation model combinations are examined and compared in this study, namely the explicit solvent model: OPLSAA/SPC model, and the implicit solvent models: OPLSAA/SGB (Surface GB), AMBER94/GBSA (GB with Solvent Accessible Surface Area), AMBER96/GBSA, and AMBER99/GBSA. Surprisingly, we find that the free energy landscapes from implicit solvent models are quite different from that of the explicit solvent model. Except for AMBER96/GBSA, all other implicit solvent models find the lowest free energy state not the native state. All implicit solvent models show erroneous salt-bridge effects between charged residues, particularly in OPLSAA/SGB model, where the overly strong salt-bridge effect results in an overweighting of a non-native structure with one hydrophobic residue F52 expelled from the hydrophobic core in order to make better salt bridges. On the other hand, both AMBER94/GBSA and AMBER99/GBSA models turn the beta-hairpin in to an alpha-helix, and the alpha-helical content is much higher than the previously reported alpha-helices in an explicit solvent simulation with AMBER94 (AMBER94/TIP3P). Only AMBER96/GBSA shows a reasonable free energy landscape with the lowest free energy structure the native one despite an erroneous salt-bridge between D47 and K50. Detailed results on free energy contour maps, lowest free energy structures, distribution of native contacts, alpha-helical content during the folding process, NOE comparison with NMR, and temperature dependences are reported and discussed for all five models.     

Protein Structure Prediction and Design in a Biologically Realistic Implicit Membrane

Protein design is a powerful tool for elucidating mechanisms of function and engineering new therapeutics and nanotechnologies. Although soluble protein design has advanced, membrane protein design remains challenging because of difficulties in modeling the lipid bilayer. In this work, we developed an implicit approach that captures the anisotropic structure, shape of water-filled pores, and nanoscale dimensions of membranes with different lipid compositions. The model improves performance in computational benchmarks against experimental targets, including prediction of protein orientations in the bilayer, ΔΔG calculations, native structure discrimination, and native sequence recovery. When applied to de novo protein design, this approach designs sequences with an amino acid distribution near the native amino acid distribution in membrane proteins, overcoming a critical flaw in previous membrane models that were prone to generating leucine-rich designs. Furthermore, the proteins designed in the new membrane model exhibit native-like features including interfacial aromatic side chains, hydrophobic lengths compatible with bilayer thickness, and polar pores. Our method advances high-resolution membrane protein structure prediction and design toward tackling key biological questions and engineering challenges.     

Predicting peptide binding to MHC pockets via molecular modeling, implicit solvation, and global optimization

Development of a computational prediction method based on molecular modeling, global optimization, and implicit solvation has produced accurate structure and relative binding affinity predictions for peptide amino acids binding to five pockets of the MHC molecule HLA-DRB1*0101. Because peptide binding to MHC molecules is essential to many immune responses, development of such a method for understanding and predicting the forces that drive binding is crucial for pharmaceutical design and disease treatment. Underlying the development of this prediction method are two hypotheses. The first is that pockets formed by the peptide binding groove of MHC molecules are independent, separating the prediction of peptide amino acids that bind within individual pockets from those that bind between pockets. The second hypothesis is that the native state of a system composed of an amino acid bound to a protein pocket corresponds to the system's lowest free energy. The prediction method developed from these hypotheses uses atomistic-level modeling, deterministic global optimization, and three methods of implicit solvation: solvent-accessible area, solvent-accessible volume, and Poisson-Boltzmann electrostatics. The method predicts relative binding affinities of peptide amino acids for pockets of HLA-DRB1*0101 by determining computationally an amino acid's global minimum energy conformation. Prediction results from the method are in agreement with X-ray crystallography data and experimental binding assays.     

Implicit solvent models for flexible protein-protein docking by molecular dynamics simulation

The suitability of three implicit solvent models for flexible protein-protein docking by procedures using molecular dynamics simulation is investigated. The three models are (i) the generalized Born (GB) model implemented in the program AMBER6.0; (ii) a distance-dependent dielectric (DDD) model; and (iii) a surface area-dependent model that we have parameterized and call the NPSA model. This is a distance-dependent dielectric model modified by neutralizing the ionizable side-chains and adding a surface area-dependent solvation term. These solvent models were first tested in molecular dynamics simulations at 300 K of the native structures of barnase, barstar, segment B1 of protein G, and three WW domains. These protein structures display a range of secondary structure contents and stabilities. Then, to investigate the performance of the implicit solvent models in protein docking, molecular dynamics simulations of barnase/barstar complexation, as well as PIN1 WW domain/peptide complexation, were conducted, starting from separated unbound structures. The simulations show that the NPSA model has significant advantages over the DDD and GB models in maintaining the native structures of the proteins and providing more accurate docked complexes.     

Comparison of binding energies of SrcSH2-phosphotyrosyl peptides with structure-based prediction using surface area based empirical parameterization

The prediction of binding energies from the three-dimensional (3D) structure of a protein-ligand complex is an important goal of biophysics and structural biology. Here, we critically assess the use of empirical, solvent-accessible surface area-based calculations for the prediction of the binding of Src-SH2 domain with a series of tyrosyl phosphopeptides based on the high-affinity ligand from the hamster middle T antigen (hmT), where the residue in the pY+ 3 position has been changed. Two other peptides based on the C-terminal regulatory site of the Src protein and the platelet-derived growth factor receptor (PDGFR) are also investigated. Here, we take into account the effects of proton linkage on binding, and test five different surface area-based models that include different treatments for the contributions to conformational change and protein solvation. These differences relate to the treatment of conformational flexibility in the peptide ligand and the inclusion of proximal ordered solvent molecules in the surface area calculations. This allowed the calculation of a range of thermodynamic state functions (deltaCp, deltaS, deltaH, and deltaG) directly from structure. Comparison with the experimentally derived data shows little agreement for the interaction of SrcSH2 domain and the range of tyrosyl phosphopeptides. Furthermore, the adoption of the different models to treat conformational change and solvation has a dramatic effect on the calculated thermodynamic functions, making the predicted binding energies highly model dependent. While empirical, solvent-accessible surface area based calculations are becoming widely adopted to interpret thermodynamic data, this study highlights potential problems with application and interpretation of this type of approach. There is undoubtedly some agreement between predicted and experimentally determined thermodynamic parameters: however, the tolerance of this approach is not sufficient to make it ubiquitously applicable.     

Amino acid encoding schemes from protein structure alignments: multi-dimensional vectors to describe residue types

Bioinformatic software has used various numerical encoding schemes to describe amino acid sequences. Orthogonal encoding, employing 20 numbers to describe the amino acid type of one protein residue, is often used with artificial neural network (ANN) models. However, this can increase the model complexity, thus leading to difficulty in implementation and poor performance. Here, we use ANNs to derive encoding schemes for the amino acid types from protein three-dimensional structure alignments. Each of the 20 amino acid types is characterized with a few real numbers. Our schemes are tested on the simulation of amino acid substitution matrices. These simplified schemes outperform the orthogonal encoding on small data sets. Using one of these encoding schemes, we generate a colouring scheme for the amino acids in which comparable amino acids are in similar colours. We expect it to be useful for visual inspection and manual editing of protein multiple sequence alignments.     

Modeling loop reorganization free energies of acetylcholinesterase: a comparison of explicit and implicit solvent models

The treatment of hydration effects in protein dynamics simulations varies in model complexity and spans the range from the computationally intensive microscopic evaluation to simple dielectric screening of charge-charge interactions. This paper compares different solvent models applied to the problem of estimating the free-energy difference between two loop conformations in acetylcholinesterase. Molecular dynamics (MD) simulations were used to sample potential energy surfaces of the two basins with solvent treated by means of explicit and implicit methods. Implicit solvent methods studied include the generalized Born (GB) model, atomic solvation potential (ASP), and the distance-dependent dieletric constant. By using the linear response approximation (LRA), the explicit solvent calculations determined a free-energy difference that is in excellent agreement with the experimental estimate, while rescoring the protein conformations with GB or the Poisson equation showed inconsistent and inferior results. While the approach of rescoring conformations from explicit water simulations with implicit solvent models is popular among many applications, it perturbs the energy landscape by changing the solvent contribution to microstates without conformational relaxation, thus leading to non-optimal solvation free energies. Calculations applying MD with a GB solvent model produced results of comparable accuracy as observed with LRA, yet the electrostatic free-energy terms were significantly different due to optimization on a potential energy surface favored by an implicit solvent reaction field. The simpler methods of ASP and the distance-dependent scaling of the dielectric constant both produced considerable distortions in the protein internal free-energy terms and are consequently unreliable.     

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.     

Free energies of protein decoys provide insight into determinants of protein stability

We have calculated the stability of decoy structures of several proteins (from the CASP3 models and the Park and Levitt decoy set) relative to the native structures. The calculations were performed with the force field-consistent ES/IS method, in which an implicit solvent (IS) model is used to calculate the average solvation free energy for snapshots from explicit simulations (ESs). The conformational free energy is obtained by adding the internal energy of the solute from the ESs and an entropic term estimated from the covariance positional fluctuation matrix. The set of atomic Born radii and the cavity-surface free energy coefficient used in the implicit model has been optimized to be consistent with the all-atom force field used in the ESs (cedar/gromos with simple point charge (SPC) water model). The decoys are found to have a consistently higher free energy than that of the native structure; the gap between the native structure and the best decoy varies between 10 and 15 kcal/mole, on the order of the free energy difference that typically separates the native state of a protein from the unfolded state. The correlation between the free energy and the extent to which the decoy structures differ from the native (as root mean square deviation) is very weak; hence, the free energy is not an accurate measure for ranking the structurally most native-like structures from among a set of models. Analysis of the energy components shows that stability is attained as a result of three major driving forces: (1) minimum size of the protein-water surface interface; (2) minimum total electrostatic energy, which includes solvent polarization; and (3) minimum protein packing energy. The detailed fit required to optimize the last term may underlie difficulties encountered in recovering the native fold from an approximate decoy or model structure.     

A combined steepest descent and genetic algorithm (SD/GA) approach for the optimization of solvation parameters

An accurate solvation model is essential for computer modeling of protein folding and other biomolecular self-assembly processes. Compared to explicit solvent models, implicit solvent models, such as the Poisson-Boltzmann (PB) with solvent accessible surface area model (PB/SA), offer a much faster speed—the most compelling reason for the popularity of these implicit solvent models. Since these implicit solvent models typically use empirical parameters, such as atomic radii and the surface tensions, an optimal fit of these parameters is crucial for the final accuracy of properties such as solvation free energy and folding free energy. In this paper, we proposed a combined approach, namely SD/GA, which takes the advantage of both local optimization with the steepest descent (SD), and global optimization with the genetic algorithm (GA), for parameters optimization in multi-dimensional space. The SD/GA method is then applied to the optimization of solvation parameters in the non-polar cavity term of the PB/SA model. The results show that the newly optimized parameters from SD/GA not only increase the accuracy in the solvation free energies for ～200 organic molecules, but also significantly improve the free energy landscape of a β-hairpin folding. The current SD/GA method can be readily applied to other multi-dimensional parameter space optimization as well.     

A "solvated rotamer" approach to modeling water-mediated hydrogen bonds at protein-protein interfaces

Water-mediated hydrogen bonds play critical roles at protein-protein and protein-nucleic acid interfaces, and the interactions formed by discrete water molecules cannot be captured using continuum solvent models. We describe a simple model for the energetics of water-mediated hydrogen bonds, and show that, together with knowledge of the positions of buried water molecules observed in X-ray crystal structures, the model improves the prediction of free-energy changes upon mutation at protein-protein interfaces, and the recovery of native amino acid sequences in protein interface design calculations. We then describe a "solvated rotamer" approach to efficiently predict the positions of water molecules, at protein-protein interfaces and in monomeric proteins, that is compatible with widely used rotamer-based side-chain packing and protein design algorithms. Finally, we examine the extent to which the predicted water molecules can be used to improve prediction of amino acid identities and protein-protein interface stability, and discuss avenues for overcoming current limitations of the approach.     

The Limited Role of Nonnative Contacts in the Folding Pathways of a Lattice Protein

Models of protein energetics that neglect interactions between amino acids that are not adjacent in the native state, such as the Gō model, encode or underlie many influential ideas on protein folding. Implicit in this simplification is a crucial assumption that has never been critically evaluated in a broad context: Detailed mechanisms of protein folding are not biased by nonnative contacts, typically argued to be a consequence of sequence design and/or topology. Here we present, using computer simulations of a well-studied lattice heteropolymer model, the first systematic test of this oft-assumed correspondence over the statistically significant range of hundreds of thousands of amino acid sequences that fold to the same native structure. Contrary to previous conjectures, we find a multiplicity of folding mechanisms, suggesting that Gō-like models cannot be justified by considerations of topology alone. Instead, we find that the crucial factor in discriminating among topological pathways is the heterogeneity of native contact energies: The order in which native contacts accumulate is profoundly insensitive to omission of nonnative interactions, provided that native contact heterogeneity is retained. This robustness holds over a surprisingly wide range of folding rates for our designed sequences. Mirroring predictions based on the principle of minimum frustration, fast-folding sequences match their Gō-like counterparts in both topological mechanism and transit times. Less optimized sequences dwell much longer in the unfolded state and/or off-pathway intermediates than do Gō-like models. For dynamics that bridge unfolded and unfolded states, however, even slow folders exhibit topological mechanisms and transit times nearly identical with those of their Gō-like counterparts. Our results do not imply a direct correspondence between folding trajectories of Gō-like models and those of real proteins, but they do help to clarify key topological and energetic assumptions that are commonly used to justify such caricatures.     

Water as a conformational editor in protein folding

As molecules approach one another in aqueous solution, desolvation free energy barriers to association are encountered. Experiments suggest these (de)solvation effects contribute to the free energy barriers separating the folded and unfolded states of protein molecules. To explore their influence on the energy landscapes of protein folding reactions, we have incorporated desolvation barriers into a semi-realistic, off-lattice protein model that uses a simplified physico-chemical force-field determined solely by the sequence of amino acids. Monte Carlo sampling techniques were used to study the effects on the thermodynamics and kinetics of folding of a number of systems, diverse in structure and sequence. In each case, desolvation barriers increase the stability of the native conformation and the cooperativity of the major folding/unfolding transition. The folding times of these systems are reduced significantly upon inclusion of desolvation barriers, demonstrating that the particulate nature of the solvent engenders a more defined route to the native fold.     

A simple physical model for the prediction and design of protein-DNA interactions

Protein-DNA interactions are crucial for many biological processes. Attempts to model these interactions have generally taken the form of amino acid-base recognition codes or purely sequence-based profile methods, which depend on the availability of extensive sequence and structural information for specific structural families, neglect side-chain conformational variability, and lack generality beyond the structural family used to train the model. Here, we take advantage of recent advances in rotamer-based protein design and the large number of structurally characterized protein-DNA complexes to develop and parameterize a simple physical model for protein-DNA interactions. The model shows considerable promise for redesigning amino acids at protein-DNA interfaces, as design calculations recover the amino acid residue identities and conformations at these interfaces with accuracies comparable to sequence recovery in globular proteins. The model shows promise also for predicting DNA-binding specificity for fixed protein sequences: native DNA sequences are selected correctly from pools of competing DNA substrates; however, incorporation of backbone movement will likely be required to improve performance in homology modeling applications. Interestingly, optimization of zinc finger protein amino acid sequences for high-affinity binding to specific DNA sequences results in proteins with little or no predicted specificity, suggesting that naturally occurring DNA-binding proteins are optimized for specificity rather than affinity. When combined with algorithms that optimize specificity directly, the simple computational model developed here should be useful for the engineering of proteins with novel DNA-binding specificities.     

Density functional calculations of 15N chemical shifts in solvated dipeptides

We performed density functional calculations to examine the effects of solvation, hydrogen bonding, backbone conformation, and the side chain on 15N chemical shielding in proteins. We used N-methylacetamide (NMA) and N-formyl-alanyl-X (with X being one of the 19 naturally occurring amino acids excluding proline) as model systems. In addition, calculations were performed for selected fragments from protein GB3. The conducting polarizable continuum model was employed to include the effect of solvent in the density functional calculations. Our calculations for NMA show that the augmentation of the polarizable continuum model with the explicit water molecules in the first solvation shell has a significant influence on isotropic 15N chemical shift but not as much on the chemical shift anisotropy. The difference in the isotropic chemical shift between the standard beta-sheet and alpha-helical conformations ranges from 0.8 to 6.2 ppm depending on the residue type, with the mean of 2.7 ppm. This is in good agreement with the experimental chemical shifts averaged over a database of 36 proteins containing >6100 amino acid residues. The orientation of the 15N chemical shielding tensor as well as its anisotropy and asymmetry are also in the range of values experimentally observed for peptides and proteins.     

Origins of Pressure-Induced Protein Transitions

The molecular mechanisms underlying pressure-induced protein denaturation can be analyzed based on the pressure-dependent differences in the apparent volume occupied by amino acids inside the protein and when they are exposed to water in an unfolded conformation. We present here an analysis for the peptide group and the 20 naturally occurring amino acid side chains based on volumetric parameters for the amino acids in the interior of the native state, the micelle-like interior of the pressure-induced denatured state, and the unfolded conformation modeled by N-acetyl amino acid amides. The transfer of peptide groups from the protein interior to water becomes increasingly favorable as pressure increases. Thus, solvation of peptide groups represents a major driving force in pressure-induced protein denaturation. Polar side chains do not appear to exhibit significant pressure-dependent changes in their preference for the protein interior or solvent. The transfer of nonpolar side chains from the protein interior to water becomes more unfavorable as pressure increases. We conclude that a sizeable population of nonpolar side chains remains buried inside a solvent-inaccessible core of the pressure-induced denatured state. At elevated pressures, this core may become packed almost as tightly as the interior of the native state. The presence and partial disappearance of large intraglobular voids is another driving force facilitating pressure-induced denaturation of individual proteins. Our data also have implications for the kinetics of protein folding and shed light on the nature of the folding transition state ensemble.