Home | About | Journals | Submit | Contact Us | Français |

**|**HHS Author Manuscripts**|**PMC2738853

Formats

Article sections

- Abstract
- I. Introduction
- II. Physical Considerations
- III. Using the RR Approach
- IV. RESULTS
- V. DISCUSSION
- CONCLUSIONS AND FUTURE APPLICATIONS
- References

Authors

Related links

J Phys Chem B. Author manuscript; available in PMC 2010 May 21.

Published in final edited form as:

PMCID: PMC2738853

NIHMSID: NIHMS114340

The publisher's final edited version of this article is available at J Phys Chem B

See other articles in PMC that cite the published article.

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, while 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 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 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.

Entropic effects play a crucial role in biological and chemical processes contributing to folding,^{1}^{–}^{6} binding processes ^{3}^{,}^{7}^{–}^{11} and sometimes, modulating catalytic effects.^{12}^{–}^{15} Therefore, it is important to be able to assess the magnitude of the entropic contributions in a quantitative way. Unfortunately, this task is extremely challenging and popular approaches such as the quasiharmonic (QH) and related treatments are quite inaccurate.^{16}^{,}^{17} Several approaches have been proposed in the recent years for more accurate calculation of configurational entropies^{1}^{,}^{6}^{,}^{12}^{,}^{13}^{,}^{17}^{–}^{27} and of these, some have been validated in the actual studies of biological systems.^{12}^{,}^{13}^{,}^{17} Furthermore, it has been shown that the entropic effects can be evaluated by brute force approaches but these require an extremely large computer time.^{15}^{,}^{17}^{,}^{28} Interestingly however, most of the attempts to study the entropic contributions to binding and catalysis conveniently ignored the effect of the solvation entropies.^{10}^{,}^{17}^{,}^{29}^{–}^{32} In some case s, it was found that the solvation entropy for the given process in water can be determined quite accurately by the phenomenological approach of Florian and Warshel^{33} who calibrated the polarization of their Langevin dipole(LD) model to reproduce the solvation entropies (see below). While the LD approach is extremely effective in solution studies, it is not clear how to implement this type of approach in the studies of the solvation entropies in the proteins. In such a case, one needs a fully microscopic approach and the development of a proper strategy which is the subject of the present work. This paper begins by considering the physical basis of the solvation entropies. It is shown that the main contributions are not associated with the solvent vibrations but with the solvent configurations. We also explore the ability of the restraint release (RR) approach to evaluate the solvation entropies. It is found that the RR approach provides an effective tool for the quantitative evaluation of the solvation energies of the charged molecules.

In studying the solvation entropies, it is important to have a working hypothesis regarding the origin of these effects. The simplest idea is that somehow the vibrational frequencies of the first solvation shell are higher near uncharged solute than near the solute in its charged form. One way to explore this idea is to look at the solvent vibrations that are associated with the solute-solvent interactions whilst the solute is charged and uncharged. This type of analysis can be done conveniently by the dispersed polaron (Spin Boson) method^{34} which provides the projection of the solvent modes along the solvation coordinates. The results of such an analysis for sodium and calcium ion in water gave the results summarized in Fig 1. Evidently, the vibrational frequencies are similar and thus, unlikely to account for the solvation entropy.

The Fourier transform of the time dependence of the energy gap (Δ*V*_{ab}) between the charged (state a) and uncharged (state b) states of **(a)** Na^{+ }ion **(b)** Ca^{2+} ions. The dimensionless coordinate displacement **(**Δ**)** values are obtained from a simple **...**

Attempts to use the calculated spectrum in the evaluation of the vibrational entropy *S*(*v*), Δ*S _{vib}* =∫

The above finding tells us to look in the same direction as explored in the LD treatment. That is, it is more likely that the solvation entropies reflect the tendency of the dipoles of the first solvation shell to be oriented toward the charged solute and to assume random conformation when the solute is uncharged. This has been in fact the idea beyond our LD treatment^{33} which will be discussed below.

According to the LD model, the hydration entropy for charged solutes is dominated by the solvent immobilization in the presence of the solute electrostatic field, *ξ*. This entropy contribution can be estimated using the concept of accessible configurational volumes introduced by Frank and Evans^{37}^{,}^{38} and expressed as

$$\mathrm{\Delta}{S}_{\mathit{immob}}=\sum _{j}{k}_{B}ln({V}_{sl}^{j}/{V}_{pl}^{j}),\phantom{\rule{0.38889em}{0ex}}j\in \text{solvent}\phantom{\rule{0.16667em}{0ex}}\text{dipoles},$$

(1)

where *k _{B}* is the Boltzmann constant while

Schematic representation of the restricted configurational volume of a grid-centered solvent dipole in the presence of an external electrostatic field, *ξ*. In the pure liquid, the rotation of the solvent dipole is assumed to be free and the accessible **...**

$${V}_{sl}=\frac{2\alpha}{\pi}{V}_{pl}$$

(2)

The angle *α* corresponds to the average angle between the direction of the solvent dip ole and the external field. In the LD solvation model, it is assumed that the average solvent polarization is determined even in solution by the Langevin function L(*x*)^{39}:

$$cos\phantom{\rule{0.38889em}{0ex}}\alpha =L(x)=coth(x)-1/x$$

(3)

$$\text{where}\phantom{\rule{0.16667em}{0ex}}x=\frac{{\mu}_{0}\xi}{3{k}_{B}T}C$$

(4)

In Eq 4, *μ*_{0}, *ξ*, *k _{B}*,

$$\mathrm{\Delta}{S}_{\mathit{immob}}^{j}={k}_{B}ln\left(\frac{\pi}{2arccos(L(x))}\right)$$

(5)

Eq 5 does not include explicitly the dipole-dipole interaction term. Thus, this equation loses its validity for solvent dipoles positioned in the regions where the magnitude of the solute electrostatic field is comparable to the field from other solvent dipoles. Consequently, we introduced an empirical correction function

$$\mathrm{\Delta}{S}_{\mathit{cor}}^{j}=-5{k}_{B}arctan\left(\frac{x}{27}\right)exp\left(\frac{{x}^{6}}{32}\right)$$

(6)

that decreases the magnitude of Δ*S _{immob}* for small fields. For +1 or −1 charged solutes, this correction affects the solvent dipoles beyond the first hydration shell, whereas for +2 and −2 charged solutes this term affects dipoles outside the first two hydration shells. Finally, a multiplicative constant,

$$\mathrm{\Delta}{\stackrel{\sim}{S}}_{\mathit{immob}}={C}_{i}\mathrm{\sum}(\mathrm{\Delta}{S}_{\mathit{immob}}+\mathrm{\Delta}{S}_{\mathit{cor}}),j\in \text{solvent}\phantom{\rule{0.16667em}{0ex}}\text{dipoles}$$

(7)

The behavior of the uncorrected (Eq 5) and the corrected (Eq 7) immobilization entropy as a function of *x* is presented in Fig. 2 of ref^{33}. Finally, the total hydration entropies are calculated as the sum of the hydrophobic and immobilization entropies

$$\mathrm{\Delta}\phantom{\rule{0.16667em}{0ex}}{S}_{\mathit{hydr}}=\mathrm{\Delta}\phantom{\rule{0.16667em}{0ex}}{S}_{\mathit{phob}}+\mathrm{\Delta}{\stackrel{\sim}{S}}_{\mathit{immob}}$$

(8)

This expression can be applied to both neutral and ionic solutes. Note that it has been suggested that the hydration entropy of large positive ions may be significantly different from that of large negative ions of the same size.^{42}^{,}^{43} In our model, this entropy does not depend on the sign of the solute charge since we believe that the magnitude of this dependence is still an open question.

It should be mentioned here that Δ*S _{phob}* could have also been expressed by considering the decrease in the number of configurations available for the solvent dipoles near the nonpolar solute surface relative to the configurations available in the bulk solvent.

In order to explore the above idea, we evaluated the polarization of the water molecules around Na^{+}, K^{+} and Ca ^{2+}ions and around the uncharged forms of these ions. The corresponding results are summarized in Figs 3–6.

The time-dependent fluctuations of three water molecules at distances of 1.96Å (a), 4.90Å (b) and 6.21Å, respectively (c) from charged (left) and uncharged (right) Na^{+} ion during a 2ns molecular dynamics simulation.

Histogram of the average orientation of three water molecules at distances of 1.96Å (a), 4.90Å (b) and 6.21Å, respectively (c) from charged (left) and uncharged (right) Ca^{2+} ion during a 2ns molecular dynamics simulation.

Figs 3 and and44 describe the fluctuations of the angle between the solvent field and the bisector of the two O-H vectors of the representative water molecules. Figs 5 and and66 provide histogram of the angle, θ. As expected, the figures indicate that the orientation of the molecules in the first solvation shell is more restricted in the charged solute than in the neutral solute. However, the trend in the case of molecules in the second and third solvation shell is less clear.

The time-dependent fluctuations of three water molecules at distances of 1.96Å (a), 4.90Å (b) and 6.21Å, respectively (c) from charged (left) and uncharged (right) Ca^{2+} ion during a 2ns molecular dynamics simulation.

Histogram of the average orientation of three water molecules at distances of 1.96Å (a), 4.90Å (b) and 6.21Å, respectively (c) from charged (left) and uncharged (right) Na^{+} ion during a 2ns molecular dynamics simulation.

While some qualitative aspects of the origin of the entropic effect are observed in Figs 3–6, we need a more quantitative way that can also be “imported” to similar studies in the proteins. Thus, we turn our attention to the RR approach which is described in the next section.

The RR approach is based on the idea that with strong harmonic Cartesian restraints we have no entropy effects so, in some way, the free energy of releasing the restraints can tell us about the magnitude of the entropic effect. This general idea has been proposed by Herman and Wang^{45} for evaluating entropic effects in limited internal coordinates and is similar in some respects to that suggested in an early proposal.^{46} However, the problem is to formulate the RR idea correctly in a general set of coordinates and also, to take into account the fact recognized in ref ^{47}, that the RR free energy can include a large enthalpic contribution. The resolution of this problem and the introduction of a practical general RR approach have been introduced in refs^{12}^{,}^{47} with a general cartesian treatment and an optimization approach that eliminates the enthalpic contribution from the RR free energy. The main points of this approach and the additional recent improvements will be discussed below.

In order to obtain the solute entropic contribution to a transfer between two potential surfaces (*U*_{I}→*U*_{II}), one has to separate the corresponding partition functions to the solute and the solvent contributions (where in the present case, the ‘solute’ includes explicit solvent molecules and the ‘solvent’ is defined as the rest of the system) by writing

$$\begin{array}{l}{Q}_{N}=\int dR\phantom{\rule{0.16667em}{0ex}}\int dr\phantom{\rule{0.16667em}{0ex}}{e}^{-{U}_{N}(R,r)\beta}=\int d\overline{R}\phantom{\rule{0.16667em}{0ex}}[\int dR\int dr\delta (R-\overline{R})\phantom{\rule{0.16667em}{0ex}}{e}^{-{U}_{N}(R,r)\beta}]\\ =\int d\overline{R}{q}_{N}\phantom{\rule{0.16667em}{0ex}}(\overline{R})=\int d\overline{R}\phantom{\rule{0.16667em}{0ex}}{e}^{-{W}_{N}(\overline{R})\beta}\end{array}$$

(9)

and

$${q}_{N}=\int dR\int dr\delta (R-\overline{R})\phantom{\rule{0.16667em}{0ex}}{e}^{-{U}_{N}\beta}$$

where *U _{N}* is the potential surface for state

$${W}_{N}(\overline{R})={U}_{N}(\overline{R})+\mathrm{\Delta}{g}_{\mathit{sol}}^{N}(\overline{R})={W}_{N}^{0}({\overline{R}}_{0}^{N})+{({K}_{W}^{N}/2(\overline{R})-{\overline{R}}_{0}^{N})}^{2}$$

(10)

where we could use *R* but instead, took as our variable to indicate that this is our restraint coordinate. Here, Δ*g _{sol}* () is the solvation free energy at ,

For the sake of simplicity, we continue by examining a case where the PMF corresponds to a one-dimensional solute coordinate, although extension to many dimensions is straightforward. Our one-dimensional case is described by the thermodynamic cycle of Fig 7. In this cycle, we divide *R* to segments that can be defined by the delta function of Eq 9 or by introducing a strong quadratic constraint ((*K*′_{1}/2)( − ′)^{2}). In this case, we have

A schematic illustration of the effect of and the enthalpic contribution to Δ*G*′. The figure considers one-dimensional potentials of mean force (*W*) for two systems (e.g. charged and uncharged solute) and the effect of **...**

$$\begin{array}{l}{{W}^{\prime}}_{N}(\overline{R})={W}_{N}^{0}({\overline{R}}_{0}^{N})+({K}_{W}/2){(\overline{R}-{\overline{R}}_{0}^{N})}^{2}+({{K}^{\prime}}_{1}/2){(\overline{R}-{{\overline{R}}^{\prime}}_{N})}^{2}\\ ={W}_{N}({\overline{R}}^{\prime})+({K}_{1}/2){(\overline{R}-{{\overline{R}}^{\prime}}_{N})}^{2}\end{array}$$

(11)

where the notation *W*′ designates the addition of the constraint.

Next, we use Eq 11 and our quadratic constraint and evaluate the relevant partition functions.

$${Q}_{N}=\underset{-\infty}{\overset{\infty}{\int}}{e}^{-{W}_{N}(\overline{R})}d\overline{R}={e}^{-{W}_{N}({\overline{R}}_{0})\beta}\underset{-\infty}{\overset{\infty}{\int}}d\overline{R}{e}^{-({{K}_{W}}^{N}/2){(\overline{R}-{\overline{R}}_{0}^{N})}^{2}\beta}={e}^{-{W}_{N}^{0}({\overline{R}}_{0})\beta}\sqrt{2\pi /\beta {K}_{W}^{N}}$$

(12)

$${{Q}^{\prime}}_{N}({\overline{R}}^{\prime})=\underset{-\infty}{\overset{\infty}{\int}}{e}^{-{{W}^{\prime}}_{N}(\overline{R})\beta}d\overline{R}={e}^{-{W}_{N}({{\overline{R}}^{\prime}}_{N})\beta}\underset{-\infty}{\overset{\infty}{\int}}d{\overline{R}}^{\prime}{e}^{-({K}_{1}/2){(\overline{R}-{{\overline{R}}^{\prime}}_{N})}^{2}\beta}={e}^{-{W}_{N}({{\overline{R}}^{\prime}}_{N})\beta}\sqrt{2\pi /\beta {K}_{1}}$$

where we use *Q*′ or the partition function evaluated with *W*′* _{N}* and where
${W}_{N}^{0}={W}_{N}({\overline{R}}_{0}^{N})$. Using the above equation we can write

$$\mathrm{\Delta}{G}_{\text{I}\to \text{II}}=-{\beta}^{-1}ln({Q}_{\text{II}}/{Q}_{\text{I}})=\mathrm{\Delta}{W}_{\text{I}\to \text{II}}^{0}-\frac{{\beta}^{-1}}{2}ln({K}_{W}^{\text{I}}/{K}_{W}^{\text{II}})$$

$$\begin{array}{l}\mathrm{\Delta}{{G}^{\prime}}_{{\text{I}}^{\prime}\to \text{I}}({\overline{R}}^{\prime})={G}_{\text{I}}-{{G}^{\prime}}_{\text{I}}({\overline{R}}^{\prime})=-{\beta}^{-1}ln({Q}_{\text{I}}/{{Q}^{\prime}}_{\text{I}}({\overline{R}}^{\prime}))\\ =-\frac{{\beta}^{-1}}{2}ln({K}_{1}/{K}_{W}^{\text{I}})+({W}_{\text{I}}({\overline{R}}_{0}^{\text{I}})-{W}_{\text{I}}({{\overline{R}}^{\prime}}_{\text{I}}))\le 0\end{array}$$

(13)

$${\overline{E}}_{I}={k}_{B}{T}^{2}{(\partial ln{Q}_{\text{I}}/\partial T)}_{v}={W}_{\text{I}}^{0}+({k}_{B}/2)T-(T\phantom{\rule{0.16667em}{0ex}}\partial {W}_{\text{I}}^{0}/\partial T)={W}_{\text{I}}^{0}+({k}_{B}/2)T-T\mathrm{\Delta}{S}_{\mathit{sol}}^{\text{I}}({\overline{R}}_{0}^{\text{I}})$$

where is the average energy and *S _{sol}* is the contribution from the surrounding of the specified system included in our explicit treatment. We also have

$${S}_{\text{I}}={k}_{B}ln{Q}_{\text{I}}+{\overline{E}}_{\text{I}}/T=({k}_{B}/2)[1+ln(2\pi /\beta {K}_{W}^{\text{I}})]+\mathrm{\Delta}{S}_{\mathit{sol}}^{\text{I}}({\overline{R}}_{0}^{\text{I}})$$

(14)

$${S}_{{\text{I}}^{\prime}}({\overline{R}}^{\prime})=({k}_{B}/2)[1+ln(2\pi /\beta {K}_{\text{I}})]+\mathrm{\Delta}{S}_{\mathit{sol}}^{\text{I}}({{\overline{R}}^{\prime}}_{\text{I}})$$

Now, we can write

$$\begin{array}{l}\mathrm{\Delta}{S}_{\text{I}\to \text{II}}=({k}_{B}/2)ln({K}_{W}^{I}/{K}_{W}^{\text{II}})+\mathrm{\Delta}{S}_{\mathit{sol}}^{\text{II}}({\overline{R}}_{0}^{\text{II}})-\mathrm{\Delta}{S}_{\mathit{sol}}^{\text{I}}({\overline{R}}_{0}^{\text{I}})\\ =({k}_{B}/2)ln({K}_{W}^{\text{I}}/{K}_{W}^{\text{II}})+\mathrm{\Delta}{S}_{\mathit{sol}}^{\text{I}\to \text{II}}\\ =\mathrm{\Delta}{S}_{\text{I}\to \text{II}}+\mathrm{\Delta}{S}_{\mathit{sol}}^{\text{I}\to \text{II}}\end{array}$$

(15)

$$\mathrm{\Delta}{{S}^{\prime}}_{{\text{I}}^{\prime}\to \text{I}}({{\overline{R}}^{\prime}}_{\text{I}})=({k}_{B}/2)ln({{K}^{\prime}}_{\text{I}}/{K}_{W}^{\text{I}})+\mathrm{\Delta}{S}_{\mathit{sol}}^{\text{I}}({{\overline{R}}^{\prime}}_{0})-\mathrm{\Delta}{S}_{\mathit{sol}}^{\text{I}}({{\overline{R}}^{\prime}}_{\text{I}})$$

(16)

Thus,

$$\mathrm{\Delta}{{S}^{\prime}}_{{\text{I}}^{\prime}\to \text{II}}({\overline{R}}_{0}^{\text{I}})=({k}_{B}/2)ln({{K}^{\prime}}_{\text{I}}/{K}_{W}^{\text{I}})$$

(17)

Finally, from Eq 13 and Eq 15 we have in the special case when ${{\overline{R}}^{\prime}}_{\text{I}}={\overline{R}}_{0}^{\text{I}}$ and ${{\overline{R}}^{\prime}}_{\text{II}}={\overline{R}}_{0}^{\text{II}}$

$$\mathrm{\Delta}{{G}^{\prime}}_{\text{I}{\text{I}}^{\prime}\to \text{II}}({\overline{R}}_{0}^{\text{II}})-\mathrm{\Delta}{{G}^{\prime}}_{{\text{I}}^{\prime}\to \text{I}}({\overline{R}}_{0}^{\text{I}})=-({k}_{B}T/2)ln({K}_{W}^{\text{I}}/{K}_{W}^{\text{II}})=-T\mathrm{\Delta}{{S}^{\prime}}_{\text{I}\to \text{II}}$$

(18)

The usage of this relationship however, requires knowing the values of
${\overline{R}}_{0}^{\text{I}}$ and
${\overline{R}}_{0}^{\text{II}}$. This can be done, in principle, by minimizing the corresponding PMF (see section III. *4*). Furthermore, using Eq 13 we can write

$$\mid \mathrm{\Delta}{{G}^{\prime}}_{{\text{I}}^{\prime}\to \text{I}}({\overline{R}}_{0})\mid \le \mid \mathrm{\Delta}{{G}^{\prime}}_{\text{I}\to {\text{I}}^{\prime}}(\overline{R})\mid $$

(19)

Thus, we estimate Δ*G*′(*R*_{0}) by trying different ’s and taking the Δ*G*′() with the smallest absolute value. Now, we can use Eqs 17 and 18 to evaluate Δ*S*′_{I→II} and
$\mathrm{\Delta}{{S}^{\prime}}_{{\text{I}}^{\prime}\to \text{II}}({\overline{R}}_{0}^{\text{I}})$. The that minimizes |Δ*G*′| is also used in evaluating the
$\mathrm{\Delta}{S}_{\mathit{sol}}^{\text{I}\to \text{II}}$ of Eq 17 and the total entropy Δ*S*_{I→II}.

The validity of the inequality of Eq. 19 has been established above for harmonic systems but probably it is also valid for an harmonic potentials. Here, we only demonstrate the validity in some limiting cases. We first consider a square well potential (Fig. 8)

In this case, we can define a unit length Δ*R* by

$$\mathrm{\Delta}R=\sqrt{2\pi /\beta {K}_{1}}$$

(20)

and then we have,

$${Q}_{I}=n\mathrm{\Delta}{Re}^{-{W}_{0}^{I}\beta}$$

(21)

where *n* = *L*/Δ*R* and *L* is the width of the potential well. This leads to

$$\mathrm{\Delta}{G}_{I}=-{\beta}^{-1}ln{Q}_{I}=-{\beta}^{-1}ln(n\mathrm{\Delta}R)+{W}_{0}^{I}$$

(22)

$${\overline{E}}_{I}={k}_{B}{T}^{2}(\partial ln{Q}_{I}/\partial T)=-{k}_{B}{T}^{2}\partial ({W}_{0}^{I}/{k}_{B}T)/\partial T={W}_{0}^{I}$$

$${S}_{I}={k}_{B}ln{Q}_{I}+{\overline{E}}_{I}=({\beta}^{-1}/T)[ln({n}_{I}\mathrm{\Delta}R)-{W}_{0}^{I}\beta )+{W}_{0}^{I}/T={k}_{B}ln({n}_{I}\mathrm{\Delta}R)$$

Now, we also have

$${{Q}^{\prime}}_{I}=\int {e}^{-({K}_{1}/2)\phantom{\rule{0.16667em}{0ex}}{(R-{R}_{0}^{1})}^{2}\beta}{e}^{-{W}_{0}^{I}\beta}={e}^{-{W}_{0}^{I}\beta}\sqrt{2\pi /\beta {K}_{1}}={e}^{-{W}_{0}^{I}\beta}\mathrm{\Delta}R$$

(23)

$$\mathrm{\Delta}{{G}^{\prime}}_{I}={W}_{0}^{I}-{\beta}^{-1}ln\mathrm{\Delta}R$$

(24)

$$\mathrm{\Delta}{G}_{{I}^{\prime}-I}=-{\beta}^{-1}(ln({n}_{1}\mathrm{\Delta}R)-ln\mathrm{\Delta}R)$$

Using the above relationship for transition from *I* to *II*, with the same *K*_{1}, we obtain

$$\mathrm{\Delta}{S}_{I\to II}\cong {k}_{B}ln({n}_{II}/{n}_{I})=\mathrm{\Delta}{G}_{I{I}^{\prime}\to II}-\mathrm{\Delta}{G}_{{I}^{\prime}\to II}$$

(25)

as long as we are in the region where *R* ≥ 0 and *W* = *W*_{0}, and under the condition that the contributions for Δ*G _{I}* from the region with
$W={W}_{1}^{I}$ are small. Once we move to the region where
$W={W}_{1}^{I}$, where
${W}_{1}^{I}\gg {W}_{0}^{I}$ (e.g. at

$$\mathrm{\Delta}{{G}^{\prime}}_{I}({{R}^{\prime}}_{1}<0)={W}_{1}^{I}-{\beta}^{-1}ln\mathrm{\Delta}R$$

(26)

$$\mathrm{\Delta}{{G}^{\prime}}_{{I}^{\prime}\to I}({{R}^{\prime}}_{1}<0)={W}_{1}^{I}-{\beta}^{-1}ln({n}_{I})$$

Hence,

$$\mid \mathrm{\Delta}{{G}^{\prime}}_{{I}^{\prime}\to I}(0<{{R}^{\prime}}_{1}<L)\mid <\mid \mathrm{\Delta}{{G}^{\prime}}_{{I}^{\prime}\to I}({{R}^{\prime}}_{1}<0)\mid $$

thus, satisfying our condition of obtaining Δ*S _{I}*

Schematic representation of a periodic system that may be utilized to obtain Δ*S*_{I→II} from *R*′_{1} that minimizes Δ*G*′.

In this case

$${W}_{I}(R)={W}_{I}^{0}+\sum _{l}(K/2){(R-{R}_{l}^{0})}^{2}\theta {(R)}_{l}$$

(27)

where *θ*(*R*)* _{l}* = 0 when

$${Q}_{I}=\sum _{l}\underset{{R}_{l}^{-}}{\overset{{R}_{l}^{-}}{\int}}exp\{-({W}_{II}^{0}+(K/2){(R-{R}_{l}^{0})}^{2})\beta \}dR$$

(28)

$$=N\underset{-\mathrm{\Delta}/2}{\overset{+\mathrm{\Delta}/2}{\int}}exp\{-({W}_{I}^{0}+(K/2){(R-{R}_{l}^{0})}^{2})\beta \}dR\cong Nexp\{-{W}_{I}^{0}\beta \}\sqrt{2\pi /K\beta}$$

where $\mathrm{\Delta}={R}_{l}^{+}-{R}_{l}^{-}$. Now, we can write

$$\mathrm{\Delta}{G}_{I}=-{\beta}^{-1}ln{Q}_{I}={W}_{I}^{0}-{\beta}^{-1}ln(\sqrt{2\pi /K\beta})-{\beta}^{-1}ln{N}_{I}$$

(29)

We also, obtain

$$ln{Q}^{\prime}\cong {W}_{I}({{R}^{\prime}}_{l})-{\beta}^{-1}ln\sqrt{2\pi /\beta {K}_{1}}$$

(30)

where *R*′* _{l}* is the center of the constraint potential. We also, have

$${\overline{E}}_{I}={k}_{B}{T}^{2}(\partial ln{Q}_{I}/\partial T)\cong {W}_{I}^{0}$$

(31)

$${S}_{I}={k}_{B}ln{Q}_{I}+\overline{E}/T\cong -{\beta}^{-1}(ln(\sqrt{2\pi /\beta K})+ln{N}_{I})$$

$${S}_{I\to II}={k}_{B}ln({N}_{II}/{N}_{I})$$

Following the same consideration as in the previous case, and utilizing the fact that ${W}_{I}({{R}^{\prime}}_{l})>{W}_{I}^{0}$ when ${{R}^{\prime}}_{l}\ne {R}_{l}^{0}$, we will obtain the same inequality.

The RR treatment is based on starting the simulations with a strong positional restraint where the entropic contributions are negligible. However, the free energy perturbation (FEP) release requires an enormous amount of computer time even where *K*_{1} is significantly large. Fortunately, we realized that the RR procedure is not really needed at the stage when *K*_{1} is large since the corresponding entropic contribution can be easily estimated by the QH approximations.^{48}^{,}^{49} In general, the QH tends to be valid when restraints are significant and only starts to become problematic when *K*_{1} is small (resulting in a range of very shallow and an harmonic potential energy surfaces).^{16}^{,}^{17}^{,}^{47} We exploit this fact by combining the QH and RR approximations by evaluating the entropic contribution at a moderate value of *K*_{1} by the QH approach followed by the RR procedure from this specific value. This treatment can be expressed by

$$-T\mathrm{\Delta}{S}_{\mathit{conf}}=-T\mathrm{\Delta}S{(K={K}_{1})}_{QH}+min[\mathrm{\Delta}{G}_{RR}(K={K}_{1}\to K=0)]$$

(32)

where the −*T*Δ*S*(*K* = *K*_{1})* _{QH}* designates the entropy computed by the QH approximation, where

Our next task is to optimize Δ*G _{RR}*. This free energy term is obtained by a FEP/umbrella sampling approach using a mapping potential

$${U}_{m}={U}^{(1)}(1-{\lambda}_{m})+{U}^{(2)}{\lambda}_{m}$$

(33)

where

$${U}^{(1)}={U}_{\mathit{syst}}+({K}^{(1)}/2)\sum _{j}{({R}_{j}-{\overline{R}}_{j})}^{2}$$

(34)

$${U}^{(2)}={U}_{\mathit{syst}}$$

where *U _{syst}* is the potential energy of the unconstraint system.

If we use the end point linear response approximation (LRA)instead of the full FEP treatment, we will have^{10}

$$\mathrm{\Delta}{G}_{RR}\approx \frac{1}{2}[{\langle {U}^{(2)}-{U}^{(1)}\rangle}_{{U}_{2}}+{\langle {U}^{(2)}-{U}^{(1)}\rangle}_{{U}_{1}}]=\frac{1}{2}[{\langle \mathrm{\Delta}U\rangle}_{{U}_{2}}+{\langle \mathrm{\Delta}U\rangle}_{{U}_{1}}]$$

(35)

Now we have to find *R*^{0} that minimizes Δ*G _{RR}*. This is done by

$$\partial \mathrm{\Delta}{G}_{RR}/\partial {\overline{R}}_{l}\cong (1/2)\partial ({\langle \mathrm{\Delta}U\rangle}_{{U}_{2}}+{\langle \mathrm{\Delta}U\rangle}_{{U}_{1}})/\partial {\overline{R}}_{l}$$

Now,

$$\begin{array}{l}\partial ({\langle \mathrm{\Delta}U\rangle}_{{U}_{1}})/\partial {\overline{R}}_{l}=-{K}^{(1)}{\langle {R}_{l}-{\overline{R}}_{l}\rangle}_{{U}_{1}}\\ \partial ({\langle \mathrm{\Delta}U\rangle}_{{U}_{2}})/\partial {\overline{R}}_{l}=0\\ {(\nabla \mathrm{\Delta}{G}_{RR})}_{l}=\partial \mathrm{\Delta}{G}_{RR}/\partial {\overline{R}}_{l}=-1/2{K}^{(1)}[{\langle ({R}_{l}-{\overline{R}}_{l})\rangle}_{{U}_{1}}]\\ {\partial}^{2}{\langle \mathrm{\Delta}U\rangle}_{{U}_{1}}/\partial {\overline{R}}_{l}\phantom{\rule{0.16667em}{0ex}}\partial {\overline{R}}_{{l}^{\prime}}={\delta}_{l{l}^{\prime}}{K}^{(1)}/2\end{array}$$

(36)

Thus,

$$\partial \mathrm{\Delta}{G}_{RR}/\partial {\overline{R}}_{l}\partial {\overline{R}}_{{l}^{\prime}}={F}_{l{l}^{\prime}}={K}^{(1)}{\delta}_{l{l}^{\prime}}/2$$

(37)

We also know that the step from the current ** R** to the optimal

$$\mathrm{\Delta}\overline{\mathit{R}}=-{\mathbf{F}}^{-1}\nabla \mathrm{\Delta}{\mathbf{G}}_{\mathit{RR}}$$

(38)

This gives

$$\mathrm{\Delta}\overline{\mathit{R}}={\langle \mathit{R}-\overline{\mathit{R}}\rangle}_{{\mathit{U}}_{1}}$$

(39)

The LRA approximation of Eq 35 works the best when *K*^{(1)} is small since this is when we get the largest contribution to Δ*G _{RR}*. Since Eq 39 is not exact (due to the use of the LRA approach), we found it useful to scale it and to use

$$\mathrm{\Delta}\overline{\mathit{R}}=\alpha {\langle \mathit{R}-\overline{\mathit{R}}\rangle}_{{\mathit{U}}_{1}}$$

(40)

where *α* is a scaling factor.

All calculations were performed using the MOLARIS software package, and the ENZYMIX force field.^{10}^{,}^{52}

The constraint release free energy terms are evaluated by a free energy perturbation (FEP)^{53}^{,}^{54} (also see ref^{55}). The FEP method evaluates the free energy associated with the change of the potential surface from *U*^{(1)} to *U*^{(2)} by gradually changing the potential surface using the relationship

$${U}_{m}{\lambda}_{m}={U}^{(1)}(1-{\lambda}_{m})+{U}^{(2)}{\lambda}_{m}$$

where
${U}^{(1)}={U}_{\mathit{syst}}+({K}^{(1)}/2){\displaystyle \sum _{j}}{R}_{j}-{\overline{R}}_{j}{)}^{2}$, *U*^{(2)} = *U _{syst}* and

$$exp\{-\delta G({\lambda}_{m}\to {\lambda}_{{m}^{\prime}})\beta \}={\langle exp\{-({U}_{{m}^{\prime}}-{U}_{m})\beta \}\rangle}_{m}$$

where * _{m}* indicates that the given average is evaluated by propagating trajectories over

$$\mathrm{\Delta}G({U}^{(1)}\to {U}^{(2)})=\sum _{m=0}^{n-1}\delta G({\lambda}_{m}\to {\lambda}_{m+1})$$

The RR calculations require imposing strong positional constraints on the Cartesian coordinates of the solvent molecules in the presence of a solute (either in its charged or uncharged state) whose position is fixed by applying a moderate positional constraint *K*=10 kcal mol^{−1}Å^{−2} and to subsequently evaluate the free energies associated with releasing these restraints using the following equation:

$$-T\mathrm{\Delta}{S}_{\mathit{conf}}=min(\mathrm{\Delta}{G}_{RR}^{\mathit{chg}})-min(\mathrm{\Delta}{G}_{RR}^{\mathit{unchg}})$$

(41)

$\mathrm{\Delta}{G}_{RR}^{\mathit{chg}}$ and
$\mathrm{\Delta}{G}_{RR}^{\mathit{unchg}}$ denotes RR free energies of the systems with a fixed, charged and uncharged solute, respectively. All RR free energies contain a residual contribution from the enthalpy of the system. However, this contribution approaches zero for restraint coordinates that give the lowest RR energy.^{45}^{,}^{47}

To calculate the optimized restraint coordinates (**) that would give the optimal Δ***G* for entropy calculations, the given system with the solute in its charged state, was subjected to MD simulations using a constraint force of 5.0 kcal mol^{−1}Å^{−2} and a scaling factor (*α*) of 0.1 while averaging ** R** −

An initial 30 ps molecular dynamics simulation was performed on the solvent molecules around the fixed charged and uncharged solute, with snapshots taken every 5 ps. MD simulations, 40 ps in length, were performed on each snapshot (with coordinates collected every 0.01 ps) to calculate the QH entropy for each snapshot, and this procedure was carried out for *K*= 3.0 kcalmol^{−1} in magnitude. This allows the calculation of the absolute entropy for this restraint. Subsequently, the free energy for the release of these restraints was evaluated using the RR approach employing Eq 32.

The RR-FEP calculations were performed with an 18Å simulation sphere of explicit water molecules subject to the surface constraint all-atom solvent (SCAAS) boundary conditions.^{56} The RR simulations, consisting of 42 windows each with a simulation time of 40ps, involved the release of the position restraints in 5 FEP stages, changing *K*_{1} from 3.0 to 1.0, from 1 to 0.3, from 0.3 to 0.03, from 0.03 to 0.003 and finally, from 0.003 to 0.0003, where all the values of *K*_{1} are given in kcal mol^{−1} Å^{−2}. These simulations were done at 300K with 1fs time steps.

Simultaneously, the LD model^{33} implemented in the program ChemSol 2.0^{57} was also used to evaluate the solvation entropy. This method evaluates the solvation entropy by estimating the restriction of the solvent motion due to the field from the solute.

As stated above, we tried to explore in this work the ability of the RR approach to reproduce solvation entropies. The most obvious benchmark is provided by the solvation entropies of ions and this challenge has been addressed here by evaluating the solvation entropies of Na^{+}, K^{+} and Ca^{2+} ions. The corresponding RR calculations were performed on the inner 29 water molecules that are closest to the ion, while the rest of the 18Å SCAAS system to move with a very small position constraint on the oxygen atoms in addition to the SCAAS polarization and distance restraints. In short, the entropic effect associated with the effect of the solute on the inner 29 water molecules was explored while representing rest of the system by the SCAAS model. In calculating the solvation entropy, one should also consider the hydrophobic energy of forming the cavity for the uncharged ion. This hydrophobic contribution is small for the ions considered here and were added to the final calculated results. This was done using the corresponding values obtained from ChemSol. In future, we hope to obtain such contributions from our RR calculations. The results of these calculations are summarized in Tables 1–3. Remarkably, we reproduce the observed trend in an almost quantitative way. The tables also give the quantitative LD results but as mentioned before, the LD treatment is empirical in nature and cannot be used within the protein active site. It should also be mentioned that the recent works of Aqvist and coworkers^{28} reproduced quantitative results for the same system studied here. However, this was done by evaluating the van’t Hoff approach that requires a major investment in terms of computer time.

The fact that the computational methodology employed in this work reproduces the observed solvation entropy indicates that the physics of this effect is captured by this approach and that the effect is associated with the change in the effective space available for the solvent molecules upon change in the solute charges.

The overall trend is clearly discussed in section II, but here we have a well-defined simulation approach that reproduces the entropic effect quantitatively. Furthermore, in order to clarify the nature of the challenge overcome by the RR calculations, we also report in Table 4 the QH results obtained with no restraint (*K*= 0 kcal mol^{−1} Å^{−2}) for Na^{+} and Ca^{2+} ions. Although QH entropic estimation for Na^{+} ion seems fairly accurate, it provides completely erroneous results for Ca^{2+} ion. Thus, we may say that the QH is simply unable to capture the solvation entropies. This is probably one of the best illustrations of the fundamental problems with the widely used QH and related approaches.

The present work focuses on the calculations of the entropic effects, despite the fact that most questions about biophysics functions are related to the corresponding free energy surfaces.^{58}^{,}^{59} Additionally, we are well aware of the fact that free energy calculations converge much faster than the entropic calculations. However, in many cases complete understanding of the entropy-enthalpy compensation is useful for the deeper understanding of the underlying physics of the system and in particular, the origin of the temperature dependence of various processes. Moreover, the evaluation of the entropic effects is one of the utmost challenges of computational biophysics and in fact, the failure to reproduce entropic effects is frequently used by the experimental community in trying to diminish the important free energy calculations (also see discussion in ref^{60}).

In any case, the accurate evaluation of the entropic contributions to the free energy changes in condensed phases is an extremely challenging task. Popular approaches such as the QH and related treatments have been widely used to extract the configurational entropy from molecular dynamics simulations.^{48}^{,}^{61}^{,}^{62} Unfortunately, although QH approximation functions quite well for simple systems with single energy wells, this approach provides a poor approximation of the configurational entropies of systems with multiple occupied energy wells, as it merges multiple narrow energy wells into one broad energy well.^{16}^{,}^{17} In fact, the validity of such approaches in systems with very shallow potentials is questionable. Therefore, the use of the QH approach in the studies of binding energies or related problems is problematic. Several approaches have been proposed for obtaining quantitative entropic contributions, but most of these approaches have not been able to provide quantitative results in general cases, except by brute force studies of the free energy change with temperature that requires enormous simulation time for convergence.^{28} In our view and experience, one of the most promising strategies is our RR approach that has been shown to give quantitative results for activation entropies.^{10}^{,}^{12}^{,}^{13}^{,}^{47}^{,}^{50} Nevertheless, this approach has not been validated in the studies of binding entropies.

One of the most serious problems in calculations of binding entropy is the evaluation of the solvation entropies in water and in the protein active site. This contribution has been generally neglected in key studies^{29}^{–}^{32} although it has been evaluated or included in calculations of activation free energies.^{13}^{,}^{15}^{,}^{47} The problem of activation entropies is also a problem of general interest both for studies of molecules in solutions and in proteins.

When one deals with molecules in solution, it is possible to use the LD models to obtain a relatively reliable estimate of the solvation entropies. Unfortunately, this approach has not been generalized to protein active sites where one must use a microscopic approach. The studies of ions in solution should provide the best benchmark for the development and validation of the corresponding approach since the corresponding values are known experimentally.

The molecular simulations of the entropy contribution in the areas of protein-ligand binding and enzyme catalysis, is an extremely challenging task and the subject of intense research. The analyses of the role of solvent entropies in ligand binding and enzyme catalysis are still rare, although knowledge of their values would be of great utility considering that the binding processes can be relatively dominated by solvent entropy effects. Therefore, there is a strong motivation to obtain the solvation entropic contributions that can lead to a better understanding of the process under consideration. Although the solvent entropies can be extracted from the temperature dependence of the explicit solvent FEP simulations, these calculations are difficult to converge and their convergence tends to be considerably slower than that of the free energy itself. The present work focuses on the development of efficient microscopic model for the studies of the solvation entropies in aqueous solution. Our preliminary studies clarified that the solvation entropies are not associated with the solvent vibrations but rather with the availability of larger configurational space upon transition from the charged to the uncharged form of the solute. This effect is easily observed in the cases of divalent ions but in the case of monovalent ions it is restricted to the first solvation shell. With this in mind, we demonstrated that the solvent entropies can be successfully and accurately evaluated by our RR approach. This has been illustrated for the simple but relevant test case of monovalent and divalent metal ions. It is clear from the results that solvation entropies depend strongly on the magnitude of the charge of the solvated ion. This is obviously related to the corresponding restriction in the solvent fluctuations.

Although convergence problems may exist for larger and more complex simulation systems, a good agreement is obtained between the calculated and the experimental data for solvation entropies of monovalent and divalent metal ions. Furthermore, we also developed a systematic approach for the minimization of *ΔG _{RR}* through the Newton-Rhapson minimization of the

The success of the present approach in the evaluation of the electrostatic components of the solvation entropies indicates that the RR method will probably be useful in the studies of hydrophobic effects where the physics (of restricted orientation) is smaller. Efforts are underway in our laboratory to extend the utility of this approach in the studies of the contribution of the solvation entropies to the binding entropies, hoping to help advance the understanding of this important effect. The corresponding strategy will be similar to the one used in the present study and will require releasing the restraints on the protein residues and the water molecules within a given cut-off distance. However, this approach will require careful refinement since releasing the restraints on all of the protein residues may lead to computational instability while releasing the restraints on a few residues may not capture the full solvation entropy contribution. In future, it should be feasible to estimate the magnitude of the solvent entropic contributions involved in ligand binding to proteins. Generally, the binding affinity of a compound can be improved by generating a favorable binding enthalpy, favorable solvation entropy, or by minimizing the unfavorable conformational entropy. In some cases, it is much easier to optimize the entropy. For example, medicinal chemists have attempted to conformationally constrain the substrate thus, restricting the entropy loss upon binding. Further improvement may be obtained by trying to optimize the solvation entropies by choosing specialized functional groups. However, the main advancement in our approach is in obtaining more reliable predictions of the entropic contributions to ligand binding. Our RR approach should therefore, be viewed as a step towards the generation of algorithms capable of thermodynamic dissection particularly, those that can predict the entropic consequences of introducing different chemical functionalities in the molecular scaffold generated either *de novo* (during lead generation) or while being optimized (during lead optimization)for improving their binding affinity towards a target protein. Thus, this study may have implications towards our understanding of protein-ligand interactions, for theoretical modeling of molecular recognition, and for rational drug design.

This work was supported by National Institutes of Health Grant R01 GM40283. All computational work was supported by the University of Southern California High Performance Computing and Communication Center (HPCC). The authors are grateful to Drs. Z. T. Chu and Pankaz Sharmafor their technical assistance and insightful suggestions throughout this work.

1. Meirovitch H. Current Opinion in Structure Biology. 2007;17:181. [PubMed]

2. Amzel LM. Proteins: Structure, Function, and Genetics. 1997;28:144. [PubMed]

3. Brady GP, Sharp KA. Current Opinion in Structural Biology. 1997;7:215. [PubMed]

4. Privalov PL, Makhatadze GI. J Mol Biol. 1993;232:660. [PubMed]

5. Creamer TP. Methods in Molecular Biology (Totowa, NJ, United States) 2001;168:117. [PubMed]

6. Li DW, Khanlarzadeh M, Wang J, Huo S, Brüschweiler R. J Phys Chem B. 2007;111:13807. [PubMed]

7. Finkelstein A, Janin J. Protein Eng. 1989;3:1. [PubMed]

8. Murray C, Verdonk M. J Comput Aided Mol Des. 2002;16:741. [PubMed]

9. Gilson MK, Given JA, Bush BL, McCammon JA. Biophys J. 1997;72:1047–1069. [PubMed]

10. Sham YY, Chu ZT, Tao H, Warshel A. Proteins: Struct Funct Genet. 2000;39:393. [PubMed]

11. Chang CE, Chen W, Gilson MK. Proc Natl Acad Sci U S A. 2007;104:1534. [PubMed]

12. Villà J, Štrajbl M, Glennon TM, Sham YY, Chu ZT, Warshel A. Proc Natl Acad Sci USA. 2000;97:11899. [PubMed]

13. Sharma PK, Xiang Y, Kato M, Warshel A. Biochemistry. 2005;44:11307. [PubMed]

14. Warshel A. Annu Rev Biophys Biomol Struct. 2003;32:425. [PubMed]

15. Bjelic S, Brandsdal B, Aqvist J. Biochemistry. 2008;47:10049. [PubMed]

16. Chang CE, Chen W, Gilson MK. J Chem Theory Comput. 2005;1:1017. [PubMed]

17. Carlsson J, Aqvist J. J Phys Chem B. 2005;109:6448. [PubMed]

18. Meirovitch H. Calculation of the free energy and entropy of macromolecular systems by computer simulation. Vol. 12 Wiley-VCH; 1998.

19. Hnizdo V, Fedorowicz A, Singh H, Demchuk E. J Comput Chem. 2003;24:1172. [PubMed]

20. Cheluvaraja S, Meirovitch H. Proc Natl Acad Sci USA. 2004;101:9241. [PubMed]

21. Darian E, Hnizdo V, Fedorowicz A, Singh H, Demchuk E. J Comput Chem. 2005;26:651. [PubMed]

22. Cheluvaraja S, Meirovitch H. J Chem Phys. 2006;125:24905. [PubMed]

23. Wang J, Bruschweiler R. J Chem Theory Comput. 2006;2:18.

24. Hnizdo V, Darian E, Fedorowicz A, Demchuk E, Li S, Singh H. J Comput Chem. 2007;28:655. [PubMed]

25. Hnizdo V, Tan J, Killian BJ, Gilson MK. J Comput Chem. 2008;29:1605. [PMC free article] [PubMed]

26. Cheluvaraja S, Mihailescu M, Meirovitch H. J Phys Chem B. 2008;112:9512. [PMC free article] [PubMed]

27. Levy RM, Karplus M, Kushick J, Perahia D. Macromolecules. 1984;17:1370.

28. Carlsson J, Aqvist J. Phys Chem Chem Phys. 2006;8:5385. [PubMed]

29. Lee MS, Olsson MA. Biophys J. 2006;90:864. [PubMed]

30. Chen W, Chang CE, Gilson MK. Biophys J. 2004;87:3035. [PubMed]

31. Ruvinsky AM. Journal of Computer-Aided Molecular Design. 2007;21:361. [PubMed]

32. Chang CE, Chen W, Gilson MK. Proc Natl Acad Sci U S A. 2007;104:1534. [PubMed]

33. Florián J, Warshel A. J Phys Chem B. 1999;103:10282.

34. Warshel A, Hwang JK. J Chem Phys. 1986;84:4938–4957.

35. Warshel A, Lifson S. Journal of Chemical Physics. 1970;53:582.

36. King G, Warshel J Chem Phys. 1990;93:8682.

37. Frank HS. J Chem Phys. 1945;13:478.

38. Frank HS, Evans MW. J Chem Phys. 1945;13:507.

39. Langevin P. Annales de Chimie Physique. 1905;8:70.

40. Warshel A, Levitt M. J Mol Biol. 1976;103:227. [PubMed]

41. Papazyan A, Warshel A. J Phys Chem B. 1998;102:5248.

42. Rashin AA, Bukatin MA. Biophys Chem. 1994;51:167. [PubMed]

43. Rashin AA, Bukatin MA. J Phys Chem. 1994;98:386.

44. Luzhkov V, Warshel A. J Comp Chem. 1992;13:199.

45. Hermans J, Wang L. J Am Chem Soc. 1997;119:2707.

46. Warshel A. Office of Naval Research, Award Number N00014-91-J-1318. 1991. Computer simulation of chemical reactions in synthetic model compounds and genetically engineered active sites.

47. Štrajbl M, Sham YY, Villà J, Chu ZT, Warshel A. J Phys Chem B. 2000;104:4578.

48. Karplus M, Kushick JN. Macromolecules. 1981;14:325.

49. DiNola A, Berendsen H, Edholm O. Macromolecules. 1984;17:2044.

50. Rosta E, Kamerlin S, Warshel A. Biochemistry. 2008;47:3725. [PubMed]

51. Warshel A. Computer Modeling of Chemical Reactions in Enzymes and Solutions. John Wiley & Sons; New York: 1991.

52. Lee FS, Chu ZT, Warshel A. J Comp Chem. 1993;14:161.

53. Zwanzig RW. J Chem Phys. 1954;22:1420.

54. Valleau JP, Torrie GM. A guide to Monte Carlo for statistical mechanics. 2. Byways. In: Berne BJ, editor. Modern Theoretical Chemistry. Plenum Press; New York: 1977. pp. 5–169.

55. Warshel A. In Encyclopedia of Molecular Biology. John Wiley & Sons; 1999.

56. Warshel A, King G. Chem Phys Lett. 1985;121:124.

57. Florian J, Warshel A. ChemSol. 2.1. University of Southern California; Los Angeles: 1999.

58. Olsson MH, Parson WW, Warshel A. Chem Rev. 2006;106:1737. [PubMed]

59. Warshel A, Parson WW. Quart Rev Biophys. 2001;34:563. [PubMed]

60. Olsson MH, Siegbahn PE, Warshel A. Journal of the American Chemical Society. 2004;126:2820. [PubMed]

61. Andricioaei I, Dinner AR, Karplus M. J Chem Phys. 2003;118:1074.

62. Schlitter J. Chem Phys Lett. 1993;215:617.

63. Marcus Y. Ion Properties. Marcel Dekker, Inc.; New York: 1997.

64. Noyes RM. J Am Chem Soc. 1962;84:513.

PubMed Central Canada is a service of the Canadian Institutes of Health Research (CIHR) working in partnership with the National Research Council's national science library in cooperation with the National Center for Biotechnology Information at the U.S. National Library of Medicine(NCBI/NLM). It includes content provided to the PubMed Central International archive by participating publishers. |