10272
J. Phys. Chem. B 2008, 112, 10272–10279
Water Inside a Hydrophobic Cavitand Molecule Jeffrey Ewell, Bruce C. Gibb, and Steven W. Rick* Department of Chemistry, UniVersity of New Orleans, New Orleans, Louisiana 70148 ReceiVed: May 19, 2008; ReVised Manuscript ReceiVed: June 18, 2008
The structure and dynamics of water inside a water-soluble, bowl-shaped cavitand molecule with a hydrophobic interior are studied using molecular dynamics computer simulations. The simulations find that the number of inside water molecules is about 4.5, but it fluctuates from being completely empty to full on a time scale of tens of nanoseconds. The transition from empty to full is energetically favorable and entropically unfavorable. The water molecules inside have fewer hydrogen bonds than the bulk and in general weaker interactions; the lower energy results from the nearest-neighbor interactions with the cavitand atoms and the water molecules at the entrance of the cavitand, interactions that are lost upon dewetting. An analysis of translational and rotational motion suggests that the lower entropy of the inside water molecules is due to decreased translational entropy, which outweighs an increased orientational entropy. The cavitand molecule acts as a host binding hydrophobic guests, and dewetting can be induced by the presence of a hydrophobic guest molecule about 3 Å above the entrance. At this position, the guest displaces the water molecules which stabilize the inside water molecules and the empty cavitand becomes more stable than the full. I. Introduction Deep-cavity cavitand 1 is a bowl-shaped host molecule that can bind a wide variety of hydrophobic guests.1 The cavitand is water-soluble by virtue of eight carboxylic acid groups on its surface and has a large hydrophobic concave binding pocket defined by twelve aromatic rings. The interior is roughly 8 Å deep, 8 Å wide at the rim, and narrow enough at the bottom as to be closed to the passage of molecules (as represented in 1).
Depending on the hydrophobicity and size of a guest, complexes of different host-guest stoichiometries can be formed by 1. If the guest is very small2 or amphiphilic3 the host forms 1:1 complexes, whereas with larger and/or less polar guests, the host dimerizes to form capsular complexes. These supramolecular species, in which one or more guest molecule is contained within the hydrophobic interior of the capsule, have been shown to bind drug molecules,4 separate hydrocarbon gases,2 and make excellent nanoscale reactors for bringing about novel chemistries.5 Influencing the binding thermodynamics is the water structure in the binding site of the host. The geometric constraints near * To whom correspondence should be addressed. E-mail:
[email protected].
such concave surfaces leads to different water properties than water near flat or convex surfaces.6–10 There has been significantly less studies of water near concave hydrophobic surfaces than there has been for convex or flat surfaces or pores. These studies have primarily been for water inside model hydrophobic spheres (with radii from 6 to 24 Å)6,7 or near hydrophobic hemispheres (with radii varying from 6.5 to 12 Å).8,9 These studies find that water molecules have fewer hydrogen bonds and a higher energy than bulk water. In contrast, the number of hydrogen bonds is the same as the bulk for water near small methane-sized hydrophobic spheres but decreases as the curvature increases, indicating a change in the mechanism of the hydrophobic effect as solute size increases.11–16 The different energetics of water near concave surfaces may have other thermodynamic consequences. In the simulations of Paschek, water molecules in concave regions between hydrophobic solute pairs have an increased heat capacity relative to the bulk and to water molecules next to convex hydrophobic surfaces.10 The water in concave hydrophobic environments may also have different dynamical properties than the bulk, although this is not clear and may depend on the details of the system. Chau finds similar translational but faster orientational dynamics for water in hydrophobic hemispheres,8 and Brovchenko et al. find slower translational motion for water in hydrophobic spheres.6 Faster orientational dynamics would be consistent with fewer hydrogen bonds. For water near both small and large spherical solutes, experiments17–19 and simulations12,16 find slower translational and rotational dynamics. So the concave hydrophobic environment limits the number of neighboring water molecules resulting in a decreased number of hydrogen bonds, which allows for faster rotational motion but, maybe, slower translational motion. The fewer number of hydrogen bonds for water molecules near concave surfaces may make them easier to dewet the surface. Simulations of these systems do show water densities near the surface intermediate between wet and dry7,8 and demonstrate substantial oscillations in the number of water molecules, with the empty, or dry, state being observed, but
10.1021/jp804429n CCC: $40.75 2008 American Chemical Society Published on Web 07/29/2008
Water Inside a Hydrophobic Cavitand Molecule rare, in long simulations.9 Simulations of water in model hydrophobic pores20 and carbon nanotubes21–23 in which the confined environment also leads to a reduced number of hydrogen bonds show fluctuations between empty and full states. In both the concave and pore studies, the water-solute interactions are critical, with purely repulsive potentials leading to dewetting but with weak attractive interactions, corresponding to more realistic interactions, leading to surface wetting.6,7,21,22 Dewetting fluctuations are seen in simulations of the folding of proteins and other polymers,24–27 although it is not established that folding proceeds through dewetting as opposed to a gradual water expulsion mechanism.28 In the simulations of Miller et al. the rate-limiting event in hydrophobic polymer collapse involves a fluctuation in the local water density in concave bend on the polymer.27 More related to host-guest binding, the study of Setny of the potential of mean force between a methane molecule and a hydrophobic hemisphere find that the surface completely dewets as the methane center is within 1 Å of the top of the cavity (as indicated by the center of the cavity atoms), with partial dewetting occurring as larger distances.29 In this study, we use molecular dynamics computer simulations to study water in the concave hydrophobic surface of the octaacid cavitand (1). We aim to characterize the water density, the energetic factors that lead to hydration versus dewetting, and the dynamical properties of the water molecules inside the cavitand. We also describe how the water is affected by the presence of a guest molecule as it approaches the host. II. Methods The molecular dynamics simulations were performed using the AMBER 7 simulation package,30 with the TIP4P-Ew potential for water.31 The TIP4P-Ew model accurately reproduces many of the important properties of water over a range of temperatures31 and also accurately reproduces the solvation properties for simple nonpolar solutes.32 Charges for the octaacid molecule were generated from a RESP charge fitting procedure with input from Hartree-Fock calculations at the 6-31G* level using the Gaussian03 program.33 Additional parameters for the octaacid and ethane molecules were generated using the gaff parameter set.34 The initial coordinates for the octaacid 1 were from an unpublished X-ray structure. Of the eight carboxylic acids, the four at the top rim of the structure, which are extended out into the solvent, were assumed to be unprotonated. Of the four at the bottom, which are in relatively close proximity, two, taken to be diagonally across from each other, were assumed to be protonated. This gives a total net charge of -6. Charge neutrality of the system was created by adding 6 sodium ions, using the Amber 99 parameter set.35 The system contained 1178 water molecules. Simulations were carried out primarily in the isothermal-isobaric (T,P,N) ensemble, for lengths of 10 ns, except for microcanonical (E,V,N) simulations performed for determining dynamical properties. A 1-fs time step was used, and long-ranged electrostatics were treated with particle mesh Ewald. A cutoff for the Lennard-Jones and the real-space part of the Ewald interactions equal to 9.0 Å is used. All analyses of the simulations were done with our own programs. The octaacid molecule was free to translate and rotate, as well as change conformation, so the spatial analysis of water structure were done by rotating the simulation coordinates onto the original octaacid structure, which has its C4 axis along the z axis. The x axis is defined to be along the direction connecting two opposite carbon atoms connected to the carboxylic acids on the bottom of the molecule. The atoms at the top of the molecule roughly form a square; the x axis connects across the
J. Phys. Chem. B, Vol. 112, No. 33, 2008 10273 middle of two sides. By symmetry, the y axis would be equivalent to the x axis. All error estimates represent one standard deviation of the data. Simulations were also carried out which excluded water from the interior of the cavitand. Water were kept out of the interior using a short ranged repulsive potential inside the cavity, E(z) ) Rz12, where z is the distance from a plane perpendicular to the C4 axis, with R equal to 2345 Å12 kcal/mol. The plane is positioned at a distance d from the bottom of the cavitand. III. Results Water Density Structure in the Cavitand. Figure 1 shows the water density as a function of x and z (see Methods). The density shows the average position of the water oxygen positions. The 4-fold symmetry of the molecule was used to generate the density. The origin of z corresponds to the position of the four methine carbon atoms from which the -CH2CH2-CO2 groups connect. A number of conclusions can be made from this figure. First, water does enter the cavitand to a significant amount. There is no water density at the bottom of the cavitand, indicating that the bottom of the cavitand is closed to the passage of water. Inside the cavitand, two separate regions of water density are apparent. The two regions are separated by a region of low density right below z ) 4 Å, corresponding to the positions of four benzal hydrogen atoms. The lower regions contains a single, considerably localized, water molecule. The upper region contains more water molecules, the number of which depends on the definition of the interior region. The choice of value of z that separates interior water molecules from the others is arbitrary. We chose to use the positions of the 8 oxygen atom on the 32-member ring of the rim, which are the uppermost heavy atoms. These have an average z position equal to 7.3 Å. Notice that the excluded volume of the cavitand molecule extends beyond this distance, to almost 10 Å. That value of z was chosen because it is in the region where the cavitand-water surface changes from convex to concave, this being one way to define interior rather than exterior solute surfaces, and also because above this position there are always water molecules, whereas below there are configurations for which there are no water molecules. In addition, because this is the position of oxygen atoms, this could be considered the upper limit of the hydrophobic part of the cavitand. With this definition, there are on average 4.34 ( 0.04 water molecules, although this number fluctuates considerably. For the 10-ns simulation, the cavitand is empty of water molecules for only a small fraction of the time, 3 out of 10000 stored configurations. (These represent distinct events, separated by configurations in which there are water molecules in the cavitand.) From this, the free energy of the filled relative to empty cavitand can be estimated from ∆G ) -kT ln[(10000 3)/3) or -4.8 ( 0.6 kcal/mol. Because the empty configurations are so rare, there is a good deal of uncertainty in this number. This is in addition to its arbitrariness, due to the ambiguity of the definition of interior. However, the simulations do show that a dewetting transition or a liquid-vapor oscillation takes place on a time scale of approximately 5 ns. Less ambiguous is the presence of a lower-most water molecule, between 0 and 4.0 Å. A water molecule is present in this position a fraction 0.88 ( 0.04 of the time, giving a free energy difference (occupied relative to empty) equal to -kT ln[0.88/(1 - 0.88)] ) -1.2 ( 0.1 kcal/mol. From the probabilities of having N water molecules in the cavitand, P(N), a free energy relative to the empty state can be found from -kT ln[P(N)/P(0)]. The probabilities and free energies are shown in Figure 2. The most
10274 J. Phys. Chem. B, Vol. 112, No. 33, 2008
Ewell et al.
Figure 1. Two-dimensional density of water molecules in the vicinity of the octaacid cavity. The density is in units relative to the bulk. The plot was prepared with the XFARBE program.36
probable, or lowest free energy, number of water molecules is 4, with 5 being very close. The highest number of water molecules seen inside is 7. One measure of the effect of the confined environment on the water structure is the number of hydrogen bonds between water molecules. We use the Mancera and Buckingham definition of a hydrogen bond: a hydrogen bond exists if the oxygen-oxygen distance is less than 3.6 Å and the angle between the O-H vector on the hydrogen bonding donor and the O-O vector is between 130 and 180°.37 The most confined water is in the lowest position; this water has on average 1.33 ( 0.02 hydrogen bonds. The lower water tends to form hydrogen bonds as a hydrogen bond acceptor. The other water molecules inside the cavitand (from around 4 to 7.3 Å) have on average 2.65 ( 0.01 hydrogen bonds. In the region at the top of the cavitand, from about 7.3 to 10.3 Å, the number of hydrogen bonds is increased to 3.33 ( 0.01, closer to the bulk value of 3.64 ( 0.02. (Bulk properties are calculated using the ten water molecules furthest from the center of the cavitand.) Because the water molecules in the concave part of the cavitand have one less hydrogen bond than bulk water, the energies of these molecules must be significantly different.
The number of hydrogen bonds with neighboring waters depends both on the geometric constraints of the cavitand, which limits the number of water molecules and the density of the water molecules within the available volume. From Figure 1 it is apparent that the density in the lower part of the cavity is structured, with regions of high density (around x ) 0 Å, z ) 2 Å, and x ) (1 Å, z ) 5 Å) and regions of low density (x ) 0 Å, z ) 6 Å). Above z ) 7 Å, the density appears similar to the bulk value (above z ) 11 Å). The structure of the water can be further characterized by F(z), the number of water molecules that are inside the cavitand and have an oxygen position between z and z +δz [Figure 3A]. Three sections of water density can be seen. A section from 0 to 4 Å integrates to give 0.88 ( 0.04 water molecules (as reported above). A middle section from 4 to 6.6 Å has 2.82 ( 0.04 water molecules, and the a top section from 6.6 to 9 Å contains 3.21 ( 0.02 molecules. To make comparisons to the bulk density, some estimates of the volume of the cavitand must be made. For a particular value of z, a volume can be found from the area, A(z), of the plane tangential to the z axis, which is within the solventaccessible surface38 of the cavitand. To calculate the solventaccessible surface, we use a water radius of 1.4 Å and radii for
Water Inside a Hydrophobic Cavitand Molecule
Figure 2. (A) The probability of having N water molecules inside the cavity and the (B) the free energy for N water molecules relative to zero water molecules. The line in part B is a fourth-order polynomial fit to the free energy data, and the line in part A is from the same polynomial.
J. Phys. Chem. B, Vol. 112, No. 33, 2008 10275 similar values. The average energy for a water molecule in this position is -15.31 ( 0.03 kcal/mol (the water-water contribution is -8.57 ( 0.04, and the water-octaacid contribution is -6.74 ( 0.03). For the other water molecules in the concave section of the cavitand from around 4 to 7.3 Å, the energy is also noticeably higher than the bulk, with an average value of -18.99 ( 0.02 kcal/mol (the water-water energy is -16.16 ( 0.02, and the water-octaacid energy is -2.83 ( 0.02 kcal/ mol). The energetically stable position in the cavitand around 4.4 Å is stabilized by favorable water-water interactions between a water molecule at that position and the water localized in the bottom position of the cavitand. This analysis shows that the energy lost from having fewer water neighbors inside the cavitand is not completely compensated for by the interaction with the cavitand atoms, which is not surprising given the nonpolar nature on the cavitand interior. It is also apparent that the regions of high water density (near 2 and 5 Å) do not correspond to regions where the energy is relatively low. Clearly, water molecules are in places that are not low in energy themselves, so they must lower the free energy of the entire system by being in these positions. Given that the energy of the water molecules inside the cavitand is higher than the bulk water, the question is what factors favor hydration of cavity. To fully understand the energetic contribution to the hydration process, the energy of the entire system with the empty and full cavitand should be compared. However, the difference from the movement of a few water molecules would be very small relative to the total energy of the system, which contains over a thousand water molecules. For this reason, we use a method which only includes the small subset of water molecules in the vicinity of the cavity. We split up the water molecules into three types. There are the n water molecules inside the cavitand, the m water molecules at the boundary between the inside and bulk water molecules, and the N - n - m remaining or bulk (even though these include many water molecules near the solute) waters. The potential energy as a function of n is ww ow 〈E(n) 〉 ) 〈Einside (n) 〉 + 〈Einside (n) 〉 + 〈Eww boundary(n,m) 〉 +
Figure 3. (A) The density of water molecules inside the cavity, expressed in terms of F(z), the number of water molecules which have an oxygen position between z and z +δz (solid line) and F(z) scaled by the area within the solvent-accessible surface, A(z) (dashed line). (B) The total (solid line), water-water (dashed line), and water-octaacid (dotted line) interaction energy of a water molecule as a function of its oxygen position along the z axis inside the cavitand. The constant dot-dashed line is the bulk total energy.
the cavitand atoms were calculated using the method of Hummer et al. in which the radii correspond to the distance where the water-solute interaction equals 1 kT.39 This gives radii equal to 1.1, 1.5, and 1.6 Å for hydrogen, oxygen, and carbon atoms, respectively. The resulting density (the dotted line in Figure 3A) shows that the lowest and second lowest position have densities slightly higher than the bulk density, and near the top the density is about equal to the bulk value (1 g/cm3). Energetics of the Cavitand Water Molecules. The energy of a molecule depends on its position in the cavitand, as illustrated in Figure 3B. The bulk energy is -22.26 ( 0.03 kcal/ mol, -21.96 ( 0.03 kcal/mol from the water-water interactions and -0.31 ( 0.02 kcal/mol from the water-octaacid interactions. The energy is less than the bulk value inside the cavitand and only reaches the bulk value at positions outside the cavitand, above 10 Å. The energy is highest for the lower-most water, where the water-water and water-octaacid energies have
ww 〈Eow boundary(n, m) 〉 + 〈Ebulk(N - n - m) 〉 + oo 〈Eow bulk(N - n - m) 〉 + 〈E (n)〉 (1)
where the superscripts ww, ow, and oo indicate water-water, octaacid-water, and octaacid-octaacid energies, respectively. ww The energy 〈Einside (n)〉 is the sum of all the water-water interactions involving the inside water molecules, 〈Eww boundary(n,m)〉 is the sum of the water-water interactions between the boundary and bulk water molecules, and 〈Eww bulk(N-n-m)〉 is the sum of interactions only between the bulk water molecules, all averaged over those configurations with n inside water molecules. Our ow simulations find that both 〈Eboundary (n, m)〉 and 〈Eoo(n)〉 are independent of n, with this assumption the energy difference between the hydrated and empty cavitand is ww ow 〈∆E(n) 〉 ) 〈E(n) 〉 -〈E(0) 〉 ) 〈Einside (n) 〉 + 〈Einside (n) 〉 + ww ww 〈Eww boundary(n, m) 〉 -〈Eboundary(0, m) 〉 -〈Ebulk/Nbulk 〉 n -
〈Eow bulk/Nbulk 〉 n (2) The change in energy involves the gain in energy from the water-octaacid interactions and the water-water interactions involving the inside water molecules, a loss of the bulk water-octaacid and water-water interactions, and a modification of the water-water interactions in the boundary. The loss
10276 J. Phys. Chem. B, Vol. 112, No. 33, 2008
Ewell et al. closest to the inside, has an average of 3.34 ( 0.02 hydrogen bonds, 2.38 ( 0.03 with other boundary water molecules and 0.96 ( 0.03 with the inside water molecules. With n ) 0, the water molecule at the bottom of the boundary region has 2.63 ( 0.03 hydrogen bonds, all, obviously, with other boundary region water molecules. So the one hydrogen bond lost with the absence of an inside water molecule is partially made up (with only 0.25 of a hydrogen bond) through more hydrogen bonds with waters in the boundary region. The energies needed as input to eq 2 are given in Table 1. The energies for each value of n are given, along with the weighted average for each energy 7
〈A 〉 ) Figure 4. The total interaction energy for a water molecule with four (solid line) and no (dashed line) water molecules inside the cavitand. The dotted line shows the total interaction energy with four inside water molecules minus the interaction with the four inside water molecules.
in the bulk water-water energy is only half the water-water energy (-21.94 kcal/mol) to avoid overcounting. This factor of a half is consistent with the fact that, as a water molecule leaves the liquid phase, the neighboring water molecules will adjust to fill the space left by the exiting water molecule and remake their interactions. This is in contrast to the water molecules inside the cavitand. When they exit the cavitand, the water-octaacid interactions and the water-water interactions, most significantly the interactions between the inside and boundary water molecules, are lost. On average, there are 3.15 ( 0.03 hydrogen bonds between the inside and boundary water molecules. How the boundary water molecules respond when these interactions are lost will affect the energy change through ww ww the terms 〈Eboundary (n, m)〉 - 〈Eboundary (0, m)〉. To answer that question, we have to define the boundary region. Figure 4 shows the total energy of a water molecule as a function of its z position within the cavitand averaged for the configurations with n ) 0 (the dashed line) and n ) 4 (the solid line). To generate more configurations with n ) 0, a plane inside the cavitand (at a position z ) 6 Å) with a short-ranged repulsive interaction on the water molecules was used to eliminate waters from the cavitand, as described in Methods. We started with the three n ) 0 configurations and four additional configurations with water molecules all above 6 Å (equilibrating for 3 ps first to empty the cavitand of waters), taken from the unbiased simulations, and ran for 0.4 ns, to give a variety of starting configurations and 1.4 ns of total simulation time. The energies become independent of the number of inside water molecules once a distance over 10 Å is reached, and so the boundary region is taken to be from 7.3 to 10.3 Å. This region contains on average 8 water molecules, so we set m ) 8. Also shown is the energy of water molecules with n ) 4 less the interaction with those four inside water molecules (the dotted line). The dotted line then shows what the energy of the boundary water molecules would be if there were no changes in structure between the full and empty cavitands. The fact that the dotted line is slightly higher in energy than the dashed line indicates that the water molecules do reorientate to partially make up for the interactions lost by the removal of the interior water molecules. This same effect is seen in the hydrogen bond averages for water molecules in boundary region. With n ) 4, a water molecule at the bottom of the boundary region, and so
∑ A(n)e-∆G(n)⁄kT
n)1
7
∑e
(3)
-∆G(n)⁄kT
n)1
For the free-energy change, the average given is the total free energy change between the empty and hydrated ww cavitand as given earlier. The energy 〈∆Eboundary 〉 is equal to ww ww ww 〈Eboundary(n)〉 - 〈Eboundary(n ) 0)〉, with 〈Eboundary (n ) 0)〉 ) -137.3 ( 1.2 kcal/mol. This energy only contributes a small amount to the overall energy change, about 1 kcal/mol, indicating that there is not a lot reorganization of water molecules in the boundary region. The amount of the total water-water energy involving the n inside water molecules, ww 〈Einside 〉, that comes from interactions with the m boundary ww water molecules, 〈Einside-boundary 〉, is also shown. The energy change, 〈∆E(n)〉, for all n except for n ) 1 (but here the error bars are large) is negative, so the process of hydrating the cavitand is exothermic, even though the interactions of water molecules inside the cavitand are less than they are in bulk. There is also a PV contribution to the enthalpy change, arising from the larger volume of the n ) 0 state from the larger number of water molecules in the bulk. The molecular volume of liquid water is 29.9 Å3, and the PV energy per molecule is only 4 × 10-4 kcal/mol. The negative energy change can be attributed to two factors. The water-water interactions for the inside water molecules which average to -54.5 ( 0.2 kcal/mol are less than the water-water interactions the water molecules have in the bulk, -47.70 ( 0.04 kcal/mol. This is largely due to the fact that the water molecules in the boundary region do not make new interactions when the interior water molecules are absent and the energy between the inside and boundary water molecules (-11.3 ( 0.2 kcal/mol) is lost. The boundary water molecules change structure to make up only about 1 kcal/mol of that energy. The empty cavitand creates a small vapor-liquid ww interface costing about -10 kcal/mol in energy (〈Einside-boundary 〉 ww - 〈∆Eboundary〉). Some (about 4.2 kcal/mol) of that energy must be compensated for by the stronger interactions in the bulk so the total contribution to the energy change from the water-water ww ww interactions is -5.7 ( 1.3 kcal/mol (〈Einside 〉 - 〈∆Eboundary 〉ww n〈Ebulk/Nbulk〉). The other factor is the lost interactions between the interior water molecules and the octaacid molecule (-15.4 ( 0.2 kcal/mol), which is significantly greater than the interaction between bulk water molecules and the octaacid (-1.35 ( 0.05 kcal/mol). This factor then contributes about -14 kcal/ mol, which together with the water-water contribution makes up the total energy change of about -20 kcal/mol. If the enthalpy change is as large as -20 kcal/mol, then there must be an entropy change, -T∆S, of about 15 kcal/mol to give ∆G ) -5 kcal/mol. The large entropy change would mean
Water Inside a Hydrophobic Cavitand Molecule
J. Phys. Chem. B, Vol. 112, No. 33, 2008 10277
TABLE 1: Energetic Contributions to Cavitand Hydration (in kcal/mol) as a Function of the Number of Inside Water Molecules na n ww 〈Einside 〉 ww 〈Einside-boundary 〉 ow 〈Einside 〉 ww 〈∆Eboundary〉 ww n〈Ebulk /Nbulk〉 ow n〈Ebulk /Nbulk〉
〈∆E(n)〉 〈∆G(n)〉 〈∆HT(n)〉 a
1
2
3
4
5
6
7
average
-12.5(2.1) -7.4(0.7) -2.3(1.9) 5.9(4.9) -10.97(0.02) -0.31(0.02) 2.5(5.6) -1.1(0.5) -14.6(7.8)
-21.0(0.6) -6.5(0.5) -7.2(0.4) 1.6(3.0) -21.94(0.03) -0.62(0.04) -4.1(3.0) -2.5(0.4) -15.9(6.2)
-34.9(0.4) -9.2(0.4) -11.5(0.4) 1.1(1.3) -32.91(0.05) -0.93(0.06) -11.4(1.3) -3.6(0.4) -18.1(6.0)
-49.5(0.2) -11.0(0.3) -14.4(0.2) 0.8(1.5) -43.88(0.06) -1.24(0.08) -18.0(1.5) -4.3(0.4) -20.3(6.0)
-63.6(0.4) -12.3(0.3) -17.5(0.4) 1.3(1.4) -54.85(0.08) -1.55(0.10) -23.4(1.4) -4.2(0.4) -21.5(6.0)
-79.8(0.9) -13.4(0.3) -19.3(0.7) 1.4(1.5) -65.82(0.09) -1.86(0.12) -30.2(1.5) -3.3(0.4) -22.4(6.0)
-101.0(2.0) -17.2(1.5) -17.6(1.4) 1.2(3.0) -76.79(0.11) -2.17(0.14) -38.5(3.5) -1.5(0.5) -21.9(6.7)
-54.5(0.2) -11.3(0.2) -15.4(0.2) 1.0(0.8) -47.70(0.04) -1.35(0.05) -19.8(0.8) -4.8(0.6) -20.5(3.4)
Numbers in parentheses give error estimates.
that there is a tendency for the cavitand to be less occupied with water as temperature increases. To test this, we ran additional simulations at 313 and 328 K. The simulations did demonstrate a trend for the cavitand to be empty with a greater probability as temperature increased. The probability of being empty increases from 0.0003 ( 0.0002 at 298 K to 0.003 ( 0.002 at 313 K to 0.007 ( 0.005 at 328 K. The entropy change, 〈∆S(n)〉, can be found from the slope (as determined from a linear fit) of 〈∆G(n)〉. An estimate of the enthalpy change can then be found from 〈∆G(n)〉 + T〈∆S(n)〉. This is shown in Table 1, as ∆HT. The agreement between the energy change from eq 2 and from the temperature dependent data is very good. Both give -20 kcal/mol. The agreement is not as good for some of the less common hydration states, particularly n ) 1, but both methods agree that not only is the average energy change negative but also that it tends to decrease as n increases. Water Structure in the Presence of a Guest Molecule. Because the inside water molecules are stabilized through interactions with water molecules just outside the cavitand, disrupting the water structure outside the cavitand could lead to less water molecules inside. To test how a guest molecule influences the water structure inside the cavitand, we performed simulations with using ethane as a potential guest. Ethane is the one (purely) hydrocarbon guest that has been shown to form a stable 1:1 complex with cavitand 1.2 (Methane does not bind to the pocket, whereas propane triggers dimerization of the host to form a 2:2 capsular complex.)2 The center-of-mass of the ethane molecule in constrained at various distances along the C4 (z) axis from the bottom of the cavitand, using a harmonic restraint, with force constant 2.0 kcal/mol/Å2. Four separate simulations of 10 ns each were performed, at distances of 10, 11, 12, and 14 Å. From these simulations, the fraction that the cavitand was empty as a function of ethane distance was calculated (Figure 5). The ethane molecule does induce dewetting of the cavitand when it is in a position to displace the boundary region water molecules. The hydration of the cavitand is not changed until the ethane gets to a point around 11 Å. Beyond 12 Å the probability of finding an empty cavitand is too low for it to show in the probability histograms. For comparison, the probability of the empty cavitand without an ethane molecules (0.0003) is shown; this is essentially the same probability as that when the ethane is constrained to be 12 Å or farther away. At a distance of 11 Å, the number of inside water molecules oscillates from being fully occupied, with 4 or more molecules, to empty, on a time scale of 1.1 ( 0.4 ns. When the ethane molecule is closer than 10 Å to the bottom of the cavitand, the probability to be empty increases to almost one. When the ethane molecule is actually inside the cavitand (unconstrained), there is not much room for water molecules,
Figure 5. The probability of an empty cavitand as a function of ethane center-of-mass distance along the z axis of the cavitand molecules. The dotted line shows the probability in the absence of an ethane molecule.
and there is on average less than one (0.3) water molecule inside. Inside the cavitand, the ethane molecule is positioned with an average center-of-mass position equal to z ) 4.26 ( 0.06, with one carbon atom in the lowermost position (less than 4 Å) a large fraction (0.77 ( 0.04) of the time. When the ethane is in the cavitand, a water molecule was never in the lower position. Dynamics of the Interior Water Molecules. The translational velocity autocorrelation function is given by N
CV(t) )
∑
1 m〈ViCM(t)ViCM(0)〉 N i)1
(4)
where m is the mass of a water molecule and ViCM(t) is the center of mass velocity of molecule i. Rotational times scales can be found from the angular velocity autocorrelation function N
cω(t) )
3
∑∑
1 I 〈ω (t)ωij(0)〉 N i)1 j)1 j ij
(5)
where Iij is the jth principle moment of inertia of molecule i and ωij is angular velocity around the principal axis j of molecule i. These correlation functions are shown in Figure 6 for bulk water, the lowermost water, and the inside water molecules (excluding the lowermost water). From these functions the translational diffusion constant can be found from40
DT ) 1/3m
∫0∞ cT(t)dt
and the rotational diffusion constant can be found from
(6)
10278 J. Phys. Chem. B, Vol. 112, No. 33, 2008
DR )
1 3
∑ Ij
∫0∞ cwtdt
Ewell et al.
(7)
j)1
The values for the translational diffusion constants are 2.6 ( 0.1 × 10 -5 cm2/s (bulk), 1.9 ( 0.2 × 10 -5 cm2/s (inside), and 0.6 ( 0.3 × 10 -5 cm2/s (lowermost). The bulk diffusion constant is in good agreement with the previously reported value for TIP4P-Ew (2.4 ( 0.00 × 10 -5 cm2/s).31 The rotational diffusion constants are 2.1 ( 0.1 × 1011 s -1 (bulk), 3.1 ( 0.1 × 10 11 s -1 (inside), and 33. ( 2 × 10 11 s -1 (lowermost). It is evident that the lower water has much different dynamics than both the bulk and the other inside water molecules. Both the translational and rotational motion resemble that of a damped harmonic oscillator with a period of about 0.5 ps for translations and 0.1 ps for rotations. The lower water moves in a relatively wide potential minimum with a low frequency until it leaves this site. The average lifetime in this site is about 0.8 ps. The rotations are relatively fast as would be consistent with a molecule that forms 1.33 ( 0.02 hydrogen bond, mostly as a hydrogen bond acceptor so the hydrogens are free to rotate. The bulk and inside water time correlation functions both show the same time scales, with peaks and minima at similar points. The main difference between the translational correlation functions is the negative region for the inside water molecules at 0.25 ps indicating that, after this amount of time, the water molecules have had collisions with their nearest neighbors, turned around and are on average moving in the opposite direction from the original point. This negative correlation is much more pronounced for the inside water molecules than the bulk. The turning round leads to a decreased translational diffusion constant. Some of this turning around must be due to collisions with the cavitand surface. For rotational dynamics, the main difference is the stronger negative correlation for bulk relative to the inside water, as given by the deeper minimum at 0.025 ps. This indicates that the rotations are more hindered, or librational, for the bulk waters. Conclusion The simulation results show that the interior of the cavitands are occupied with water molecules, which average about 4 or 5,
Figure 6. The translational velocity autocorrelation function (A) and rotational translational function (B) for bulk (solid line), inside (dashed line), and bottom most water (dotted line).
but fluctuate between 0 and 7. The water molecules have a higher energy than bulk waters (Figure 3B), forming fewer hydrogen bonds than bulk water. The higher energy and few hydrogen bonds is in agreement with previous studies of water in concave hydrophobic environments.6–9 Despite the higher energy of the water molecules inside, the process of adding all the water molecules to the interior is exothermic for two reasons. First, even though the water-octaacid interactions are weak individually, they sum up to about -15 kcal/mol over the 4 or 5 water molecules, and these interactions are lost when the cavitand is empty (see Table 1). The octaacid-water interactions for the bulk water is only -1.35 kcal/mol, so emptying the cavity results in a loss of -14 kcal/mol in the water-octaacid interactions. Second, the interactions between the water molecules in cavitand and water molecules in the boundary between the interior and bulk regions are also lost, as a small liquid-vapor interface is created. This energy is about -11 kcal/mol. Only 1 kcal/mol of this energy is recovered by the water molecules in the boundary reorganizing in the absence of interior water molecules to form more hydrogen bonds with other boundary region water molecules. The small degree of reorganization means there is only a small entropy change associated with the creation of the interface. An additional amount of energy is gained in moving from the interior to bulk water from stronger water-water interaction in the bulk, which compensates for the energy lost in creating the liquid-vapor interface. In total, the ∆H is about -20 kcal/mol, ∆G is about -5 kcal/mol, and -T∆S is 15 kcal/mol. These values, found from an analysis of water populations and energetics at 298 K are consistent with enthalpy changes found from the temperature dependence of ∆G. Therefore, it is energetically unfavorable for the water molecules that transfer from the bulk to the cavitand interior but energetically favorable for the neighboring water and octaacid atoms, resulting in a total energy change that is favorable. On the other hand, because there is not much structural rearrangement of either the boundary water molecules or the octaacid molecule, the entropy change must be mostly due to the entropy difference of the transferring water molecules. Because the overall entropy change is negative, the entropy of the cavitand water molecules must be lower than that of bulk water. An analysis of the translational and orientational velocity time correlation functions reveals that translational motion is slower and that orientational motion is faster for water molecules in the cavitand relative to the bulk. Slower translational motion makes sense for a molecule in a concave environment in which movement is restricted along some directions and the water density also shows increased structure inside the cavitand (Figure 1). Faster orientational motion is consistent with fewer hydrogen bonds because the hydrogens not involved in a hydrogen bond should be freer to rotate. The dynamical results suggest that the decreased entropy of the interior water molecules is due to decreased translational entropy, which must overcompensate any gain in rotational entropy. A correlation between the configurational entropy and the translational diffusion constant, as suggested by this data, is consistent with a number of previous studies.40–43 That the contribution from the octaacid-water interactions are so important for the stability of the inside water molecules agrees with earlier studies of water in concave surface and hydrophobic pores. As mentioned in the Introduction, these studies find that water will only not be found near these surface if the solvent-solute interactions are purely repulsive, adding weak attractions brings the water molecules in contact with the surface.7,6,21,22 The other contribution to the negative energy of filling the cavitand, the interaction with the water molecules at the boundary between the interior and the bulk, can be
Water Inside a Hydrophobic Cavitand Molecule
J. Phys. Chem. B, Vol. 112, No. 33, 2008 10279
eliminated as a potential guest molecule approaches. Our simulations find that as an ethane molecule gets to about 11 Å from the bottom of the cavity, or 4 Å from the top, the water structure is sufficiently disrupted, inducing dewetting of the cavitand (Figure 5). With the ethane molecule in the boundary region, the probability of having an empty cavitand is close to 1. This is similar to the Setny study of a methane molecule near a hydrophobic concave hemisphere, in which the surface becomes empty of water as the methane molecule get close to top of the cavity.29 For the octaacid/ethane system, with the ethane near the entry to the cavitand, transitions between empty and filled, or dry and wet, take place on a 1-ns time scale. These results suggest a mechanism for host-guest binding, in which the guest approaches the guest binding site, promoting fluctuations in which the cavitand becomes empty of water. Once this happens the guest could presumably move to the binding site quickly, without a large energetic barrier. Our results can tell us about the role of water in the thermodynamics of host-guest binding process
host{water}(aq) + guest(aq) h host{guest}(aq)
(8)
where the occupant of the host binding site (water or guest) is indicated inside the curly brackets. The process can be split up into steps
host{water}(aq) h host{}(aq) guest(aq) h guest(gp) host{}(aq) + guest(gp) h host{guest}(aq)
(9) (10) (11)
where eq 9 removes the water from the guest binding site, eq 10 removes the guest from the solvent, and eq 11 places the gas-phase guest into the empty host. For eq 9, our calculations find ∆G ) 5, ∆H ) 20, and -T∆S ) -15 kcal/mol. For the desolvation step for the guest, in the case of ethane, ∆G ) -2, ∆H ) 4, and -T∆S ) -6 kcal/mol.44 For larger hydrophobic guests, the entropy change will be larger. By assumption that the guest occupies the same volume as that created by the water that left in eq 9, the final step should not involve substantial rearrangement of the water and so the contribution to ∆S from the solvent is about zero. There will be solvent contribution to ∆H from the interactions between the guest and the water molecules bordering the binding site. Overall, for the host-guest assembly process there should be a large contribution to the entropy from the solvent, -T∆S ) -21 kcal/mol or ∆S ) 70 cal/mol/K. The enthalpic contribution from the solvent should be positive. These calculations provide valuable information regarding the 1:1 complexes formed by cavitand 1. In the formation of capsular complexes formed by 1, for example, 2:2 host-guest complexes, this is only the first step toward the final assembly. The capular complex would have different thermodynamics from the 1:1 complex because it eliminates the water in contact with the guest. The formation of this complex would make for some interesting future studies. Acknowledgment. S.W.R. gratefully acknowledges the financial support of the National Science Foundation (CHE0611679), and B.C.G. gratefully acknowledges the financial support of the National Science Foundation (CHE-0718461). References and Notes (1) Gibb, B. G. In Organic Nano-Structures; Atwood, J. L., Steed, J. W., Ed.; John Wiley and Sons: 2007. (2) Gibb, C. L. D.; Gibb, B. C. J. Am. Chem. Soc. 2006, 128, 16498. (3) Sun, H.; Gibb, C. L. D.; Gibb, B. C. Supramolecular Chem 2008, 20, 141.
(4) Gibb, C. L. D.; Gibb, B. C. J. Am. Chem. Soc. 2004, 126, 11408. (5) Gibb, C. L. D.; Sundaresan, A. K.; Ramamurthy, V.; Gibb, B. C. J. Am. Chem. Soc. 2008, 130, 4069. (6) Brovchenko, I.; Paschek, D.; Geiger, A. J. Chem. Phys. 2000, 113, 5026. (7) Wallqvist, A.; Gallicchio, E.; Levy, R. M. J. Phys. Chem. B 2001, 105, 6745. (8) Chau, P. L. Mol. Phys. 2001, 99, 1289. (9) Setny, P.; Geller, M. J. Chem. Phys. 2006, 125, 144717. (10) Paschek, D. J. Chem. Phys. 2004, 120, 10605. (11) Lee, C. Y.; McCammon, J. A.; Rossky, P. J. J. Chem. Phys. 1984, 80, 4448. (12) Chau, P. L.; Forester, T. R.; Smith, W. Mol. Phys. 1996, 89, 1033. (13) Cheng, Y.; Rossky, P. J. Nature 1998, 392, 696. (14) Lum, K.; Chandler, D.; Weeks, J. D. J. Phys. Chem. B 1999, 103, 4570. (15) Southall, N. T.; Dill, K. A. J. Chem. Phys. B 2000, 104, 1326. (16) Chau, P. L. Mol. Phys. 2003, 101, 3121. (17) Laaksonen, A.; Stilbs, P. Mol. Phys. 1991, 74, 747. (18) Ludwig, R.; Weinhold, F.; Farrar, T. C. J. Chem. Phys. 1997, 107, 499. (19) Haselmeier, R.; Holz, M.; Marbach, W.; Weingartner, H. J. Phys. Chem. 1995, 99, 2243. (20) Beckstein, O.; Sansom, M. S. Proc. Natl. Acad. Sci. U.S.A. 2003, 100, 7063. (21) Hummer, G.; Rasaiah, J.; Noworyta, J. Nature 2001, 414, 188. (22) Wagne, A.; Rasaiah, J.; Hummer, G. J. Chem. Phys. 2002, 117, 10789. (23) Reichman, S. A. D.; Hummer, G. J. Chem. Phys. 2005, 123, 194502. (24) ten Wolde, P. R.; Chandler, D. Proc. Natl. Acad. Sci. USA 2002, 99, 6539. (25) Liu, P.; Huang, X.; Zhou, R.; Berne, B. J. Nature 2005, 437, 159. (26) Zhou, R.; Huang, X.; Margulis, C. J.; Berne, B. J. Science 2004, 305, 1605. (27) Miller, III, T. F.; Vanden-Eijnden, E.; Chandler, D. Proc. Natl. Acad. Sci. USA 2007, 104, 14559. (28) Levy, Y.; Onuchic, J. N. Annu. ReV. Biophys. Struc. 2006, 35, 389. (29) Setny, P. J. Chem. Phys. 2007, 127, 054505. (30) Case, D. A.; Pearlman, D. A.; Caldwell, J. W.; Cheatham, T. E.,. III; Wang, J.; Ross, W. S.; Simmerling, C. L.; Darden, T. A.; Merz, K. M.; Stanton, R. V.; Cheng, A. L.; Vincent, J. J.; Crowlet, M.; Tsui, V.; Gohlke, H.; Radmer, R. J.; Duan, Y.; Pitera, J.; Massova, I.; Seibel, G. L.; Singh, U. C.; Weiner, P. K.; Kollman, P. A. AMBER; University of California, San Francisco; San Francisco, 2002. (31) Horn, H. W.; Swope, W. C.; Pitera, J. W.; Madera, J. D.; Dick, T. J.; Hura, G. L.; Head-Gordon, T. J. Chem. Phys. 2004, 120, 9665. (32) Krouskop, P. E.; Madura, J. D.; Paschek, D.; Krukau, A. J. Chem. Phys. 2006, 124, 016102. (33) Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Montgomery, J. A., Jr.; Vreven, T.; Kudin, K. N.; Burant, J. C.; Millam, J. M.; Iyengar, S. S.; Tomasi, J.; Barone, V.; Mennucci, B.; Cossi, M.; Scalmani, G.; Rega, N.; Petersson, G. A.; Nakatsuji, H.; Hada, M.; Ehara, M.; Toyota, K.; Fukuda, R.; Hasegawa, J.; Ishida, M.; Nakajima, T.; Honda, Y.; Kitao, O.; Nakai, H.; Klene, M.; Li, X.; Knox, J. E.; Hratchian, H. P.; Cross, J. B.; Bakken, V.; Adamo, C.; Jaramillo, J.; Gomperts, R.; Stratmann, R. E.; Yazyev, O.; Austin, A. J.; Cammi, R.; Pomelli, C.; Ochterski, J. W.; Ayala, P. Y.; Morokuma, K.; Voth, G. A.; Salvador, P.; Dannenberg, J. J.; Zakrzewski, V. G.; Dapprich, S.; Daniels, A. D.; Strain, M. C.; Farkas, O.; Malick, D. K.; Rabuck, A. D.; Raghavachari, K.; Foresman, J. B.; Ortiz, J. V.; Cui, Q.; Baboul, A. G.; Clifford, S.; Cioslowski, J.; Stefanov, B. B.; Liu, G.; Liashenko, A.; Piskorz, P.; Komaromi, I.; Martin, R. L.; Fox, D. J.; Keith, T.; Al-Laham, M. A.; Peng, C. Y.; Nanayakkara, A.; Challacombe, M.; Gill, P. M. W.; Johnson, B.; Chen, W.; Wong, M. W.; Gonzalez, C.; Pople, J. A. Gaussian 03, revision C.02; Gaussian, Inc.: Wallingford, CT, 2004. (34) Wang, J.; Wolf, R. M.; Caldwell, J. M.; Kollman, P. A.; Case, D. A. J. Comput. Chem. 2004, 25, 1157. (35) Wang, J.; Cieplak, P.; Kollman, P. A. J. Comput. Chem. 2000, 21, 1049. (36) Preusser, A. ACM Trans. Math. Software 1989, 15, 79. (37) Mancera, R. L.; Buckingham, A. D. J. Chem. Phys. 1995, 99, 14632. (38) Lee, B.; Richards, F. M. J. Mol. Biol. 1971, 55, 379. (39) Vaitheeswaran, S.; Yin, H.; Rasaiah, J. C.; Hummer, G.; Sengers, J. Proc. Natl. Acad. Sci. U.S.A. 2004, 101, 17002. (40) Jana, B.; Pal, S.; Maiti, P. K.; Lin, S. T.; Hynes, J. T.; Bagchi, B. J. Phys. Chem. B 2006, 110, 19611. (41) Adam, G.; Gibbs, J. H. J. Chem. Phys. 1965, 43, 139. (42) Dzugutov, M. Nature 1996, 381, 137. (43) Lin, S.; Blanco, M.; Goddard, W. A., III J. Chem. Phys. 2003, 119, 11792. (44) Ben-Naim, A.; Marcus, Y. J. Chem. Phys. 1984, 81, 2016.
JP804429N
