Article pubs.acs.org/JPCB
Challenge of Representing Entropy at Different Levels of Resolution in Molecular Simulation Wei Huang and Wilfred F. van Gunsteren* Laboratory of Physical Chemistry, Swiss Federal Institute of Technology, ETH, 8093 Zürich, Switzerland S Supporting Information *
ABSTRACT: The role of entropic contributions in processes involving biomolecules is illustrated using the process of vaporization or condensation of the solvents water and methanol and the process of polypeptide folding in solution using molecular models at different levels of resolution: subatomic, atomic, supra-atomic, and supramolecular. For the folding process, a β-hexapeptide that adopts, as inferred from NMR experiments, both a right-handed 2.710/12-helical fold and a left-handed 314-helical fold in methanol, is used to illustrate the challenge of modeling thermodynamically driven processes at different levels of resolution.
■
INTRODUCTION Four decades ago, reports on the first molecular dynamics (MD) simulations of liquid water at ambient conditions of temperature and pressure had appeared in the literature.1,2 They showed a highly mobile picture of water molecules forming and breaking hydrogen bonds and emphasized the need to include motion or configurational variability when aiming at an accurate description of the properties of compounds in the liquid phase. In other words, when proposing a molecular model not only energetic contributions but also entropic contributions to the free energy of a liquid are to be taken into account. For example, the experimental value of the excess free energy of liquid water at room temperature and pressure is with 24 kJ mol−1 about half the size of its heat of vaporization of 44 kJ mol−1. Thus, TΔS, where T is the temperature and ΔS the difference in entropy between gas and liquid phase, is about 20 kJ mol−1 (see Table 1). This is why molecular models for liquid water and other biomolecular solvents are to be developed and their parameters calibrated based on measured thermodynamic data for such compounds and configurational ensembles generated by MD or Monte Carlo (MC) simulation.3 Entropic contributions are also non-negligible when considering biomolecular processes such as polypeptide folding4 or protein−ligand association.5 The enthalpy and entropy of solvation of biomolecules may vary depending on the cosolvents present in an aqueous solution, thereby leading to different, enthalpy- versus entropy-dominated, mechanisms of partitioning6 or association.7 In order to properly model such processes, the molecular models involved should not only © XXXX American Chemical Society
represent energetic but also represent entropic effects, the distribution of molecular configurations. This was recognized by Bill Jorgensen and Herman Berendsen during the late seventies and early eighties, when they derived their simple models for liquid water, SPC,8 and TIP3P,9 and applied multistep perturbation and thermodynamic integration techniques to calculate free-energy differences.10,11 At the end of the seventies, Bill Jorgensen and one of us (van Gunsteren) published their first papers on molecular simulation12,13 and have since then contributed to the field of biomolecular simulation through force field or model development based on a thermodynamic basis. Thirty-six years later, the field of biomolecular simulation is still facing difficult challenges, one of which is the accurate representation of entropic effects by supra-atomic or supramolecular models used in simulations of molecular processes. It is therefore a pleasure to dedicate this article to Bill Jorgensen on the occasion of his 65th birthday in honor of his contributions to the field. The free energy, F, of a molecular system can be expressed as the difference between the energy, U, and the temperature, T, times the entropy, S, of a system, (1)
F = U − TS
Since the free energy, F, is an integral of an integrand exp(−U(p⃗ N ,r ⃗ N )/kBT) over the whole phase space, its value Special Issue: William L. Jorgensen Festschrift Received: May 22, 2014 Revised: August 18, 2014
A
dx.doi.org/10.1021/jp505045m | J. Phys. Chem. B XXXX, XXX, XXX−XXX
The Journal of Physical Chemistry B
Article
During the past decades, the use of so-called supra-atomic or supramolecular models for biomolecules20−22 such as proteins or lipids in membranes has become popular, this because of the enhanced computational efficiency of such coarse-grained (CG) molecular models compared to the atomic-level of resolution, fine-grained (FG) models. However, by eliminating atomic or molecular degrees of freedom in the process of coarse-graining, the energy and entropy of the system are generally reduced. Yet, their relative sizes should be preserved from the FG to the CG model or level of resolution. For example, comparing two FG atomic-resolution models for water and methanol with two CG supramolecular resolution models for these compounds, the relative contributions of the energy and entropy to the excess free energy are similar to the experimental values, both at the FG and the CG levels of modeling (see Table 1). The CG models for water23 and methanol24 were developed by calibration to thermodynamic quantities with this property in mind. However, most CG models for proteins and lipids and their solvents currently in use have not been calibrated or evaluated against thermodynamic energy and entropy data. This leaves the question of their suitability for the prediction of thermodynamic properties and thermodynamically driven processes wide open. As the data in Table 1 illustrate, it is possible to develop supramolecular models for simple liquids that approximately preserve the relative contributions of entropy and energy to the free energy at a particular thermodynamic state point. The next challenge is to investigate whether this is also possible for the process of secondary-structure formation of polypeptides using an atomic level of resolution FG model for the polypeptide solute and a supramolecular CG model for the solvent. In the present study, the entropy difference for a β-hexapeptide in methanol between its two alternative, helically folded conformations and its unfolded conformations is evaluated in three different ways from MD simulations at different temperatures using either an FG or a CG model for the methanol solvent. The results of the different models and ways to compute entropy differences are compared and illustrate the challenge of representing entropic contributions to polypeptide folding.
Table 1. Excess Free Energy (ΔFexc) and Heat of Vaporization (ΔHvap) for Liquid Water and Liquid Methanol at Room Temperature and Pressure as Obtained from Experiment and from MD Simulations Based on FineGrained (FG) United-Atom Level of Resolution Models and Coarse-Grained (CG) Supra-Molecular Level of Resolution Models name of model experiment polarizable FG model nonpolarizable FG model polarizable CG model experiment polarizable FG model nonpolarizable FG model polarizable CG model
COS/ G228 SPC8 CGW23
COS/ M29 B327 CGM24
ΔFexc (kJ mol−1)
ΔHvap (kJ mol−1)
ΔHvap − ΔFexc (kJ mol−1)
Water 24.041 21.825
44.042 43.725
20.0 21.9
23.6
43.9
20.3
11.0
25.9
14.9
38.1a,b 37.8a,b
20.3 23.1
16.7
37.9
21.2
7.3
18.0
10.7
Methanol 17.829 14.729
a,b
Values were calculated using eq 3, the potential energies are from refs 43 and 29, respectively.
can only be calculated for very simple systems, such as an ideal gas or a set of independent harmonic oscillators. Here ⎯p , → ⎯p , ..., ⎯→ ⎯r , → ⎯ ⎯→ p⃗ N = (→ pN ) and r ⃗ N = (→ 1 r2, ..., rN ) are the 1 2 momenta and the positions of the N particles of the system, and kB is the Boltzmann constant. Fortunately, one is generally interested in free-energy differences that drive particular processes, and these can be obtained from molecular simulation trajectories or ensembles.14 Since the energy U (p⃗ N , r ⃗ N ) = K (p⃗ N ) + V ( r ⃗ N )
(2)
of a system is a function of the phase-space variables, it can be evaluated as an average over a trajectory or the phasespace configurations of the system generated by MD or MC simulation. The kinetic energy (K) is usually only a function of p⃗ N and the potential energy V is usually only a function of r ⃗ N . Like free-energy differences or relative free energies, ΔF between two systems, two states of a system or two parts of the configurational space of a system, the corresponding entropy differences can be calculated using techniques such as thermodynamic integration15 or perturbation16 and variations thereof. Yet, obtaining converged values for entropy differences requires much longer simulations than for free-energy differences.17 In the present study, we use three different ways to compute the entropic contribution to secondary-structure formation of a β-hexapeptide in methanol that adopts, according to NMR experiments,18 in equilibrium a left-handed 314-helix and a right-handed 2.710/12-helix (e.g., about 2% 314helix and 16% 2.710/12-helix at 340 K and 1 atm using the GROMOS 43A1 force field).19 Since the equilibrium between the two different helical folds and the unfolded conformations is sensitive to the temperature and to the force field or interaction model used in the simulation, it offers the opportunity to investigate energy−entropy compensation or cooperation in regard to secondary-structure stability.
■
METHOD Theoretical Background. For the pure liquids, the heat of vaporization, ΔHvap, was obtained as the ensemble average of the potential energy, EV (or V), from MD simulations at constant temperature and pressure, ΔH vap = ( −⟨EV ⟩ + pV )/Nmol
(3)
where Nmol is the number of molecules in the system, p is the pressure, and V the volume of the system. The excess free energy, ΔFexc, is obtained through thermodynamic integration by switching off the nonbonded interaction V as a function of a coupling parameter λ in MD simulations at constant temperature and volume, −1 ΔFexc = Nmol
∫0
1
⟨∂V /∂λ⟩λ dλ
(4)
where the λ dependence of the molecular interaction function is given in ref 25. In order to obtain an estimate of the free enthalpy, enthalpy, and entropy of folding of a polypeptide in solution at equilibrium, the configurational space of the peptide is B
dx.doi.org/10.1021/jp505045m | J. Phys. Chem. B XXXX, XXX, XXX−XXX
The Journal of Physical Chemistry B
Article
Figure 1. (a) Structural formula of the β-hexapeptide studied. In the simulations, both end groups were protonated, in accordance with previous studies and experimental data. (b) The right-handed 2.710/12-helical model structure. (c) The left-handed 314-helical model structure.
separated into folded configurations F (e.g., L-helical and Rhelical ones) and unfolded configurations U, using a geometric criterion.26 The free enthalpy of folding ΔGTFU at temperature T can then be obtained from T ΔG FU
=
G FT
−
G UT
= −RT
ln(PFT /PUT )
Cp =
(5)
T2 T1 ΔSFU − ΔSFU =
∫T
T2−1
−1 1
⟨ΔHFU⟩T dT −1
T2 ΔSFU
−
T1 ΔSFU
=
T2
−
T1 ΔHFU
T1
ST2 − ST1 =
∫T
1
T
dT
T
dT
(11)
(12)
and this difference at temperature T between two not too different temperatures T1 and T2 (T1