Thermodynamics of amide hydrogen bond formation in polar and apolar solvents

Molecular Simulations On Supercomputers

Work from our laboratory using molecular dynamics applied to the study of peptide/protein folding dynam ics and thermodynamics in aqueous solution, and the exploration and development of molecular solvent in teraction models including electronic polarizability is described. The basic computational techniques we have used in developing algorithms for vector/parallel, shared memory, multiple instruction, multiple data (MIMD) and massively parallel single instruction, multi ple data (SIMD) architectures are discussed in the con text of these application areas. The technique we have used to parallelize the molecular dynamics/molecular mechanics code CHARMM is detailed, with special fo cus given to the performance issues faced in develop ing algorithms for the CRAY Y-MP. New developments of a molecular dynamics kernel for simulation of elec tronically polarizable solvent systems using a parallel- extended Fortran language (Fortran-9X) on the mas sively parallel Connection Machine (CM-2) are out lined. We provide examples of applications using both of these platforms and discuss the development of software for such simulations in a heterogeneous com putational environment.

Exploring Multiple Binding Modes Using Confined Replica Exchange Molecular Dynamics

Molecular docking is extensively applied to determine the position of a ligand on its receptor despite the rather poor correspondence between docking scores and experimental binding affinities found in several studies, especially for systems structurally unrelated with those used in the scoring functions' training sets. Here, we present a method for the prediction of binding modes and binding free energies, which uses replica exchange molecular dynamics in combination with a receptor-shaped piecewise potential, confining the ligand in the proximity of the receptor surface and limiting the accessible conformational space of interest. We assess our methodology with a set of protein receptor-ligand test cases. In every case studied, the method is able to locate the ligand on the experimentally known receptor binding site, and it gives as output the binding free energy. The added value of our approach with respect to other available methods is that it quickly performs a conformational space search, providing a set of bound (or unbound) configurations, which can be used to determine phenomenological structural and energetic properties of an experimental binding state as a result of contributions provided by diversified multiple binding poses.

Evaluations of the Absolute and Relative Free Energies for Antidepressant Binding to the Amino Acid Membrane Transporter LeuT with Free Energy Simulations

The binding of ligands to protein receptors with high affinity and specificity is central to many cellular processes. The quest for the development of computational models capable of accurately evaluating binding affinity remains one of the main goals of modern computational biophysics. In this work, free energy perturbation/molecular dynamics simulations were used to evaluate absolute and relative binding affinity for three different antidepressants to a sodium-dependent membrane transporter, LeuT, a bacterial homologue of human serotonin and dopamine transporters. Dysfunction of these membrane transporters in mammals has been implicated in multiple diseases of the nervous system, including bipolar disorder and depression. Furthermore, these proteins are key targets for antidepressants including fluoxetine (aka Prozac) and tricyclic antidepressants known to block transport activity. In addition to being clinically relevant, this system, where multiple crystal structures are readily available, represents an ideal testing ground for methods used to study the molecular mechanisms of ligand binding to membrane proteins. We discuss possible pitfalls and different levels of approximation required to evaluate binding affinity, such as the dependence of the computed affinities on the strength of constraints and the sensitivity of the computed affinities to the particular partial charges derived from restrained electrostatic potential fitting of quantum mechanics electrostatic potential. Finally, we compare the effects of different constraint schemes on the absolute and relative binding affinities obtained from free energy simulations.

Strength of a bifurcated H bond

Macromolecules are characterized by their particular arrangement of H bonds. Many of these interactions involve a single donor and acceptor pair, such as the regular H-bonding pattern between carbonyl oxygens and amide H(+)s four residues apart in α-helices. The H-bonding potential of some acceptors, however, leads to the phenomenon of overcoordination between two donors and one acceptor. Herein, using isotope-edited Fourier transform infrared measurements and density functional theory (DFT) calculations, we measured the strength of such bifurcated H bonds in a transmembrane α-helix. Frequency shifts of the (13)C=(18)O amide I mode were used as a reporter of the strength of the bifurcated H bond from a thiol and hydroxyl H(+) at residue i + 4. DFT calculations yielded very similar frequency shifts and an energy of -2.6 and -3.4 kcal/mol for the thiol and hydroxyl bifurcated H bonds, respectively. The strength of the intrahelical bifurcated H bond is consistent with its prevalence in hydrophobic environments and is shown to significantly impact side-chain rotamer distribution.

Structural stability of proteins in aqueous and nonpolar environments

A protein folds into its native structure with the α-helix and∕or β-sheet in aqueous solution under the physiological condition. The relative content of these secondary structures largely varies from protein to protein. However, such structural variability is not exhibited in nonaqueous environment. For example, there is a strong trend that alcohol induces a protein to form α-helices, and many of the membrane proteins within the lipid bilayer consists of α-helices. Here we investigate the structural stability of proteins in aqueous and nonpolar environments using our recently developed free-energy function F = (Λ - TS)∕(k(B)T(0)) = Λ∕(k(B)T(0)) - S∕k(B) (T(0) = 298 K and the absolute temperature T is set at T(0)) which is based on statistical thermodynamics. Λ∕(k(B)T(0)) and S∕k(B) are the energetic and entropic components, respectively, and k(B) is Boltzmann's constant. A smaller value of the positive quantity, -S, represents higher efficiency of the backbone and side-chain packing promoted by the entropic effect arising from the translational displacement of solvent molecules or the CH(2), CH(3), and CH groups which constitute nonpolar chains of lipid molecules. As for Λ, in aqueous solution, a transition to a more compact structure of a protein accompanies the break of protein-solvent hydrogen bonds: As the number of donors and acceptors buried without protein intramolecular hydrogen bonding increases, Λ becomes higher. In nonpolar solvent, lower Λ simply implies more intramolecular hydrogen bonds formed. We find the following. The α-helix and β-sheet are advantageous with respect to -S as well as Λ and to be formed as much as possible. In aqueous solution, the solvent-entropy effect on the structural stability is so strong that the close packing of side chains is dominantly important, and the α-helix and β-sheet contents are judiciously adjusted to accomplish it. In nonpolar solvent, the solvent-entropy effect is substantially weaker than in aqueous solution. Λ is crucial and the α-helix is more stable than the β-sheet in terms of Λ, which develops a tendency that α-helices are exclusively chosen. For a membrane protein, α-helices are stabilized as fundamental structural units for the same reason, but their arrangement is performed through the entropic effect mentioned above.

Study of base pair mutations in proline-rich homeodomain (PRH)-DNA complexes using molecular dynamics

Proline-rich homeodomain (PRH) is a regulatory protein controlling transcription and gene expression processes by binding to the specific sequence of DNA, especially to the sequence 5'-TAATNN-3'. The impact of base pair mutations on the binding between the PRH protein and DNA is investigated using molecular dynamics and free energy simulations to identify DNA sequences that form stable complexes with PRH. Three 20-ns molecular dynamics simulations (PRH-TAATTG, PRH-TAATTA and PRH-TAATGG complexes) in explicit solvent water were performed to investigate three complexes structurally. Structural analysis shows that the native TAATTG sequence forms a complex that is more stable than complexes with base pair mutations. It is also observed that upon mutation, the number and occupancy of the direct and water-mediated hydrogen bonds decrease. Free energy calculations performed with the thermodynamic integration method predict relative binding free energies of 0.64 and 2 kcal/mol for GC to AT and TA to GC mutations, respectively, suggesting that among the three DNA sequences, the PRH-TAATTG complex is more stable than the two mutated complexes. In addition, it is demonstrated that the stability of the PRH-TAATTA complex is greater than that of the PRH-TAATGG complex.

Characterization of experimentally determined native-structure models of a protein using energetic and entropic components of free-energy function

We show how to characterize the native-structure models of a protein using our free-energy function F which is based on hydration thermodynamics. Ubiquitin is adopted as an example protein. We consider models determined by the X-ray crystallography and two types of NMR model sets. A model set of type 1 comprises candidate models for a fixed native structure, and that of type 2 forms an ensemble of structures representing the structural variability of the native state. In general, the X-ray models give lower F than the NMR models. There is a trend that, as a model deviates more from the model with the lowest F among the X-ray models, its F becomes higher. Model sets of type 1 and those of type 2, respectively, exhibit two different characteristics with respect to the correlation between the deviation and F. It is argued that the total amount of constraints such as NOEs effectively taken into account in constructing the NMR models can be examined by analyzing the behavior of F. We investigate structural characteristics of the models in terms of the energetic and entropic components of F which are relevant to intramolecular hydrogen bonding and to backbone and side-chain packing, respectively.

A scoring function based on solvation thermodynamics for protein structure prediction

We predict protein structure using our recently developed free energy function for describing protein stability, which is focused on solvation thermodynamics. The function is combined with the current most reliable sampling methods, i.e., fragment assembly (FA) and comparative modeling (CM). The prediction is tested using 11 small proteins for which high-resolution crystal structures are available. For 8 of these proteins, sequence similarities are found in the database, and the prediction is performed with CM. Fairly accurate models with average Cα root mean square deviation (RMSD) ∼ 2.0 Å are successfully obtained for all cases. For the rest of the target proteins, we perform the prediction following FA protocols. For 2 cases, we obtain predicted models with an RMSD ∼ 3.0 Å as the best-scored structures. For the other case, the RMSD remains larger than 7 Å. For all the 11 target proteins, our scoring function identifies the experimentally determined native structure as the best structure. Starting from the predicted structure, replica exchange molecular dynamics is performed to further refine the structures. However, we are unable to improve its RMSD toward the experimental structure. The exhaustive sampling by coarse-grained normal mode analysis around the native structures reveals that our function has a linear correlation with RMSDs < 3.0 Å. These results suggest that the function is quite reliable for the protein structure prediction while the sampling method remains one of the major limiting factors in it. The aspects through which the methodology could further be improved are discussed.

Free-energy-based methods for binding profile determination in a congeneric series of CDK2 inhibitors

Free-energy pathway methods show great promise in computing the mode of action and the free energy profile associated with the binding of small molecules with proteins, but are generally very computationally demanding. Here we apply a novel approach based on metadynamics and path collective variables. We show that this combination is able to find an optimal reaction coordinate and the free energy profile of binding with explicit solvent and full flexibility, while minimizing human intervention and computational costs. We apply it to predict the binding affinity of a congeneric series of 5 CDK2 inhibitors. The predicted binding free energy profiles are in accordance with experiment.

CHARMM: the biomolecular simulation program

CHARMM (Chemistry at HARvard Molecular Mechanics) is a highly versatile and widely used molecular simulation program. It has been developed over the last three decades with a primary focus on molecules of biological interest, including proteins, peptides, lipids, nucleic acids, carbohydrates, and small molecule ligands, as they occur in solution, crystals, and membrane environments. For the study of such systems, the program provides a large suite of computational tools that include numerous conformational and path sampling methods, free energy estimators, molecular minimization, dynamics, and analysis techniques, and model-building capabilities. The CHARMM program is applicable to problems involving a much broader class of many-particle systems. Calculations with CHARMM can be performed using a number of different energy functions and models, from mixed quantum mechanical-molecular mechanical force fields, to all-atom classical potential energy functions with explicit solvent and various boundary conditions, to implicit solvent and membrane models. The program has been ported to numerous platforms in both serial and parallel architectures. This article provides an overview of the program as it exists today with an emphasis on developments since the publication of the original CHARMM article in 1983.

Packaging effects on site-specific DNA-protein interactions

We show that the rate of site-specific association of a protein molecule of interest with the DNA chain can be approximately 10(2) times higher than that of the three-dimensional diffusion-controlled collision rate limit approximately 10(8) mol(-1) s(-1) only when the protein molecule of interest searches for its specific site on the DNA chain in a reduced dimensional space with a dimensionality dr of dr<1. Upon considering the concurrent dynamics of the linear DNA chain that is embedded in a d-dimensional space along with the one-dimensional diffusion dynamics of the nonspecifically bound protein molecule on the DNA chain, we derive the generalized scaling law epsilon approximately 2(3(2-d)+3), where epsilon is the number of times by which the rate of site-specific association of the protein molecule with the DNA chain can be enhanced over the three-dimensional diffusion-controlled collision rate limit and d is the dimensionality of the reduced search space. Using the analogy between the self-intersection loop length in the theory of random walks and the ring-closure events in the theory of site specific interactions of a protein molecule with the DNA chain, we further show that the extent of packaging and volume compression of the genomic DNA inside the living cell is designed in such a way that the efficiency of the protein molecule in the process of searching for its specific site on the genomic DNA is a maximum. Our simulation results suggest that the volume compression factor theta which is the ratio between the total volume of the living cell and the volume occupied by the DNA chain along with all the other bound protein molecules should be such that theta>or=100 for an efficient site specific interaction of a protein molecule of interest with the linear DNA chain that is embedded in a three-dimensional space. Our theoretical and simulation results agree well with the E. coli cellular system.

Importance of translational entropy of water in biological self-assembly processes like protein folding

We briefly review our studies on the folding/unfolding mechanisms of proteins. In biological self-assembly processes such as protein folding, the number of accessible translational configurations of water in the system increases greatly, leading to a large gain in the water entropy. The usual view looking at only the water in the close vicinity of the protein surface is capable of elucidating neither the large entropic gain upon apoplastocyanin folding, which has recently been found in a novel experimental study, nor the pressure and cold denaturation. With the emphasis on the translational entropy of water, we are presently constructing a reliable method for predicting the native structure of a protein from its amino-acid sequence.

Computational evaluation of protein-small molecule binding

Determining protein-small molecule binding affinity is a key component of present-day rational drug discovery. To circumvent the time, labor, and materials costs associated with experimental protein-small molecule binding assays, a variety of structure-based computational methods have been developed for determining protein-small molecule binding affinities. These methods can be placed in one of two classes: accurate but slow (Class 1), and fast but approximate (Class 2). Class 1 methods, which explicitly take into account protein flexibility and include an atomic-level description of solvation, are capable of quantitatively reproducing experimental protein-small molecule absolute binding free energies. However, Class 1 computational requirements make screening thousands to millions of small molecules against a protein, as required for rational drug design, infeasible for the foreseeable future. Class 2 methods, on the contrary, are sufficiently fast to perform such inhibitor screening, yet they suffer from limited descriptions of protein flexibility and solvation, which in turn limit their ability to select and rank-order small molecules by computed binding affinities. This review presents an overview of Class 1 and Class 2 methods, and avenues of research in Class 2 methods aimed at bringing them closer to Class 1 accuracy.

Generalized theory of site-specific DNA-protein interactions

We develop a generalized theory of the site-specific DNA-protein interactions, which includes both the static as well as the dynamical factors influencing the one-dimensional diffusion of the nonspecifically bound protein molecule which is in the process of searching for the specific site on the DNA lattice. We argue that the chemically driven condensation of the DNA molecule introduces a static distribution in the one-dimensional phenomenological diffusion coefficient associated with the protein molecule and the conformational dynamics of the DNA introduces temporal fluctuations in the one-dimensional diffusion coefficient over the static distribution. We further derive the generalized inequality conditions and the scaling laws which are required to enhance the three-dimensional diffusion controlled site-specific association rate to an arbitrary order. Our model predicts that when the degree of condensation of the DNA molecule under consideration is very high, then the probability distribution associated with the stationary state one-dimensional diffusion coefficient variable as well as the stationary state one-dimensional diffusion length variable will be a flat one. Further analysis reveals that to achieve a site-specific association rate which is higher than that of the three-dimensional diffusion controlled rate, the one-dimensional diffusion length associated with the dynamics of the nonspecifically bound protein molecule on the DNA lattice should fall in certain critical ranges. Comparison of our theoretical results with the recent experimental observations reveals that when the DNA molecule is under a stretched condition, then the static distribution of the one-dimensional diffusion coefficient associated with the dynamics of the protein molecule on the DNA lattice is a Gaussian and therefore the fluctuations in the one-dimensional diffusion coefficient generated by the dynamical factors are confined in a harmonic type potential.

Flexible docking in solution using metadynamics

We apply our recently developed metadynamics method to the docking of ligands on flexible receptors in water solution. This method mimics the real dynamics of a ligand exiting or entering an enzyme and in so doing reconstructs the free energy surface. We apply it to four docking cases: beta-trypsin/benzamidine, beta-trypsin/chlorobenzamidine, immunoglobulin McPC-603/phosphocholine, and cyclin-dependent kinase 2/staurosporine. In every case studied, the method is able to predict the docked geometry and the free energy of docking. Its added value with respect to many other available methods is that it reconstructs the complete free energy surface, including all the relevant minima and the barriers between them.

Peptide and protein folding and conformational equilibria: theoretical treatment of electrostatics and hydrogen bonding with implicit solvent models

Since biomolecules exist in aqueous and membrane environments, the accurate modeling of solvation, and hydrogen bonding interactions in particular, is essential for the exploration of structure and function in theoretical and computational studies. In this chapter, we focus on alternatives to explicit solvent models and discuss recent advances in generalized Born (GB) implicit solvent theories. We present a brief review of the successes and shortcomings of the application of these theories to biomolecular problems that are strongly linked to backbone H-bonding and electrostatics. This discussion naturally leads us to explore existing areas for improvement in current GB theories and our approach towards addressing a number of the key issues that remain in the refinement of these models. Specifically, the critical importance of balancing solvation forces and intramolecular forces in GB models is illustrated by examining the influence of backbone hydrogen bond strength and backbone dihedral energetics on conformational equilibria of small peptids.

Standard free energy of releasing a localized water molecule from the binding pockets of proteins: double-decoupling method

Localized water molecules in the binding pockets of proteins play an important role in noncovalent association of proteins and small drug compounds. At times, the dominant contribution to the binding free energy comes from the release of localized water molecules in the binding pockets of biomolecules. Therefore, to quantify the energetic importance of these water molecules for drug design purposes, we have used the double-decoupling approach to calculate the standard free energy of tying up a water molecule in the binding pockets of two protein complexes. The double-decoupling approach is based on the underlying principle of statistical thermodynamics. We have calculated the standard free energies of tying up the water molecule in the binding pockets of these complexes to be favorable. These water molecules stabilize the protein-drug complexes by interacting with the ligands and binding pockets. Our results offer ideas that could be used in optimizing protein-drug interactions, by designing ligands that are capable of targeting localized water molecules in protein binding sites. The resulting free energy of ligand binding could benefit from the potential free energy gain accompanying the release of these water molecules. Furthermore, we have examined the theoretical background of the double-decoupling method and its connection to the molecular dynamics thermodynamic integration techniques.

A Sequential Molecular Mechanics/Quantum Mechanics Study of the Electronic Spectra of Amides

We report gas-phase electronic spectra of formamide, N-methyformamide, acetamide, and N-methylacetamide at 300 K calculated using a combination of classical molecular dynamics and time-dependent density functional theory (TDDFT). In comparison to excitation energies computed using the global minima structures, the valence npi* and pi(nb)pi* states show a significant red-shift of 0.1-0.35 eV, while smaller shifts are found for the n3s and pi(nb)3s Rydberg states. In this work, we have identified the physical origin of these shifts arising from variations of the molecular structure. We present simple relationships between key geometrical parameters and spectral shifts. Consequently, electronic spectra can be generated directly from ground-state structures, without additional quantum chemical calculations. The electronic spectrum of formamide in aqueous solution is computed using TDDFT using an explicit solvent model. This provides a quantitative determination of the condensed-phase spectrum. In general, this study shows that temperature effects can change the predicted excitation energies significantly and demonstrates how electronic spectra at elevated temperatures can be computed in a computationally efficient way.

Protein and peptide folding explored with molecular simulations

Molecular simulations, comprising models with atomic details of polypeptide and solvent as well as minimalist models employing only C alpha atoms, are being used with specialized simulation methods from statistical mechanics to examine fundamental questions in peptide and protein folding mechanism, kinetics, and thermodynamics. Detailed calculations of free energy changes along coordinates describing the formation of hydrogen-bonding interactions in helical, turn, and beta-sheet models provide insights into the time scale and mechanism of secondary structure formation. Potential roles for these processes in directing protein folding are also elucidated by such calculations. Analogous methodologies extended to more complex polypeptides with tertiary structures (proteins) are used to explore global questions about protein folding landscapes, to delineate atomic details of folding mechanism, and to elucidate putative roles for solvent in the late stages of folding.

From folding theories to folding proteins: a review and assessment of simulation studies of protein folding and unfolding

Beginning with simplified lattice and continuum "minimalist" models and progressing to detailed atomic models, simulation studies have augmented and directed development of the modern landscape perspective of protein folding. In this review we discuss aspects of detailed atomic simulation methods applied to studies of protein folding free energy surfaces, using biased-sampling free energy methods and temperature-induced protein unfolding. We review studies from each on systems of particular experimental interest and assess the strengths and weaknesses of each approach in the context of "exact" results for both free energies and kinetics of a minimalist model for a beta-barrel protein. We illustrate in detail how each approach is implemented and discuss analysis methods that have been developed as components of these studies. We describe key insights into the relationship between protein topology and the folding mechanism emerging from folding free energy surface calculations. We further describe the determination of detailed "pathways" and models of folding transition states that have resulted from unfolding studies. Our assessment of the two methods suggests that both can provide, often complementary, details of folding mechanism and thermodynamics, but this success relies on (a) adequate sampling of diverse conformational regions for the biased-sampling free energy approach and (b) many trajectories at multiple temperatures for unfolding studies. Furthermore, we find that temperature-induced unfolding provides representatives of folding trajectories only when the topology and sequence (energy) provide a relatively funneled landscape and "off-pathway" intermediates do not exist.

Thermodynamic basis for sequence-specific recognition of ssDNA by an autoantibody

11F8 is a sequence-specific DNA binding monoclonal autoantibody previously isolated from an autoimmune lupus-prone mouse [Stevens, S. Y., and Glick, G. D. (1999) Biochemistry 38, 560-568]. This antibody, like many other lupus anti-DNAs, localizes to kidney tissue and eventually leads to renal damage through a process that first involves the binding of DNA antigens. A series of experiments were conducted to investigate the thermodynamic and structural basis by which this antibody discriminates between specific, noncognate, and nonspecific sequences. Sequence-specific binding occurs with a minimal dependence on the polyelectrolyte effect along with a favorable binding enthalpy reflecting the presence of base stacking and contacts to DNA bases. This favorable binding enthalpy apparently is derived from desolvation at the binding interface and is consistent with recent models of the nonclassical hydrophobic effect. Noncognate recognition is also driven by the nonclassical hydrophobic effect, but is accompanied by highly unfavorable entropies that are responsible for reduced affinity relative to the high-affinity consensus sequence. Nonspecific recognition is driven completely by the polyelectrolyte effect involving extensive electrostatic interactions with the phosphate backbone. Collectively, the data demonstrate the ability of 11F8 to adapt its mode of binding to the available DNA surface and provide a thermodynamic model for sequence-specific recognition of single-stranded DNA. The salient features of this model employ the paradigms invoked to explain protein.dsDNA, protein.RNA, and antibody.antigen binding.

Energetics of the interaction between water and the helical peptide group and its role in determining helix propensities

The alanine helix provides a model system for studying the energetics of interaction between water and the helical peptide group, a possible major factor in the energetics of protein folding. Helix formation is enthalpy-driven (-1.0 kcal/mol per residue). Experimental transfer data (vapor phase to aqueous) for amides give the enthalpy of interaction with water of the amide group as approximately -11.5 kcal/mol. The enthalpy of the helical peptide hydrogen bond, computed for the gas phase by quantum mechanics, is -4.9 kcal/mol. These numbers give an enthalpy deficit for helix formation of -7.6 kcal/mol. To study this problem, we calculate the electrostatic solvation free energy (ESF) of the peptide groups in the helical and beta-strand conformations, by using the delphi program and parse parameter set. Experimental data show that the ESF values of amides are almost entirely enthalpic. Two key results are: in the beta-strand conformation, the ESF value of an interior alanine peptide group is -7.9 kcal/mol, substantially less than that of N-methylacetamide (-12.2 kcal/mol), and the helical peptide group is solvated with an ESF of -2.5 kcal/mol. These results reduce the enthalpy deficit to -1.5 kcal/mol, and desolvation of peptide groups through partial burial in the random coil may account for the remainder. Mutant peptides in the helical conformation show ESF differences among nonpolar amino acids that are comparable to observed helix propensity differences, but the ESF differences in the random coil conformation still must be subtracted.

Nucleic acid base-pairing and N-methylacetamide self-association in chloroform: affinity and conformation

A recently developed computational method, 'mining minima', is used to examine the hydrogen-bonding interactions of nucleic acid base-pairs and of the N-methylacetamide homodimer in chloroform. The mining minima algorithm aggressively samples molecular conformations, identifies the most important local minima, and computes their contributions to the overall free energy of the system. Here, the CHARMM 98 parameter set is used for the potential energy and the generalized Born/surface area solvent model is used to account for the influence of the solvent. Good agreement with experiment is obtained for the non-covalent binding affinities of a series of complexes. The computational approach used here is applicable to a range of molecular systems.

Thermodynamics of a reverse turn motif. Solvent effects and side-chain packing

The linear pentapeptide, Ala-Tyr-cis-Pro-Tyr-Asp-NMA (AYPYD) is known to have a significant population of type VI turn conformers in aqueous solvent. We have carried out theoretical studies of the conformational energetics of this peptide using a potential of mean force (PMF) consisting of the AMBER/OPLS empirical potential energy function, a macroscopic electrostatic model of polar solvation, and a surface area-based model of non-polar solvation. Conformers were taken from molecular dynamics simulations reported elsewhere, or generated by a random search method reported here. The chain entropy of folding was calculated by a systematic search of accessible dihedral angle space. The intra-peptide component was found to strongly favor folding and was nearly cancelled by the polar solvation term which disfavored folding. The non-polar solvation term had little effect. Fluctuations about the average value of the PMF were small and in accord with estimates from a simple harmonic model. When applied to conformers generated by a random search, the PMF selected a conformer close to the NMR-determined structure as the lowest energy conformer. The conformer with the second-lowest energy was extended, but was found to fold rapidly to the turn state in a subsequent molecular dynamics study, and may be an important state on the folding-unfolding pathway. Averages of the PMF were combined with the entropy estimates to provide an estimate of the free energy of folding that is in reasonable agreement with experimental results. In terms of the interplay between backbone electrostatic interactions and the packing of apolar side-chains, this peptide provides a model for the energetics of protein folding, and therefore makes a useful test case for calculations.