Entropy and free energy of a mobile protein loop in explicit water

Free-energy calculations for semi-flexible macromolecules: applications to DNA knotting and looping

We present a method to obtain numerically accurate values of configurational free energies of semiflexible macromolecular systems, based on the technique of thermodynamic integration combined with normal-mode analysis of a reference system subject to harmonic constraints. Compared with previous free-energy calculations that depend on a reference state, our approach introduces two innovations, namely, the use of internal coordinates to constrain the reference states and the ability to freely select these reference states. As a consequence, it is possible to explore systems that undergo substantially larger fluctuations than those considered in previous calculations, including semiflexible biopolymers having arbitrary ratios of contour length L to persistence length P. To validate the method, high accuracy is demonstrated for free energies of prime DNA knots with L/P = 20 and L/P = 40, corresponding to DNA lengths of 3000 and 6000 base pairs, respectively. We then apply the method to study the free-energy landscape for a model of a synaptic nucleoprotein complex containing a pair of looped domains, revealing a bifurcation in the location of optimal synapse (crossover) sites. This transition is relevant to target-site selection by DNA-binding proteins that occupy multiple DNA sites separated by large linear distances along the genome, a problem that arises naturally in gene regulation, DNA recombination, and the action of type-II topoisomerases.

Absolute Free Energy of Binding and Entropy of the FKBP12-FK506 Complex: Effects of the Force Field

The hypothetical scanning molecular dynamics (HSMD) method combined with thermodynamic integration (HSMD-TI) has been extended recently for calculating ΔA(0)-the absolute free energy of binding of a ligand to a protein. With HSMD-TI, ΔA(0) is obtained in a new way as a sum of several components, among them is ΔSligand-the change in the conformational entropy as the ligand is transferred from the bulk solvent to the active site-this entropy is obtained by a specific reconstruction procedure. This unique aspect of HSMD (which is useful in rational drug design) is in particular important for treating large ligands, where ΔSligand might be significant. Technically, one should verify that the results for ΔSligand converge-a property that might become more difficult for large ligands; therefore, studying ligands of increasing size would define the range of applicability of HSMD-TI for binding. In this paper, we check the performance of HSMD-TI by applying it to the relatively large ligand FK506 (126 atoms) complexed with the protein FKBP12, where ΔA(0) = -12.8 kcal/mol is known experimentally as well as the crystal structure of the complex. This structure was initially equilibrated by carrying out a 100 ns molecular dynamics trajectory, where the system is modeled by the AMBER force field, TIP3P water, and Particle Mesh Ewald. HSMD-TI calculations were carried out in three conformational regions defined by the intervals [0.2,2], [2,5], and [5,100] ns along the trajectory, where local equilibration of the total energy has been observed; we obtained ΔA(0) = -13.6 ± 1.1, -16.6 ± 1.4, and -16.7 ± 1.4 kcal/mol, respectively indicating the following: (1) The second and third regions belong to the same conformational subspace of the complex, which is different from the [0.2,2] ns subspace. (2) The unsatisfactory result for ΔA(0) obtained in the well equilibrated (hence theoretically preferred) latter regions reflects the nonperfect modeling used, which however (3) has led to the experimental ΔA(0) in the [0.2,2] ns region close to the crystal structure. Keeping the complex near its crystal structure has been a successful approach in the literature. To check this avenue further, we applied harmonic restraints on backbone atoms and obtained unsatisfactory results for ΔA(0), suggesting that implementation of this approach is not straightforward. Converging results for ΔSligand were obtained in all regions, where the result ΔSligand([0.2,2]) = 7.1 ± 1.2 kcal/mol is less region dependent than ΔA(0) and is relatively large probably due to the large ligand.

Absolute free energy of binding of avidin/biotin, revisited

The binding of biotin to avidin is one of the strongest in nature with absolute free energy of binding, ΔA(0) = -20.4 kcal/mol. Therefore, this complex became a target for a large number of computational studies, which all, however, are based on approximate techniques or simplified models and have led to a wide range of results Therefore, ΔA(0) is calculated here by rigorous statistical mechanical methods and models that consider long-range electrostatics. (1) We apply our method, "hypothetical scanning molecular dynamics with thermodynamic integration" (HSMD-TI) to avidin-biotin modeled by periodic boundary conditions with particle mesh ewald (PME). (2) We apply the double decoupling method (DDM) to this system modeled by the spherical solvent boundary potential (SSBP) and the generalized solvent boundary potential (GSBP). The corresponding results for neutral biotin, ΔA(0) = -29.1 ± 0.8 and -25.2 ± 0.5 kcal/mol are significantly lower than the experimental value; we also provide the result for a charged biotin, ΔA(0) = -33.3 ± 0.8 kcal/mol. It is plausible to suggest that this disagreement with the experiment may stem from ignoring the (positive) contribution of a mobile loop that changes its structure upon ligand binding.

Calculation of the Absolute Free Energy of Binding and Related Entropies with the HSMD-TI Method: The FKBP12-L8 Complex

The hypothetical scanning molecular dynamics (HSMD) method is used here for calculating the absolute free energy of binding, ΔA(0) of the complex of the protein FKBP12 with the ligand SB2 (also denoted L8) - a system that has been studied previously for comparing the performance of different methods. Our preliminary study suggests that considering long-range electrostatics is imperative even for a hydrophobic ligand such as L8. Therefore the system is modeled by the AMBER force field using Particle Mesh Ewald (PME). HSMD consists of three stages applied to both the ligand-solvent and ligand-protein systems. (1) A small set of system configurations (frames) is extracted from an MD trajectory. (2) The entropy of the ligand in each frame is calculated by a reconstruction procedure. (3) The contribution of water and protein to ΔA(0) is calculated for each frame by gradually increasing the ligand-environment interactions from zero to their full value using thermodynamic integration (TI). Unlike the conventional methods, the structure of the ligand is kept fixed during TI, and HSMD is thus free from the end-point problem encountered with the double annihilation method (DAM); therefore, the need for applying restraints is avoided. Furthermore, unlike the conventional methods, the entropy of the ligand and water is obtained directly as a byproduct of the simulation. In this paper, in addition to the difference in the internal entropies of the ligand in the two environments, we calculate for the first time the external entropy of the ligand, which provides a measure for the size of the active site. We obtain ΔA(0) = -10.7 ±1.0 as compared to the experimental values -10.9 and -10.6 kcal/mol. However, a protein/water system treated by periodic boundary conditions grows significantly with increasing protein size and the computation of ΔA(0) would become expensive by all methods. Therefore, we also apply HSMD to FKBP12-L8 described by the GSBP/SSBP model of Roux's group (implemented in the software CHARMM) where only part of the protein and water around the active site are considered and long-range electrostatic effects are taken into account. For comparison this model was also treated by the double decoupling method (DDM). The two methods have led to comparable results for ΔA(0) which are somewhat lower than the experimental value. The ligand was found to be more confined in the active site described by GSBP/SSBP than by PME where its entropy in solvent is larger than in the active site by 1.7 and by 5.5 kcal/mol, respectively.

Relative stability of the open and closed conformations of the active site loop of streptavidin

The eight-residue surface loop, 45-52 (Ser, Ala, Val, Gly, Asn, Ala, Glu, Ser), of the homotetrameric protein streptavidin has a "closed" conformation in the streptavidin-biotin complex, where the corresponding binding affinity is one of the strongest found in nature (ΔG ∼ -18 kcal∕mol). However, in most of the crystal structures of apo (unbound) streptavidin, the loop conformation is "open" and typically exhibits partial disorder and high B-factors. Thus, it is plausible to assume that the loop structure is changed from open to closed upon binding of biotin, and the corresponding difference in free energy, ΔF = F(open) - F(closed) in the unbound protein, should therefore be considered in the total absolute free energy of binding. ΔF (which has generally been neglected) is calculated here using our "hypothetical scanning molecular-dynamics" (HSMD) method. We use a protein model in which only the atoms closest to the loop are considered (the "template") and they are fixed in the x-ray coordinates of the free protein; the x-ray conformation of the closed loop is attached to the same (unbound) template and both systems are capped with the same sphere of TIP3P water. Using the force field of the assisted model building with energy refinement (AMBER), we carry out two separate MD simulations (at temperature T = 300 K), starting from the open and closed conformations, where only the atoms of the loop and water are allowed to move (the template-water and template-loop interactions are considered). The absolute F(open) and F(closed) (of loop + water) are calculated from these trajectories, where the loop and water contributions are obtained by HSMD and a thermodynamic integration (TI) process, respectively. The combined HSMD-TI procedure leads to total (loop + water) ΔF = -27.1 ± 2.0 kcal∕mol, where the entropy TΔS constitutes 34% of ΔF, meaning that the effect of S is significant and should not be ignored. Also, ΔS is positive, in accord with the high flexibility of the open loop observed in crystal structures, while the energy ΔE is unexpectedly negative, thus also adding to the stability of the open loop. The loop and the 250 capped water molecules are the largest system studied thus far, which constitutes a test for the efficiency of HSMD-TI; this efficiency and technical issues related to the implementation of the method are also discussed. Finally, the result for ΔF is a prediction that will be considered in the calculation of the absolute free energy of binding of biotin to streptavidin, which constitutes our next project.

New method for calculating the absolute free energy of binding: the effect of a mobile loop on the avidin/biotin complex

Hypothetical scanning molecular dynamics (HSMD) is a relatively new method for calculating the absolute free energy and entropy. HSMD is extended here for the first time for calculating the absolute free energy of binding, ΔA(0), as applied to the avidin-biotin complex. With HSMD the ligand is built (more accurately reconstructed) from nothing in solvent and in the protein, in contrast to the commonly used methods where the ligand is annihilated (by thermodynamic integration) in these environments. Therefore, the end-point problem encountered with the latter methods does not exist with HSMD and the need for restraints is avoided. Also, the entropy of the ligand and water in both environments is obtained directly as a byproduct of the simulation. The binding mechanism of biotin to avidin involves a mobile loop that is expected to be in an open conformation in unbound avidin, which is changed to a closed one upon binding, that is, the loop moves to cover biotin in the active site. The contribution of the loop's conformational change to the total free energy of binding is calculated here for the first time. Our result, ΔA(0) = -24.9 ± 7 covers the experimental value -20.7 kcal/mol within the error bars.

Entropy and Free Energy of a Mobile Loop Based on the Crystal Structures of the Free and Bound Proteins

A mobile loop changes its conformation from "open" (free enzyme) to "closed" upon ligand binding. The difference in the Helmholtz free energy, ΔF(loop) between these states sheds light on the mechanism of binding. With our "hypothetical scanning molecular dynamics" (HSMD-TI) method ΔF(loop) = F(free) - F(bound) where F(free) and F(bound) are calculated from two MD samples of the free and bound loop states; the contribution of water is obtained by a thermodynamic integration (TI) procedure. In previous work the free and bound loop structures were both attached to the same "template" which was "cut" from the crystal structure of the free protein. Our results for loop 287-290 of AcetylCholineEsterase agree with the experiment, ΔF(loop)~ -4 kcal/mol if the density of the TIP3P water molecules capping the loop is close to that of bulk water, i.e., N(water) = 140 - 180 waters in a sphere of a 18 Å radius. Here we calculate ΔF(loop) for the more realistic case, where two templates are "cut" from the crystal structures, 2dfp.pdb (bound) and 2ace.pdb (free), where N(water) = 40 - 160; this requires adding a computationally more demanding (second) TI procedure. While the results for N(water) ≤ 140 are computationally sound, ΔF(loop) is always positive (18 ± 2 kcal/mol for N(water) = 140). These (disagreeing) results are attributed to the large average B-factor, 41.6 of 2dfp (23.4 Å(2) for 2ace). While this conformational uncertainty is an inherent difficulty, the (unstable) results for N(water) = 160 suggest that it might be alleviated by applying different (initial) structural optimizations to each template.

Methods for calculating the absolute entropy and free energy of biological systems based on ideas from polymer physics

The commonly used simulation techniques, Metropolis Monte Carlo (MC) and molecular dynamics (MD) are of a dynamical type which enables one to sample system configurations i correctly with the Boltzmann probability, P(i)(B), while the value of P(i)(B) is not provided directly; therefore, it is difficult to obtain the absolute entropy, S approximately -ln P(i)(B), and the Helmholtz free energy, F. With a different simulation approach developed in polymer physics, a chain is grown step-by-step with transition probabilities (TPs), and thus their product is the value of the construction probability; therefore, the entropy is known. Because all exact simulation methods are equivalent, i.e. they lead to the same averages and fluctuations of physical properties, one can treat an MC or MD sample as if its members have rather been generated step-by-step. Thus, each configuration i of the sample can be reconstructed (from nothing) by calculating the TPs with which it could have been constructed. This idea applies also to bulk systems such as fluids or magnets. This approach has led earlier to the "local states" (LS) and the "hypothetical scanning" (HS) methods, which are approximate in nature. A recent development is the hypothetical scanning Monte Carlo (HSMC) (or molecular dynamics, HSMD) method which is based on stochastic TPs where all interactions are taken into account. In this respect, HSMC(D) can be viewed as exact and the only approximation involved is due to insufficient MC(MD) sampling for calculating the TPs. The validity of HSMC has been established by applying it first to liquid argon, TIP3P water, self-avoiding walks (SAW), and polyglycine models, where the results for F were found to agree with those obtained by other methods. Subsequently, HSMD was applied to mobile loops of the enzymes porcine pancreatic alpha-amylase and acetylcholinesterase in explicit water, where the difference in F between the bound and free states of the loop was calculated. Currently, HSMD is being extended for calculating the absolute and relative free energies of ligand-enzyme binding. We describe the whole approach and discuss future directions.

Dynamic properties of a psychrophilic alpha-amylase in comparison with a mesophilic homologue

The cold-active, chloride-dependent alpha-amylase from Pseudoalteromonas haloplanktis (AHA) is one of the best characterized psychrophilic enzymes, and shares high sequence and structural similarity with its mesophilic porcine counterpart (PPA). An atomic detail comparative analysis was carried out by performing more than 60 ns of multiple-replica explicit-solvent molecular dynamics simulations on the two enzymes in order to characterize the differences in ensemble properties and dynamics in solution between the two homologues. We find in both enzymes high flexibility clusters in the surroundings of the substrate-binding groove, primarily involving the long loops that protrude from the main domain's barrel structure. These loops are longer in PPA and extend further away from the core of the barrel, where the active site is located: essential fluctuations in PPA mainly affect the highly solvent-accessible portions of these loops, whereas AHA is characterized by greater flexibility in the immediate surroundings of the active site. Furthermore, detailed analysis of active-site dynamics has revealed that elements previously identified through X-ray crystallography as involved in substrate binding in both enzymes undergo concerted motions that may be linked to catalysis.

Toward accurate microscopic calculation of solvation entropies: extending the restraint release approach to studies of solvation effects

The evaluation of the solvation entropies is a major conceptual and practical challenge. On the one hand, it is interesting to quantify the factors that are responsible for the solvation entropies in solutions, whereas on the other, it is essential to be able to assess the contributions of the solvation entropies to the binding free energies and related properties. In fact, the solvation entropies are neglected in almost all of the studies of the binding entropies. The main problem is that widely used approaches, such as the quasiharmonic (QH) approximation, do not provide reliable results particularly in cases of shallow potential and multidimensional surfaces while brute force evaluations of the entropic effects by simulating temperature dependence of the free energy converges very slowly. This paper addresses the above issue by starting with an analysis of the factors that are responsible for the negative solvation entropy of ions, showing that it is not due to the change in the solvent vibration modes or to the solvent force constant but to the changes in the solvent configurational space upon change in the solute charges. We begin by clarifying that when one deals with aqueous solutions, it is easy to evaluate the corresponding entropic effect by the Langevin dipole (LD) treatment. However, in this work we are interested in developing a general microscopic tool that can be used to study similar effects in the proteins. To this end, we explore the ability of our restraint release (RR) approach to evaluate the solvation entropy. We start this analysis by reviewing the foundation of this approach and in particular, the requirements of minimizing the enthalpy contribution to the RR free energy. We then establish that our approach is not a specialized harmonic treatment but a rather powerful general approach. Moving to the main topic of this work, we demonstrate that the RR approach provides quantitative results for the solvation entropies of monovalent and divalent ions and effectively captures the physics of these entropic effects. The success of the current approach indicates that it should be applicable to the studies of the solvation entropies in the proteins and also, in examining hydrophobic effects. Thus, we believe that the RR approach provides a powerful tool for evaluating the corresponding contributions to the binding entropies and, eventually, to the binding free energies. This holds promise for extending the information theory modeling to proteins and protein-ligand complexes in aqueous solutions and consequently, facilitating computer-aided drug design.

Absolute free energy and entropy of a mobile loop of the enzyme acetylcholinesterase

The loop 287-290 (Ile, Phe, Arg, and Phe) of the protein acetylcholinesterase (AChE) changes its structure upon interaction of AChE with diisopropylphosphorofluoridate (DFP). Reversible dissociation measurements suggest that the free-energy (F) penalty for the loop displacement is DeltaF=Ffree-Fbound approximately -4 kcal/mol. Therefore, this loop has been the target of two studies by Olson's group for testing the efficiency of procedures for calculating F. In this paper, we test for the first time the performance of our "hypothetical scanning molecular dynamics" (HSMD) method and the validity of the related modeling for a loop with bulky side chains in explicit water. Thus, we consider only atoms of the protein that are the closest to the loop (they constitute the "template"), where the rest of the atoms are ignored. The template's atoms are fixed in the X-ray coordinates of the free protein, and the loop is capped with a sphere of TIP3P water molecules; also, the X-ray structure of the bound loop is attached to the free template. We carry out two separate MD simulations starting from the free and bound X-ray structures, where only the atoms of the loop and water are allowed to move while the template-water and template-loop (AMBER) interactions are considered. The absolute Ffree and Fbound (of the loop and water) are calculated from the corresponding trajectories. A main objective of this paper is to assess the reliability of this model, and for this several template sizes are studied capped with 80-220 water molecules. We find that consistent results for the free energy (which also agree with the experimental data above) require a template larger than a minimal size and a number of water molecules approximately equal to the experimental density of bulk water. For example, we obtain DeltaFtotal=DeltaFwater+DeltaFloop=-3.1+/-2.5 and -3.6+/-4 kcal/mol for a template consisting of 944 atoms and a sphere containing 160 and 180 waters, respectively. Our calculations demonstrate the important contribution of water to the total free energy. Namely, for water densities close to the experimental value, DeltaFwater is always negative leading thereby to a negative DeltaFtotal (while DeltaFloop is always positive). Also, the contribution of the water entropy TDeltaSwater to DeltaFtotal is significant. Various aspects related to the efficiency of HSMD are tested and improved, and plans for future studies are discussed.

Methods for calculating the entropy and free energy and their application to problems involving protein flexibility and ligand binding

The Helmholtz free energy, F and the entropy, S are related thermodynamic quantities with a special importance in structural biology. We describe the difficulties in calculating these quantities and review recent methodological developments. Because protein flexibility is essential for function and ligand binding, we discuss the related problems involved in the definition, simulation, and free energy calculation of microstates (such as the alpha-helical region of a peptide). While the review is broad, a special emphasize is given to methods for calculating the absolute F (S), where our HSMC(D) method is described in some detail.