Letter pubs.acs.org/JPCL
Surface Polarity and Nanoscale Solvation Allan D. Friesen and Dmitry V. Matyushov* Center for Biological Physics, Arizona State University, P.O. Box 871504, Tempe, Arizona 85287-1504, United States S Supporting Information *
ABSTRACT: Linear solvation theories are well established to describe electrostatic hydration of small solutes when the hydration free energy is dominated by the electrostatic free energy of the solute multipole. In contrast, hydration of nanometer solutes is driven by surface hydration. We address the question of whether the linearresponse thermodynamics established for small multipolar solutes applies to surface hydration. To this end, molecular dynamics simulations are carried out on a model C180 solute that carries no global multipole, but the surface of which is decorated with radially pointing dipoles. Linear response is dramatically violated in this case. Further, two crossovers in the solvation thermodynamics are discovered as the surface polarity is increased. Both transformations produce strongly nonlinear solvation response. The second, more collective, crossover leads to a dramatic slowing down of the interfacial dynamics, reaching the time-scales of nanoseconds. Our picture offers the possibility of flipping water domains at interfaces of nanoparticles and biomolecules. SECTION: Liquids; Chemical and Dynamical Processes in Solution
F
are a good example. While they often carry a net charge and a large overall dipole moment (200−300 D),13 these multipoles tell little about protein solvation. Charge at the surface of a folded globular protein is fairly uniformly distributed following charge complementarity, i.e., the oppositely charged residues are neighbors.14,15 The protein multipole is largely a result of an incomplete compensation of charges within the patchwork of positive and negative sites and polar groups,16 which are locally hydrated and contribute to the overall hydration free energy.4 In other words, local interactions of hydration waters with the surface residues are more important than interaction of waters with the global multipole of the protein. While a quantitative description of the solvation thermodynamics of large solutes requires full account of its charge distribution, simplified models can help to understand whether surface solvation specific to such problems will produce measurable qualitative deviations from the standard picture of multipole solvation. This is the motivation of this study. We study the thermodynamic signatures of hydration of a simple solute with its surface populated with polar groups to address the question of whether the standard approximations developed for hydration of small multipolar solutes apply to problems of this type. We approach this problem by studying a model system based on the rigid structure of C180 fullerene.17 The polarity of its surface is artificially modified by placing dipoles at each surface carbon pointing along the radial direction specified by the vector connecting the center of mass of the fullerene to the carbon atom. The dipoles are created by placing two opposite
avorable hydration of large solutes (nanoparticles, proteins, etc.) is achieved in most cases by polar groups and charges placed at the solute−water interface.1−3 The thermodynamics of solvation is determined in such cases by surface solvation4,5 and by the relevant alteration of the surface tension, rather than by the set of global multipoles of the solute. The solute’s electric multipoles provide the main electrostatic input to the classical theories of solvation going back to Born (ions) and Onsager (dipoles). In these models, the solvation free energy is obtained by integrating the electrostatic energy density, produced by the electric field of the solute multipoles and scaled with an appropriate solvent response function, over the volume occupied by the solvent. These models predict the free energy of solvation to scale quadratically with the global solute multipole. The proportionality coefficient, or solvation susceptibility, has been the focus of many models and computations distinguishing themselves by the extent of molecular detail.6 The quadratic scaling with the solute multipole leads to a number of relations between experimentally accessible thermodynamic solvation functions7,8 and/or between cumulants of the solute−solvent interaction energy accessible by spectroscopy.9 Together, these relations establish the frameworks of the linear response theory of solvation phenomena.8,10 Note that ionic hydration also involves a term that is linear in the charge, due to water’s molecular asymmetry (molecular quadrupole is large). This component accounts for a constant electrostatic potential, independent of the solute charge, within the solute.8,11,12 While the framework of classical solvation theories and of the linear response approximation are well established, one wonders to what extent they apply to solvation of large solutes for which the surface polarity, instead of a global multipole, might be the dominant property to consider. Hydrated proteins © 2012 American Chemical Society
Received: October 17, 2012 Accepted: November 28, 2012 Published: November 28, 2012 3685
dx.doi.org/10.1021/jz301672e | J. Phys. Chem. Lett. 2012, 3, 3685−3689
The Journal of Physical Chemistry Letters
Letter
ε is given by one universal Gaussian function for all values of the surface dipole
charges, separated by 0.2 Å, at the center of mass of each carbon atom of C180 fullerene. The resulting solute carries a potentially highly polar surface, while having a negligible overall multipole. The dipole moment of the solute is essentially zero. Further, when the dipole moment of 1 D is placed on each carbon atom, the quadrupole moment of fullerene Q = [(2/3) Q:Q]1/2, determined by rotationally invariant contraction of the quadrupole matrix Q,18 is 0.65 D × Å, which is below the quadrupole moment of water ≃3 D × Å. The solute as a whole carries no traceable multipole moment, and all solvation results from surface polarity, associated here with the density of surface dipoles. The surface polarity affects the surface tension γ. This connection can be established by modifying the well-known Lippmann equation of colloid science.3 The change of the free energy of the system at constant temperature and volume is affected by the variation of the surface area dS and the electrostatic energy of the surface carbon dipoles mC dF = γ dS − NCECdmC
P(x) ∝ exp[−(δx)2 /2]
NCe0s)/σε, σ2ε
(6)
where δx = (ε − = ⟨(δε) ⟩. We indeed find that the distributions obtained from simulation trajectories at various mC can all be superimposed on this function. Deviations, however, clearly exist in the tails of the distribution. Those are often neglected, particularly when simulations suffer from insufficient sampling. They, however, hide some significant thermodynamic features, which become more prominent when higher cumulants of ε and thermodynamic potentials other than the solvation Gibbs energy are studied.21 We find here significant deviations from the linear response thermodynamics when surface solvation dominates. We leave the description of the simulation protocol to the Supporting Information (SI) and proceed to the results. Figure 1 shows the thermodynamic potentials calculated from eqs 2, 4, and 5 versus the average energy of the solute−
(1)
2
Here, NC = 180 is the number of dipoles and EC is the average field acting on the dipole. From this equation, (∂γ/∂mC)S = −(∂NCEC/∂S)m, and the variation of the surface tension caused by increasing the carbon dipole from zero to a given magnitude is Δγ =
∂ NC ∂S
∫0
mC
e0s(m)(dm /m) =
∂ μNC ∂S
(2)
where e0s = ⟨ε⟩/NC and μ = G/NC are the average solute− solvent electrostatic energy and the electrostatic chemical potential per individual dipole, respectively. Here, ε is the instantaneous solute−solvent interaction energy and G is the electrostatic hydration Gibbs energy of the solute. Further, the integral in eq 2 is the standard equation for thermodynamic integration18 connecting the solvation chemical potential to the average solute−solvent interaction energy e0s(m) (per carbon atom) as a function of the carbon dipole. Since the interaction potential scales linearly with the carbon dipole, the first order expansion of the distribution function in the solute−solvent potential,18 that is, the linear response approximation, leads to the identity7,19−21 2μ = e0s = −β⟨(δε)2 ⟩/NC ∝ (mC)2
Figure 1. Thermodynamic potentials of electrostatic solvation of C180 versus negative of the average solute−solvent interaction energy per carbon dipole, −e0s. The points are MD results, and the dashed line shows the expectation of the linear response approximation, −μ = −0.5e0s, eq 3. The dotted line of the unity slope is drawn to guide the eye.
solvent interaction e0s. The dashed line in the plot indicates the prediction of the linear response approximation, μ = 0.5e0s (eq 3). This approximation clearly does not describe the results in the whole range of surface dipoles studied here, as we discuss in more detail below. However, all thermodynamic potentials are remarkably linear in the interaction energy, with varying slopes. We find that the results are fitted by μ = 0.29e0s, ess = −0.28e0s, and e = 0.72e0s. Several features seem to be particularly noteworthy. First, the negative of chemical potential of solvation is nearly equal to the energy component ess: −μ ≃ ess. Second, we observe a substantial negative entropy of solvation per carbon, −Ts > −μ; the entropy is calculated from the standard equation μ = e − Ts. The dashed line in Figure 1 also specifies the linear response prediction for the solute−solvent entropy Ts0s = Ts − Tsss = (1/2)e0s, where the second equality is valid only in the case of linear response. The second component, the entropy associated with the change of solvent−solvent interactions, is required by thermodynamics to be equal to the corresponding energy change,22 Tsss = ess. The actual entropy Ts0s calculated from the simulations is significantly more negative, indicating that linear response underestimates the reduction of the configuration space of the shell waters imposed by the solute−solvent electrostatic interaction.
(3)
where β = 1/(kBT) is the inverse temperature. It is clear, however, that, when linear response does not hold, the chemical potential is determined by an infinite expansion in the cumulants of ε. The same statement applies to the solvation energy per carbon of C180
e = e0s + ess
(4)
which is a sum of the first solvation moment, e0s, and an infinite expansion in the cumulants of ε contained in the alteration of the solvent−solvent interaction energy produced by the solute19−22 ess = −(β /NC)
∫0
mC
⟨δεδUss⟩(dm /m)
(5)
Here, δUss is the fluctuation of the total solvent−solvent interaction energy. The linear response approximation is equivalent to approximating the statistics of the stochastic variable ε by a Gaussian distribution.7,21,23 This implies that the distribution of 3686
dx.doi.org/10.1021/jz301672e | J. Phys. Chem. Lett. 2012, 3, 3685−3689
The Journal of Physical Chemistry Letters
Letter
examine the first-order orientational parameter ⟨pI1⟩ of the hydration shell, which is equal to the average cosine of the surface water dipoles projected on the surface normal
The linear response approximation predicts quadratic dependence of the thermodynamic potentials on the magnitude of the surface dipole mC (eq 3). This prediction provides an additional opportunity to test the consistency of the linear response picture. This analysis is shown in Figure 2, where we
⟨p1I ⟩ = (N I)−1 ∑ ⟨cos θi⟩ i
(7)
I
Here, N is the number of waters with their oxygens within 9.5 Å from the solute center of mass, and subscript “1” refers to S = 1 (first order) in Legendre polynomials PS(x) defining the order parameter. This parameter is indicated by the solid squares in the inset of Figure 3a. The value of ⟨pI1⟩ is small and positive for small
Figure 2. Panel (a) shows the electrostatic chemical potential −μ/mC (diamonds), hydration energy −e/mC (circles), and the average solute−solvent energy −e0s/2mC (squares) plotted versus the carbon dipole moment mC. Linear response (dashed line labeled “L.R.”) requires all energies to scale linearly with mC, with equal slopes for −μ/mC and −e0s/2mC (eq 3). The inset in panel (a) magnifies the initial portion of the plot. Panel (b) shows the log−log plot of −μ/mC (diamonds) and ess/mC (crosses) versus mC. The slopes indicate how −μ/mC scales with mC. A unity slope corresponds to linear response. Panel (c) shows χG = −β⟨(δε)2⟩/(NCe0s) versus mC.
plot the thermodynamic potentials divided by mC, e.g., −μ/mC. These scaled functions are expected to be linear in mC if linear response holds. The response is approximately linear only up to about mC = 2 D. After that point, all thermodynamic potentials increase in magnitude significantly faster than m2C. The scaling for larger dipoles is approximately μ ∝ (mC)3.5, as is seen from Figure 2b, which shows the log−log plot of −μ/mC versus mC. The equality between the first and second (multiplied with β) cumulants of the solute−solvent interaction energy predicted by the linear response approximation (eq 3) can be tested by looking at the parameter χG = −β⟨(δε)2⟩/(NCe0s) quantifying the deviation from the linear response. This is shown in Figure 2c. As the surface dipoles are increased, χG passes through two distinct peaks, marking crossovers in the structure of the hydration layer. The first crossover was observed in simulations of a single dipole placed near the surface of a spherical solute.24 In that study, an observed maximum in χG marked the appearance of defects involving one surface water molecule reoriented by the field of the solute, with some strain introduced in the nework of interfacial hydrogen bonds. The second peak of χG, which is a new feature found for C180, marks a more dramatic crossover when singlemolecule excitations merge into an interfacial structure, as will be clear from the analysis of the solvation dynamics, below. The thermodynamic potentials presented in Figures 1 and 2 paint a global picture lacking molecular detail. To examine the molecular-level changes associated with these observations, we
Figure 3. Slow exponential relaxation time, τE, of the interaction energy ε(t) (open circles) and cosine of shell dipoles pI1(t) (solid squares) versus the average solute−solvent interaction energy per carbon dipole e0s. The inset in (a) shows the orientational parameter ⟨pI1⟩, equal to the average cosine of surface water dipoles with the interface normal (eq 7). Panel (b) presents the data from the upper panel in the logarithmic scale.
mC. With increasing dipole moment, ⟨pI1⟩ decreases, with two kinks (see also Figure 4 in the SI), before saturating at ⟨pI1⟩ ≈ −0.25. The negative values of ⟨pI1⟩ for large solute dipoles indicate interfacial water dipoles oriented somewhat against the dipoles at the carbon atoms. This seemingly counterintuitive result is easily understood from the arrangement of water molecules around fullerenes. Greater contact area and more favorable packing are achieved when waters occupy the area near the centers of the faces of a fullerene.25 Since the waters at the faces are affected by equatorial dipolar fields, they tend to orient opposite to the solute dipoles. This highlights the effect of the solute topology on the orientational properties of the interfacial waters. The second spike in the plot of χG versus mC has a dynamic signature to it. We consider first the long-time decay of the time 3687
dx.doi.org/10.1021/jz301672e | J. Phys. Chem. Lett. 2012, 3, 3685−3689
The Journal of Physical Chemistry Letters
Letter
self-correlation function, Cε(t) = ⟨ε(t)ε(0)⟩. When fitted to a sum of two exponentials, the faster relaxation occurs on 0.2−0.5 ps, consistent with the dynamics of bulk SPC/E water.26 The slower exponential relaxation time τE is plotted in Figure 3 (open circles). This relaxation time spikes to nanosecond values in the range of solute−solvent energies corresponding to the second maximum of χG (Figure 2c); the spike is followed by a sudden drop of the relaxation time. We have also calculated the time autocorrelation function of pI1(t), Cp(t) = ⟨pI1(t)pI1(0)⟩. The results are shown by the solid squares in Figure 3. The dynamics of the solute−solvent energy and the dynamics of establishing the collective orientational order in the hydration layer follow the same basic pattern. A large manifold of possible surface orientational states are available to water at hydrophobic interfaces.27 These states are reached by repopulating the configurations that in principle exist, but are less favorable, in the homogeneous liquid.28 Polar surface groups further repopulate these orientational states, restricting the interfacial water mobility by aligning water dipoles along the local electric field. The result is a large negative hydration entropy (Figure 1). The water dipoles are therefore increasingly locked in some energetically preferential positions, eventually leading to the dynamical arrest of the hydration layer and a dramatic slowing down of its dynamics (Figure 3). As the strength of the solute−solvent interactions is further increased, the weight of slow dynamics in the correlation function drops (Figure 2 in the SI) and eventually disappears. The slow dynamics of water rearrangement in the hydration layer is not resolved anymore on the length of the simulation trajectories and contributes little to the overall statistics. At this point, the hydration layer is dynamically arrested on the time-scale of simulations, and only fast librations and vibration near equilibrium positions are resolved. The picture advocated by the present simulations offers interesting extensions to problems involving heterogeneously polar surfaces, such as those of biopolymers.29 One can anticipate that a patchy interface of interchanging polarity will lock the corresponding surface waters in different preferential states, thus producing a complex landscape of surface water configurations. These states will not be completely locked at intermediate surface polarity of the surface groups, and the hydration layer will potentially flip between different energetic states, similarly to collective relaxation events leading to the slow relaxation shown in Figure 3. These collective domain motions will occur on time-scales significantly exceeding relaxation times of one-particle orientational or translational dynamics. Relaxation of polar domains will affect observables sensitive to dipolar interfacial polarization, such as Stokes-shift dynamics. Slow relaxation tails, ranging from tens of picoseconds to nanoseconds are indeed observed in the Stokes-shift dynamics of proteins.30−32 These relaxation functions are, however, strongly affected by the protein motions33−35 and it is currently not clear if the scenario of flipping water domains is realized in the hydration shells of biomolecules. The present simulations show that surface hydration carries thermodynamic and dynamic signatures distinct from those of hydration of small multipolar solutes. Surface hydration provides a context in which the electrostatic hydration can become extremely nonlinear, and in which the dynamics of the collective dipolar response of the hydration water can slow dramatically. Linear relations between solvation thermodynamic potentials found here also imply the proportionality between the energy (enthalpy) and entropy of solvation, a
phenomenon known as enthalpy−entropy compensation.36,37 If this type of proportionality is general for surface solvation, it might contribute to the understanding of this broadly observed phenomenon.37 From a general perspective, surface solvation is obviously important for the stability of solutions of large (colloidal) solutes that are nominally nonpolar, and would be predicted to be insoluble, when judged from their global multipole moments.
■
ASSOCIATED CONTENT
S Supporting Information *
Description of the simulation protocol and additional simulation results. This material is available free of charge via the Internet at http://pubs.acs.org.
■
AUTHOR INFORMATION
Corresponding Author
*E-mail:
[email protected]. Notes
The authors declare no competing financial interest.
■
ACKNOWLEDGMENTS This research was supported by the National Science Foundation (CHE-1213288). CPU time was provided by the National Science Foundation through TeraGrid resources (TGMCB080116N). The authors are grateful to Daniel Martin for his help with the setup of the simulations.
■
REFERENCES
(1) Dill, K. A. Dominant Forces in Protein Folding. Biochemistry 1990, 29, 7133−7155. (2) Ball, P. Water as an Active Constituent in Cell Biology. Chem. Rev. 2008, 108, 74. (3) Hunter, R. J. Foundations of Colloid Science, 2nd ed.; Oxford University Press: Oxford, 2001. (4) Fennell, C. J.; Dill, K. A. Physical Modeling of Aqueous Solvation. J. Stat. Phys. 2011, 145, 209−226. (5) Wang, J.; Bratko, D.; Luzar, A. Length-Scale Dependence of Hydration Free Energy: Effect of Solute Charge. J. Stat. Phys. 2011, 145, 253−264. (6) Hummer, G.; Pratt, L. R.; Garcia, A. E. Molecular Theories and Simulation of Ions and Polar Molecules in Water. J. Phys. Chem. A 1998, 102, 7885. (7) Ben-Amotz, D.; Raineri, F. O.; Stell, G. Solvation Thermodynamics: Theory and Applications. J. Phys. Chem. B 2005, 109, 6866. (8) Rajamani, S.; Ghosh, T.; Garde, S. Size Dependent Ion Hydration, Its Asymmetry, and Convergence to Macroscopic Behavior. J. Chem. Phys. 2004, 120, 4457−4466. (9) Matyushov, D. V.; Newton, M. D. Understanding the Optical Band Shape: Coumarin-153 Steady-State Spectroscopy. J. Phys. Chem. A 2001, 105, 8516−8532. (10) Lynden-Bell, R. M.; Rasaiah, J. C. From Hydrophobic to Hydrophilic Behaviour: A Simulation Study of Solvation Entropy and Free Energy of Simple Solutes. J. Chem. Phys. 1997, 107, 1981−1991. (11) Ashbaugh, H. S. Convergence of Molecular and Macroscopic Continuum Descriptions of Ion Hydration. J. Phys. Chem. B 2000, 104, 7235. (12) Beck, T. L. The Influence of Water Interfacial Potentials on Ion Hydration Free Energies and Density Profiles Near the Surface. arXiv:1209.4516. (13) Takashima, S. Electrical Properties of Biopolymers and Membranes; Adam Hilger: Bristol, U.K., 1989. (14) Barlow, D. J.; Thornton, J. M. Charge Distribution in Proteins. Biopolymers 1986, 25, 1717−1733. (15) Heringa, J.; Argos, P. Side-Chain Clusters in Protein Structures and Their Role in Protein Folding. J. Mol. Biol. 1991, 220, 151−171. 3688
dx.doi.org/10.1021/jz301672e | J. Phys. Chem. Lett. 2012, 3, 3685−3689
The Journal of Physical Chemistry Letters
Letter
(16) Takashima, S. Electric Dipole Moment of Globular Proteins: Measurement and Calculation with NMR and X-ray Databases. J. NonCryst. Solids 2002, 305, 303−310. (17) Athawale, M. V.; Jamadagni, S. N.; Garde, S. How Hydrophobic Hydration Responds to Solute Size and Attractions: Theory and Simulations. J. Chem. Phys. 2009, 131, 115102−9. (18) Gray, C. G.; Gubbins, K. E. Theory of Molecular Liquids; Clarendon Press: Oxford, 1984. (19) Matyushov, D. V.; Ladanyi, B. M. A Perturbation Theory and Simulations of the Dipole Solvation Thermodynamics: Dipolar Hard Spheres. J. Chem. Phys. 1999, 110, 994. (20) Ben-Amotz, D.; Underwood, R. Unraveling Water’s Entropic Mysteries: A Unified View of Non-polar, Polar and Ionic Hydration. Acc. Chem. Res. 2008, 41, 957−967. (21) Beck, T. L. Hydration Free Energies by Energetic Partitioning of the Potential Distribution Theorem. J. Stat. Phys. 2011, 145, 335−354. (22) Yu, H.-A.; Karplus, M. A Thermodynamic Analysis of Solvation. J. Chem. Phys. 1988, 89, 2366. (23) Asthagiri, D.; Merchant, S.; Pratt, L. R. Role of Attractive Methane−Water Interactions in the Potential of Mean Force Between Methane Molecules in Water. J. Chem. Phys. 2008, 128, 244512−7. (24) Friesen, A. D.; Matyushov, D. V. Non-Gaussian Statistics of Electrostatic Fluctuations of Hydration Shells. J. Chem. Phys. 2011, 135, 104501. (25) Chopra, G.; Levitt, M. Remarkable Patterns of Surface Water Ordering Around Polarized Buckminsterfullerene. Proc. Natl. Acad. Sci. U.S.A. 2011, 108, 14455−14460. (26) Gallo, P.; Rovere, M.; Spohr, E. Supercooled Confined Water and the Mode Coupling Crossover Temperature. Phys. Rev. Lett. 2000, 85, 4317−4320. (27) Lee, C. Y.; McCammon, J. A.; Rossky, P. J. The Structure of Liquid Water at an Extended Hydrophobic Surface. J. Chem. Phys. 1984, 80, 4448−4455. (28) Sharp, K. A.; Vanderkooi, J. M. Water in the Half Shell: Structure of Water, Focusing on Angular Structure and Solvation. Acc. Chem. Res. 2010, 43, 231. (29) Giovambattista, N.; Lopez, C. F.; Rossky, P. J.; Debenedetti, P. G. Hydrophobicity of Protein Surfaces: Separating Geometry from Chemistry. Proc. Natl. Acad. Sci. 2008, 105, 2274−2279. (30) Pal, S. K.; Zewail, A. H. Dynamics of Water in Biological Recognition. Chem. Rev. 2004, 104, 2099−2123. (31) Toptygin, D.; Gronenborn, A. M.; Brand, L. Nanosecond Relaxation Dynamics of Protein GB1 Identified by the TimeDependent Red Shift in the Fluorescence of Tryptophan and 5Fluorotryptophan. J. Phys. Chem. B 2006, 110, 26292−26302. (32) Abbyad, P.; Shi, X.; Childs, W.; McAnaney, T. B.; Cohen, B. E.; Boxer, S. G. Measurement of Solvation Responses at Multiple Sites in a Globular Protein. J. Phys. Chem. B 2007, 111, 8269−8276. (33) Li, T.; Hassanali, A. A.; Singer, S. J. Origin of Slow Relaxation Following Photoexcitation of W7 in Myoglobin and the Dynamics of its Hydration Layer. J. Phys. Chem. B 2008, 112, 16121−16134. (34) Matyushov, D. V. Nanosecond Stokes Shift Dynamics, Dynamical Transition, and Gigantic Reorganization Energy of Hydrated Heme Proteins. J. Phys. Chem. B 2011, 115, 10715−10724. (35) Martin, D. R.; Matyushov, D. V. Non-Gaussian Statistics and Nanosecond Dynamics of Electrostatic Fluctuations Affecting Optical transitions in a Green Fluorescent Protein. J. Phys. Chem. B 2012, 116, 10294−10300. (36) Grunwald, E.; Steel, C. Solvent Reorganization and Thermodynamic Enthalpy-Entropy Compensation. J. Am. Chem. Soc. 1995, 117, 5687−5692. (37) Sharp, K. Entropy−Enthalpy Compensation: Fact or Artifact? Protein Sci. 2001, 10, 661−667.
3689
dx.doi.org/10.1021/jz301672e | J. Phys. Chem. Lett. 2012, 3, 3685−3689
