Evaluation of protein-protein association energies by free energy perturbation calculations

The association energy upon binding of different amino acids in the specificity pocket of trypsin was evaluated by free energy perturbation calculations on complexes between bovine trypsin (BT) and bovine pancreatic trypsin inhibitor (BPTI). Three simulations of mutations of the primary binding residue (P(1)) were performed (P(1)-Ala to Gly, P(1)-Met to Gly and P(1)-Met to Ala) and the resulting differences in association energy (DeltaDeltaG(a)) are 2. 28, 5.08 and 2.93 kcal/mol for P(1)-Ala to Gly, P(1)-Met to Gly and to Ala with experimental values of 1.71, 4.62 and 2.91 kcal/mol, respectively. The calculated binding free energy differences are hence in excellent agreement with the experimental binding free energies. The binding free energies, however, were shown to be highly dependent on water molecules at the protein-protein interface and could only be quantitatively estimated if the correct number of such water molecules was included. Furthermore, the cavities that were formed when a large amino acid side-chain is perturbed to a smaller one seem to create instabilities in the systems and had to be refilled with water molecules in order to obtain reliable results. In addition, if the protein atoms that were perturbed away were not replaced by water molecules, the simulations dramatically overestimated the initial state of the free energy perturbations.     

Crystal structural analysis of protein-protein interactions drastically destabilized by a single mutation

The complex of barnase (bn) and barstar (bs), which has been widely studied as a model for quantitative analysis of protein-protein interactions, is significantly destabilized by a single mutation, namely, bs Asp39 --> Ala, which corresponds to a change of 7.7 kcal x mol(-1) in the free energy of binding. However, there has been no structural information available to explain such a drastic destabilization. In the present study, we determined the structure of the mutant complex at 1.58 A resolution by X-ray crystallography. The complex was similar to the wild-type complex in terms of overall and interface structures; however, the hydrogen bond network mediated by water molecules at the interface was significantly different. Several water molecules filled the cavity created by the mutation and consequently caused rearrangement of the hydrated water molecules at the interface. The water molecules were redistributed into a channel-like structure that penetrated into the complex. Furthermore, molecular dynamics simulations showed that the mutation increased the mobility of water molecules at the interface. Since such a drastic change in hydration was not observed in other mutant complexes of bn and bs, the significant destabilization of the interaction may be due to this channel-like structure of hydrated water molecules.     

Quantification of helix-helix binding affinities in micelles and lipid bilayers

A theoretical approach for estimating association free energies of alpha-helices in nonpolar media has been developed. The parameters of energy functions have been derived from DeltaDeltaG values of mutants in water-soluble proteins and partitioning of organic solutes between water and nonpolar solvents. The proposed approach was verified successfully against three sets of published data: (1) dissociation constants of alpha-helical oligomers formed by 27 hydrophobic peptides; (2) stabilities of 22 bacteriorhodopsin mutants, and (3) protein-ligand binding affinities in aqueous solution. It has been found that coalescence of helices is driven exclusively by van der Waals interactions and H-bonds, whereas the principal destabilizing contributions are represented by side-chain conformational entropy and transfer energy of atoms from a detergent or lipid to the protein interior. Electrostatic interactions of alpha-helices were relatively weak but important for reproducing the experimental data. Immobilization free energy, which originates from restricting rotational and translational rigid-body movements of molecules during their association, was found to be less than 1 kcal/mole. The energetics of amino acid substitutions in bacteriorhodopsin was complicated by specific binding of lipid and water molecules to cavities created in certain mutants.     

Thermodynamic stability of water molecules in the bacteriorhodopsin proton channel: a molecular dynamics free energy perturbation study.

The proton transfer activity of the light-driven proton pump, bacteriorhodopsin (bR) in the photochemical cycle might imply internal water molecules. The free energy of inserting water molecules in specific sites along the bR transmembrane channel has been calculated using molecular dynamics simulations based on a microscopic model. The existence of internal hydration is related to the free energy change on transfer of a water molecule from bulk solvent into a specific binding site. Thermodynamic integration and perturbation methods were used to calculate free energies of hydration for each hydrated model from molecular dynamics simulations of the creation of water molecules into specific protein-binding sites. A rigorous statistical mechanical formulation allowing the calculation of the free energy of transfer of water molecules from the bulk to a protein cavity is used to estimate the probabilities of occupancy in the putative bR proton channel. The channel contains a region lined primarily by nonpolar side-chains. Nevertheless, the results indicate that the transfer of four water molecules from bulk water to this apparently hydrophobic region is thermodynamically permitted. The column forms a continuous hydrogen-bonded chain over 12 A between a proton donor, Asp 96, and the retinal Schiff base acceptor. The presence of two water molecules in direct hydrogen-bonding association with the Schiff base is found to be strongly favorable thermodynamically. The implications of these results for the mechanism of proton transfer in bR are discussed.     

Study of the insulin dimerization: binding free energy calculations and per-residue free energy decomposition

A calculation of the binding free energy for the dimerization of insulin has been performed using the molecular mechanics-generalized Born surface area approach. The calculated absolute binding free energy is -11.9 kcal/mol, in approximate agreement with the experimental value of -7.2 kcal/mol. The results show that the dimerization is mainly due to nonpolar interactions. The role of the hydrogen bonds between the 2 monomers appears to give the direction of the interactions. A per-atom decomposition of the binding free energy has been performed to identify the residues contributing most to the self association free energy. Residues B24-B26 are found to make the largest favorable contributions to the dimerization. Other residues situated at the interface between the 2 monomers were found to make favorable but smaller contributions to the dimerization: Tyr B16, Val B12, and Pro B28, and to an even lesser extent, Gly B23. The energy decomposition on a per-residue basis is in agreement with experimental alanine scanning data. The results obtained from a single trajectory (i.e., the dimer trajectory is also used for the monomer analysis) and 2 trajectories (i.e., separate trajectories are used for the monomer and dimer) are similar.     

Computational Analysis of Phosphopeptide Binding to the Polo-Box Domain of the Mitotic Kinase PLK1 Using Molecular Dynamics Simulation

The Polo-Like Kinase 1 (PLK1) acts as a central regulator of mitosis and is over-expressed in a wide range of human tumours where high levels of expression correlate with a poor prognosis. PLK1 comprises two structural elements, a kinase domain and a polo-box domain (PBD). The PBD binds phosphorylated substrates to control substrate phosphorylation by the kinase domain. Although the PBD preferentially binds to phosphopeptides, it has a relatively broad sequence specificity in comparison with other phosphopeptide binding domains. We analysed the molecular determinants of recognition by performing molecular dynamics simulations of the PBD with one of its natural substrates, CDC25c. Predicted binding free energies were calculated using a molecular mechanics, Poisson-Boltzmann surface area approach. We calculated the per-residue contributions to the binding free energy change, showing that the phosphothreonine residue and the mainchain account for the vast majority of the interaction energy. This explains the very broad sequence specificity with respect to other sidechain residues. Finally, we considered the key role of bridging water molecules at the binding interface. We employed inhomogeneous fluid solvation theory to consider the free energy of water molecules on the protein surface with respect to bulk water molecules. Such an analysis highlights binding hotspots created by elimination of water molecules from hydrophobic surfaces. It also predicts that a number of water molecules are stabilized by the presence of the charged phosphate group, and that this will have a significant effect on the binding affinity. Our findings suggest a molecular rationale for the promiscuous binding of the PBD and highlight a role for bridging water molecules at the interface. We expect that this method of analysis will be very useful for probing other protein surfaces to identify binding hotspots for natural binding partners and small molecule inhibitors.     

Stability of macromolecular complexes

Macromolecular interactions are crucial in numerous biologic processes, yet few general principles are available that establish firm expectations for the strength of these interactions or the expected contribution of specific forces. The simplest principle would be a monotonic increase in interactions as the size of the interface grows. The exact relationship might be linear or nonlinear depending on the nature of the forces involved. Simple "linear-free energy" relationships based on atomic properties have been well documented, for example, additivity for the interaction of small molecules with solvent, and, recently, have been explored for ligand-receptor interactions. Horton and Lewis propose such additivity based on buried surface area for protein-protein complexes. We investigated macromolecular interactions and found that the highest-affinity complexes do not fulfill this simple expectation. Instead, binding free energies of the tightest macromolecular complexes are roughly constant, independent of interface size, with the notable exception of DNA duplexes. By comparing these results to an earlier study of protein-ligand interactions we find that: (1) The maximum affinity is approximately 1.5 kcal/mol per nonhydrogen atom or 120 cal/mol A(2) of buried surface area, comparable to results of our earlier work; (2) the lack of an increase in affinity with interface size is likely due to nonthermodynamic factors, such as functional and evolutionary constraints rather than some fundamental physical limitation. The implication of these results have some importance for molecular design because they suggest that: (1) The stability of any given complex can be increased significantly if desired; (2) small molecule inhibitors of macromolecular interactions are feasible; and (3) different functional classes of protein-protein complexes exhibit differences in maximal stability, perhaps in response to differing evolutionary pressures. These results are consistent with the widespread observation that proteins have not evolved to maximize thermodynamic stability, but are only marginally stable.     

Continuum electrostatic analysis of preferred solvation sites around proteins in solution

To understand water-protein interactions in solution, the electrostatic field is calculated by solving the Poisson-Boltzmann equation, and the free energy surface of water is mapped by translating and rotating an explicit water molecule around the protein. The calculation is applied to T4 lysozyme with data available on the conservation of solvent binding sites in 18 crystallographically independent molecules. The free energy maps around the ordered water sites provide information on the relationship between water positions in crystal structure and in solution. Results show that almost all conserved sites and the majority of nonconserved sites are within 1.3 A of local free energy minima. This finding is in sharp contrast to the behavior of randomly placed water molecules in the boundary layer, which, on the average, must travel more than 3 A to the nearest free energy minimum. Thus, the solvation sites are at least partially determined by protein-water interactions rather than by crystal packing alone. The characteristic water residence times, obtained from the free energies at the local minima, are in good agreement with nuclear magnetic resonance experiments. Only about half of the potential sites show up as ordered water in the 1.7 A resolution X-ray structure. Crystal packing interactions can stabilize weak or mobile potential sites (in fact, some ordered water positions are not close to free energy minima) or can prevent water from occupying certain sites. Apart from a few buried water molecules that are strong binders, the free energies are not very different for conserved and nonconserved sites. We show that conservation of a water site between two crystals occurs if the positions of protein atoms, primarily contributing to the free energy at the local minimum, do not substantially change from one structure to the other. This requirement can be correlated with the nature of the side chain contacting the water molecule in the site.     

Effect of the external electrostatic field on the conformation of a macromolecule containing charged groups

A method for calculating the free energy of a macromolecule containing charged groups in electrostatic field in aqueous solution was proposed. The non-electrostatic component of free energy was calculated with consideration of van der Waals interactions between uncharged parts of the macromolecule. The electrostatic component of free energy was calculated with regard for the interactions of charged groups of the macromolecule with each other and with water molecules. It was found that, depending on the strength of external electric field, the free energy of the system passes through a minimum, whereas the internal energy passes through a maximum. By minimizing the free energy, relative changes in the mean radius 'r' and the distance between the termini of the macromolecule 'h' were calculated. It was found that, at some values of field strength, both 'r' and 'h' decrease. An increase in strength led to an increase in 'r' and 'h'. The regularities observed depend on the charge of the macromolecule and the spatial redistribution of macromolecules and counterions.     

Binding of the Bacteriophage P22 N-Peptide to the boxB RNA Motif Studied by Molecular Dynamics Simulations

Protein-RNA interactions are important for many cellular processes. The Nut-utilization site (N)-protein of bacteriophages contains an N-terminal arginine-rich motif that undergoes a folding transition upon binding to the boxB RNA hairpin loop target structure. Molecular dynamics simulations were used to investigate the dynamics of the P22 N-peptide-boxB complex and to elucidate the energetic contributions to binding. In addition, the free-energy changes of RNA and peptide conformational adaptation to the bound forms, as well as the role of strongly bound water molecules at the peptide-RNA interface, were studied. The influence of peptide amino acid substitutions and the salt dependence of interaction were investigated and showed good agreement with available experimental results. Several tightly bound water molecules were found at the RNA-binding interface in both the presence and absence of N-peptide. Explicit consideration of the waters resulted in shifts of calculated contributions during the energetic analysis, but overall similar binding energy contributions were found. Of interest, it was found that the electrostatic field of the RNA has a favorable influence on the coil-to-α-helix transition of the N-peptide already outside of the peptide-binding site. This result may have important implications for understanding peptide-RNA complex formation, which often involves coupled folding and association processes. It indicates that electrostatic interactions near RNA molecules can lead to a shift in the equilibrium toward the bound form of an interacting partner before it enters the binding pocket.     

Intermolecular interactions between protein and other molecules including hydration effects

The structural aspects of protein functions, e.g., molecular recognition such as enzyme-substrate and antibody-antigen interactions, are elucidated in terms of dehydration and atomic interactions. When a protein interacts with some target molecule, water molecules at the interacting regions of both molecules are removed, with loss of the hydration free energy, but gaining atomic interactions between atoms of the contact sites in both molecules. The free energies of association originating from the dehydration and interactions between the atoms can be computed from changes in the accessible surface areas of the atoms involved. The free energy due to interactions between atomic groups at the contact sites is estimated as the sum of those estimated from the changes in the accessible surface area of 7 atomic groups, assuming that the interactions are proportional to the change of the area. The chain enthalpies and entropies evaluated from experimental thermodynamic properties and hydration quantities at the standard temperature for 10 proteins were available to determine the proportional constants for the atomic groups. This method was applied to the evaluation of association constants for the dimerization of proteins and the formation of proteolytic enzyme-inhibitor complexes, and the computed constants were in agreement with the experimental ones. However, the method is not accurate enough to account quantitatively for the change in the thermal stability of mutants of T4 lysozyme. Nevertheless, this method provides a way to elucidate the interactions between molecules in solution.     

Computational methodology for estimating changes in free energies of biomolecular association upon mutation. The importance of bound water in dimer-tetramer assembly for beta 37 mutant hemoglobins

The computational modeling program HINT (Hydropathic INTeractions), an empirical hydropathic force field that includes hydrogen bonding, Coulombic, and hydrophobic terms, was used to model the free energy of dimer-tetramer association in a series of deoxy hemoglobin beta 37 double mutants. Five of the analyzed mutants (beta 37W --> Y, beta 37W --> A, beta 37W --> G, beta 37W --> E, and beta 37W --> R) have been solved crystallographically and characterized thermodynamically and subsequently made a good test set for the calibration of our method as a tool for free energy prediction. Initial free energy estimates for these mutants were conducted without the inclusion of crystallographically conserved water molecules and systematically underestimated the experimentally calculated loss in free energy observed for each mutant dimer-tetramer association. However, the inclusion of crystallographic waters, interacting at the dimer-dimer interface of each mutant, resulted in HINT free energy estimates that were more accurate with respect to experimental data. To evaluate the ability of our method to predict free energies for de novo protein models, the same beta 37 mutants were computationally generated from native deoxy hemoglobin and similarly analyzed. Our theoretical models were sufficiently robust to accurately predict free energy changes in a localized region around the mutated residue. However, our method did not possess the capacity to generate the long-range secondary structural effects observed in crystallographically solved mutant structures. Final method analysis involved the computational generation of structurally and/or thermodynamically uncharacterized beta 37 deoxy hemoglobin mutants. HINT analysis of these structures revealed that free energy predictions for dimer-tetramer association in these models agreed well with previously observed energy predictions for structurally and thermodynamically characterized beta 37 deoxy hemoglobin mutants.     

Salting-in of neopentane in the aqueous solutions of urea and glycine-betaine

The effects of urea and glycine-betaine (GB) osmolytes on the hydrophobic interactions of neopentane in water have been studied using molecular dynamics simulations. From the study of the potentials of mean force, it is observed that both urea and GB decrease the association and solvation of neopentane. The calculated equilibrium constants show that urea and GB decrease the population of solvent-separated minima of neopentane. The hydrophobic association as well as solvation of neopentane molecules are stabilised by entropy and enthalpy in the mixtures. The radial distribution functions (RDFs) and running coordination numbers of water, urea and GB molecules show that neopentane shows salting-in behaviour in aqueous-GB, aqueous-urea and aqueous-urea-GB mixtures. Neopentane is preferentially solvated by GB in aqueous-GB and preferentially solvated by urea in aqueous-urea-GB solutions. The preferential solvation of neopentane by GB suggests that GB decreases the interaction between neopentane molecules i.e. salting-in of neopentane. The calculated solvation free energies and radial density profiles of neopentane also support the salting-in behaviour of neopentane in the mixtures of these osmolytes.     

Molecular dynamics simulation of the human U2B" protein complex with U2 snRNA hairpin IV in aqueous solution.

A 2200-ps molecular dynamics (MD) simulation of the U2 snRNA hairpin IV/U2B" complex was performed in aqueous solution using the particle mesh Ewald method to consider long-range electrostatic interactions. To investigate the interaction and recognition process between the RNA and protein, the free energy contributions resulting from individual amino acids of the protein component of the RNA/protein complex were calculated using the recently developed glycine-scanning method. The results revealed that the loop region of the U2 snRNA hairpin IV interacted mainly with three regions of the U2B" protein: 1) beta 1-helix A, 2) beta 2-beta 3, and 3) beta 4-helix C. U2 snRNA hairpin IV bound U2B" in a similar orientation as that previously described for U1 snRNA with the U1A' protein; however, the details of the interaction differed in several aspects. In particular, beta 1-helix A and beta 4-helix C in U2B" were not observed to interact with RNA in the U1A' protein complex. Most of the polar and charged residues in the interacting regions had larger mutant free energies than the nonpolar residues, indicating that electrostatic interactions were important for stabilizing the RNA/protein complex. The interaction was further stabilized by a network of hydrogen bonds and salt bridges formed between RNA and protein that was maintained throughout the MD trajectory. In addition to the direct interactions between RNA and the protein, solvent-mediated interactions also contributed significantly to complex stability. A detailed analysis of the ordered water molecules in the hydration of the RNA/protein complex revealed that bridged water molecules reside at the interface of RNA and protein as long as 2100 ps in the 2200-ps trajectory. At least 20 bridged water molecules, on average, contributed to the instantaneous stability of the RNA/protein complex. The stabilizing interaction energy due to bridging water molecules was obtained from ab initio Hartree-Fock and density functional theory calculations.     

Molecular simulation study of cooperativity in hydrophobic association

To investigate the cooperativity of hydrophobic interactions, the potential of mean force of two- and three-molecule methane clusters in water was determined by molecular dynamics simulations using two methods: umbrella-sampling with the weighted histogram analysis method and thermodynamic integration. Two water models, TIP3P and TIP4P, were used, while each methane molecule was modeled as a united atom. It was found that the three-body potential of mean force is not additive, i.e., it cannot be calculated as a sum of two-body contributions, but requires an additional three-body cooperative term. The cooperative term, which amounts to only about 10% of the total hydrophobic association free energy, was found to increase the strength of hydrophobic association; this finding differs from the results of earlier Monte Carlo studies with the free energy perturbation method of Rank and Baker (1997). As in the work of Rank and Baker, the solvent contribution to the potential of mean force was found to be well approximated by the molecular surface of two methane molecules. Moreover, we also found that the cooperative term is well represented by the difference between the molecular surface of the three-methane cluster and those of all three pairs of methane molecules. In addition, it was found that, while there is a cooperative contribution to the hydrophobic association free energy albeit a small one, the errors associated with the use of pairwise potentials are comparable to or larger than this contribution.     

Lipoplex thermodynamics: determination of DNA-cationic lipoid interaction energies

An experimental study of the cationic lipid-DNA binding affinity is presented. The binding free energy was determined by monitoring lipoplex dissociation under conditions of increasing salt concentration. The primary procedure was based on the extent of quenching by energy transfer of fluorophores on DNA molecules by fluorophore on a lipid as these molecules came into close association in the lipoplex. Titration calorimetry on the Dickerson dodecamer was also done, with results that were in agreement with the fluorescence data. Measurements on short oligonucleotides allowed estimation of the binding energy per nucleotide. The binding free energy is approximately 0.6 kcal/mole nucleotide for the Dickerson dodecamer and declines for longer oligonucleotides. The entropy gained upon complex formation is approximately 1 entropy unit per released counterion. The method was applied to long DNA molecules (herring and lambda-phage DNA) and revealed that complete dissociation occurs at 750 mM NaCl. Likely contributions of macromolecular desolvation and DNA flexibility to the binding energy are discussed.     

Antifreeze proteins at the ice/water interface: three calculated discriminating properties for orientation of type I proteins

Antifreeze proteins (AFPs) protect many plants and organisms from freezing in low temperatures. Of the different AFPs, the most studied AFP Type I from winter flounder is used in the current computational studies to gain molecular insight into its adsorption at the ice/water interface. Employing molecular dynamics simulations, we calculate the free energy difference between the hydrophilic and hydrophobic faces of the protein interacting with ice. Furthermore, we identify three properties of Type I "antifreeze" proteins that discriminate among these two orientations of the protein at the ice/water interface. The three properties are: the "surface area" of the protein; a measure of the interaction of the protein with neighboring water molecules as determined by the number of hydrogen bond count, for example; and the side-chain orientation angles of the threonine residues. All three discriminants are consistent with our free energy results, which clearly show that the hydrophilic protein face orientations toward the ice/water interface, as hypothesized from experimental and ice/vacuum simulations, are incorrect and support the hypothesis that the hydrophobic face is oriented toward the ice/water interface. The adsorption free energy is calculated to be 2-3 kJ/mol.