Hydration Interaction between Phospholipid Membranes: Insight into

Article pubs.acs.org/Langmuir

Hydration Interaction between Phospholipid Membranes: Insight into Different Measurement Ensembles from Atomistic Molecular Dynamics Simulations Matej Kanduč,*,† Emanuel Schneck,§ and Roland R. Netz† †

Department of Physics, Free University Berlin, D-14195 Berlin, Germany Department of Theoretical Physics, J. Stefan Institute, SI-1000 Ljubljana, Slovenia § Institut Laue-Langevin, 6 Rue Jules Horowitz, F-38042 Grenoble, France ‡

ABSTRACT: Using the novel thermodynamic extrapolation technique in molecular dynamics simulations, we investigate the interaction between phospholipid bilayers subject to various boundary conditions that correspond to established experimental methods for the determination of pressure− distance curves: the osmotic stress method, the hydrostatic method, and the surface force apparatus method. We discuss the roles of van der Waals and Helfrich undulation pressures in the force balance and find that they do not play a major role in the distance range below 28 water molecules per lipid as considered by us. We address the influence of experimental boundary conditions on bilayer structural changes as well as the consequences on interaction pressures. Significant discrepancies are observed between pressures obtained in osmotic stress and hydration methods on one hand and the surface force apparatus method on the other hand. We quantify the contribution of lipid volume compressibility to the total work of dehydration and find it to be substantial for high pressures. In a wide hydration range, the interaction pressure is mostly determined by the area per lipid molecule. This means that the influence of fatty acid chemistry on experimental pressure−distance curves is indirect and mediated by the area per lipid.

■

INTRODUCTION Experimental studies over the last four decades revealed a strong repulsion acting between polar surfaces in water.1,2 This force is very short-ranged, acting typically on the nanometer scale, and roughly obeys an exponential-decay law. Measured decay lengths range from 0.08 to 0.64 nm,3,4 whereas the repulsive pressure can easily reach several kilobars at small separations and thus dominates over electric and dispersion interactions, which are jointly referred to as the Derjaguin− Landau−Verwey−Overbeek (DLVO) interactions. In accord with standard notation, this repulsive force is called hydration repulsion and refers to the total interaction that acts between polar surfaces at short separation. It is particularly important in regulating biological processes such as membrane fusion, interaction of proteins with surfaces, proteins aggregation, and bacterial adhesion.5−7 From the 1970s on, the hydration repulsion between phospholipid bilayers, which are models of biological membranes,1,3,8−13 has been extensively studied using various measurement techniques. A straightforward (albeit technically demanding) way of measuring pressures between lipid membranes is using the surface force apparatus (SFA) established by Israelachvili and co-workers.14 Here, the membranes are immobilized on two cylindrical surfaces and brought into contact;15,16 see Figure 1c. The force acting between the cylinders can be translated into an interacting pressure between the planar surfaces via the Derjaguin approximation.17 In a different experimental approach, pressure−distance curves are obtained by X-ray or neutron © 2013 American Chemical Society

diffraction from membrane multilayers. Here, the lamellar periodicity Lz is determined from the diffraction angles, while the multilayers are compressed by exerting hydrostatic pressures (Figure 1b) or osmotic stress.1 The osmotic stress acting on a membrane stack is conveniently expressed by the equivalent pressure. Low equivalent pressures (100 bar) are reached when the multilayers are no longer in physical contact with bulk water but interact via vapor exchange with a water reservoir with reduced vapor pressure (e.g., saturated salt solutions18), Figure 1a. The water layer thickness Dw is worked out from the measured lamellar periodicity, either gravimetrically19 or structurally11 from scattering length density profiles. Although the various experimental approaches for determining pressure− distance curves have been routinely used in a complementary manner, they are not strictly equivalent. Membranes are subject to different boundary conditions: in an SFA experiment, the area per lipid is kept approximately constant, while at the same time a mechanical force perpendicular to the membrane surfaces is exerted. In osmotic stress measurements at Received: April 2, 2013 Revised: June 17, 2013 Published: June 26, 2013 9126

dx.doi.org/10.1021/la401147b | Langmuir 2013, 29, 9126−9137

Langmuir

Article

Figure 1. Schematics of the three methods for measuring the interaction pressure between lipid bilayers. (a) In the osmotic stress method, the lipid is in equilibrium with water vapor of known chemical potential μ. (b) In the hydrostatic method, the lipid is squeezed under pressure p in a cell from which water can escape through a semipermeable membrane. (c) In the SFA method lipids are grafted on two crossed cylindrical surfaces and brought into contact under external force F at fixed water chemical potential as determined by the surrounding bulk solvent.

Figure 2. (a) Structural parameters of the lipid bilayers. (b) Schematic drawing of an amphiphilic DLPC molecule with zwitterionic phosphatidylcholine head group and saturated C12 alkyl chains. (c) Simulation snapshot of a phospholipid bilayer composed of Nlip = 72 DLPC lipids and N = 2016 water molecules (nw = 28). The simulation box is shown as a red rectangle and is repeated infinitely many times in all three directions via periodic boundary conditions. Water is shown only in the red rectangle for clarity. (d) Density profiles of water and DLPC bilayer at 1 bar and full hydration nw = 28 (solid lines), dehydrated to nw = 4 in the osmotic stress scenario at p′ = 1 bar (dashed lines), and dehydrated to nw = 4 in the hydrostatic scenario at bulk-water chemical potential, that is, Δμ = 0, which corresponds to about p′ = 1600 bar (dotted lines).

of water when the membranes are hydrated with varying numbers of water molecules.25 Simulating finite-sized membrane surfaces in explicit bulk-water reservoirs26 is very timeconsuming and does not allow for the description of extended lipid bilayers with realistic behavior. Alternatively, a constant chemical potential can be assured through stochastic deletions and insertions of water molecules,27 in a method termed grand canonical Monte Carlo simulation (GCMC). Here, the bottleneck is the pressure resolution as determined by the precision with which the chemical potential of water can be controlled.28 In the present work, we employ a novel approach that we have introduced recently, termed Thermodynamic Extrapolation (TE) method.23 The method enables us to perform atomistic simulations of interacting lipid bilayers at prescribed chemical potential under various boundary conditions, and allows to infer the interacting pressure with high accuracy. This puts us in the position to scrutinize in detail the role of membrane reorganizations under various boundary conditions. We address three sets of boundary conditions, corresponding to the three above-mentioned experimental techniques, that is, (i) osmotic stress method, (ii) hydrostatic

atmospheric pressure conditions, on the other hand, only the chemical potential of water is prescribed by the reservoir, and membranes are free to rearrange their molecular organization upon dehydration. Finally, in experiments exerting hydrostatic pressures, the membranes experience isotropic compression as they get dehydrated. So far, the implications of these distinct differences in the boundary conditions have never been studied theoretically on a microscopic level. Hydration repulsion has also been the subject of many theoretical studies. But despite enormous efforts, its theoretical description turned out to be difficult.2 The orientation of water molecules20 as well as the configurational entropy of lipid molecules21 have been discussed as possible sources for the repulsion. It has also been suggested that hydration repulsion exhibits various regimes that depend on separations.22,23 It was concluded that, in order to explain the hydration interaction on a quantitative level, the structure of the solvent has to be taken into account explicitly.3 In recent years, this insight has drawn the attention toward the application of atomistic simulations, where water molecules are modeled explicitly including all relevant degrees of freedom.24 The main challenge in today’s atomistic simulations is to provide constant chemical potential 9127


Langmuir

Article

Table 1. Summary of the Different MD Simulation Setupsa

method, and (iii) SFA measurements, and examine their influence on the interacting pressure.

■

METHODS

Computer Model of Interacting Lipid Membranes. We simulate a single periodically replicated hydrated planar lipid bilayer in the liquid-crystalline Lα phase. A snapshot of the computer model is shown in Figure 2c. The bilayer is composed of dilinoleoylphosphatidyl-choline (DLPC), a phospholipid which consists of two 12-carbon atom long aliphatic chains attached to a zwitter-ionic head group. Experimentally, PC-lipids with saturated chains exhibit a transition into the gel-like L′β phase upon dehydration.29 However, the phase transition from Lα to Lβ′ phase is not directly accessible by atomistic MD simulations.30,31 In order to be able to compare the simulations to experimental data on liquid-crystalline membranes in a wide hydration range, we chose DLPC which assumes the Lα phase down to very low hydration levels due to its rather short chains. Since the simulated DLPC does not exhibit the transition into L′β phase on the simulation time scales at all, while the experimental system does, the analysis and interpretation of the simulation data at low hydration must keep in mind that this state is in fact metastable. The simulation box with dimensions Lx, Ly, and Lz contains one lipid bilayer composed of Nlip = 72 DLPC molecules (36 in each layer) hydrated with various numbers of water molecules N. The level of hydration is defined as the number of water molecules per lipid molecule, nw = N/Nlip. The membrane is arranged parallel to the xyplane and by periodic boundary conditions in all three directions a stack of lipid bilayers is represented. In order to mimic the three experimental situations described in the introduction (osmotic stress method, hydrostatic method, and SFA method), we carry out simulations at T = 300 K subject to three different sets of boundary conditions. Simulations in the Np′T ensemble, where the simulation pressure p′ = 1 bar is isotropic, correspond to the osmotic stress method and are referred to as osmotic stress simulations in the following. Under these boundary conditions, the water chemical potential in the simulations, μ′, is the main thermodynamic control parameter. Hydrostatic measurements are also represented by simulations in the NpT ensemble, but at a hydration-dependent isotropic pressure p. Here, in analogy with the experimental setup, p is the interaction pressure and is determined by the reference water chemical potential. These simulations are referred to as hydrostatic simulations. In both osmotic stress and hydrostatic simulations the pressure is isotropic, enforced via independent box scaling in lateral and perpendicular directions, respectively. In the SFA method, the lipids are grafted to the substrate and the surface area of lipids is approximately constant upon changing the surface separations.14,15 Note that it is difficult to estimate how much the area per lipid can vary in an actual experiment; for convenience, we therefore mimic an ideal SFA with fixed area per lipid in our simulations, but also discuss the effect of deviations from this ideal behavior. In the simulations, the ideal SFA is realized using the NApzT ensemble, where the area per lipid is kept at 0.65 nm2, which corresponds to the value obtained in Np′T simulations with p′ = 1 bar at full hydration (nw = 28), and the perpendicular pressure pz is controlled through adjustments in the z-extension of the simulation box, Lz. These simulations are referred to as SFA simulations in the following. See Table 1 for short overview of all three ensembles. In each of the three simulation types, interaction pressures at the reference chemical potential of water, μ, are calculated via thermodynamic extrapolation (see next section). As the membrane interfaces are not perfectly sharp, there is no unique way of defining the surface−surface separation. In accordance with experimental practice,1,9 we define the separation as the thickness of an “equivalent water layer″ with bulk density, that is,

vN Dw = w A

actual pressure

actual chemical potential

T (K)

MD type

ensemble

DLPC (osmotic stress) DLPC (hydrostatic)

Np′T

p′ = 1 bar

μ′(nw)

300

NpT

μ′ = μ

300

DLPC (SFA), first step DLPC (SFA), second step

NAp′T

p(nw) from TE p′ = 1 bar

μ′ = μ(nw)

300

μ′ = μ

300

NApT

p(nw) from TE

The first MD setup represents the osmotic stress method and is done at p′ = 1 bar; therefore, the actual chemical potential μ′ depends on the hydration level nw. The second simulation setup corresponds to the hydrostatic method and is performed at a pressure p that is obtained via eq 3 from the osmotic stress MD. Consequently, the chemical potential μ′ matches the bulk-water chemical potential μ at all hydration levels. The third setup mimics SFA measurements where the area per lipid is kept constant and the pressure is controlled in perpendicular direction. The SFA simulations are done in two subsequent steps: the first at constant actual pressure p′ = 1 bar, and the second at extrapolated pressures p(nw) determined by the bulkwater chemical potential. a

water molecules per lipid, the water layer thickness Dw depends on the area per lipid and thus on the boundary conditions. We present our results as a function of Dw rather than of nw, consistent with most experimental studies. For more technical details about the simulations see the Appendix. From the known box size Lz, which corresponds to the bilayer repeat distance, and the water-slab thickness Dw, we define the bilayer thickness as

Dlip = Lz − Dw

(2)

Interaction Pressures at Prescribed Chemical Potential: Thermodynamic Extrapolation. Molecular dynamics simulations of biological systems are usually performed in either the canonical constant-volume (NVT) ensemble or the canonical constant-pressure (NPT) ensemble. The relevant ensemble for the study of phospholipid bilayers interacting in an aqueous environment is the grand canonical ensemble, where water molecules are allowed to exchange with an external reservoir with prescribed chemical potential μ. Here, we use the recently introduced Thermodynamic Extrapolation (TE) method,23 which provides fast numerical convergence and, thus, high pressure resolution. The idea of this method is the following: simulations are performed in the canonical ensemble, that is, with a constant number of water molecules N, at a chemical potential μ′ that typically deviates significantly from the reference potential μ. A thermodynamic relation is then used to deduce the interaction pressure p from the simulation pressure p′ and the deviation in the chemical potential μ − μ′, which is precisely measured as described in the Appendix,

p ≃ p′ +

μ − μ′ vw

(3)

Here, vw is the volume occupied by a single water molecule in bulk. In our case of the SPC/E water model, the water volume is vw = 0.030 nm3 at a pressure of 1 bar. For a detailed derivation see the Appendix, where we show that the linear relationship between p and μ in eq 3 works fairly well up to the kilobar range as water is quite incompressible and therefore vw remains almost constant.

■

RESULTS AND DISCUSSION Influence of Boundary Conditions on DLPC Bilayer Structure. Figure 2d shows simulated mass density profiles of DLPC membranes and water. At a hydration level of nw = 28 water molecules per lipid (solid lines), which corresponds to a separation of around Dw = 2.7 nm, the water density reaches

(1)

where N is the number of water molecules in the simulation box with lateral surface area A = LxLy, and vw = 0.030 nm3 is the volume per water molecule in bulk at atmospheric pressure. For a given number of 9128


Langmuir

Article

Figure 3. Structural properties of DLPC bilayers from the three simulations methods (osmotic stress simulations, hydrostatic simulations, and the second step SFA simulations) compared to data from osmotic stress experiments9 as a function of membrane distance: (a) area per lipid Alip, (b) bilayer thickness Dlip, and (c) volume per lipid (1/2)DlipAlip.

Influence of Boundary Conditions on DLPC Bilayer Interaction. In the next step, we investigate to what extent the structural parameters, which significantly depend on the set of boundary conditions, influence the interaction pressure between DLPC bilayers. For all simulation types (osmotic stress simulations, hydrostatic simulations, and SFA simulations), interaction pressures are determined by Thermodynamic Extrapolation (see above). For SFA simulations, this is done in simulations with atmospheric pressure conditions, p′ = 1 bar. Figure 4 shows pressure−distance curves of various phosphatidylcholine membranes in the liquid-crystalline Lα

the bulk value in the center. At this high hydration level, the interaction pressure is negligible, so that all three sets of boundary conditions fully coincide (not shown). When the hydration level is decreased to nw = 4 in osmotic stress simulations (i.e., while keeping the pressure at p′ = 1 bar), the separation decreases to Dw = 0.4 nm and the chemical potential drops by Δμ = −2.8 kJ/mol. This significantly affects the bilayer profile (see dashed line). The membrane thickness Dlip increases from 3.1 to 3.6 nm, while the area per lipid Alip = 2A/ Nlip decreases from 0.65 down to 0.55 nm2, which is quantified in Figure 3. However, the volume per lipid, defined as Vlip = (1/2)DlipAlip, is nearly unaffected in the osmotic stress scenario. The same trend is also observed experimentally9 and can be interpreted in the following way: as the hydration level drops, the head groups fill up the voids, resulting in a reduced area per lipid. When the hydration level is decreased to nw = 4 in hydrostatic simulations (i.e., while keeping the chemical potential constant), the pressure rises to 1600 bar. The changes in membrane thickness and area per lipid are qualitatively similar to those in osmotic stress simulations, but the change in area is more pronounced while the change in thickness is less pronounced, Figure 3a and b. The overall profile shape (see dotted line in Figure 2d) does not differ much from that found in osmotic stress simulations at the same hydration, although the membrane density is significantly increased, corresponding to a reduction in the volume per lipid from 1 down to around 0.92 nm3. Surprisingly, the density minimum in the middle of the bilayer, known as “methyl dip″, is still very pronounced at 1600 bars. In contrast, when the hydration level is decreased in SFA simulations (i.e., while keeping the area per lipid Alip constant), only the membrane thickness can decrease. The corresponding decrease in the volume per lipid at high pressures exerted normal to the membrane plane (green diamond symbols in Figure 3c) is similar to the hydrostatic case (red squares). We conclude that the area per lipid mostly depends on the level of hydration (except for SFA simulations, where Alip is fixed), while it is less sensitive to the exerted pressure. The volume per lipid, on the other hand, basically only depends on the exerted pressure, whereas the membrane thickness is sensitive to both hydration level and exerted pressure. This difference is confirmed by the comparison with experimental data obtained by the osmotic stress method,9 black circles in Figure 3, which reasonably agree with the osmotic stress simulations and exhibit distinct deviations from the hydrostatic and SFA simulations, particularly for the volume per lipid Vlip.

Figure 4. Pressure−distance plots evaluated by TE from all simulation sets and compared to experimental values. The abbreviations in the experimental references stand for the used techniques: osm = osmotic stress method, hyd = hydrostatic method, and SFA = surface force apparatus. The repeat distances in data from refs 12 and 32 were converted into separations Dw by using the relationship in ref 9. The blue dashed line is the pressure derived from a fit of the free energy in the SFA experiment.10

phase obtained experimentally by various methods in a semilogarithmic plot. We notice significant differences between the different experimental curves from different sources, even when the same lipid type was used. This comparison provides us with a rough estimate of the (systematic) uncertainties inherent to the experimental data. All experimental data, however, exhibit roughly exponential distance dependence with decay lengths in a range from λ = 0.25 to 0.42 nm.1,9,10,12,13,32 The plot also contains pressure−distance curves from osmotic stress simulations, hydrostatic simulations, and SFA simulations of DLPC bilayers, together with SFA-type simulations of DPPC bilayers at 320 K from ref 23. Note that the simulation data extend to very low hydration levels, which is possible, since the 9129


Langmuir

Article

Figure 5. (a) Difference in extrapolated pressures between DLPC bilayers in osmotic stress and hydrostatic simulations. (b) Difference in extrapolated pressures between the SFA and osmotic-stress simulations.

experience various levels of lateral mobility,33−35 which influences the area per lipid at close-contact point of crossed cylinders in SFA experiments.36 The extreme case of lipids that are freely mobile in lateral directions would correspond to a constant chemical potential of lipid molecules, where the surrounding lipid layer would serve as a reservoir. In principle, such a scenario can be theoretically realized, but is beyond the scope of our current work. We remark that, under such circumstances, the further increased area per lipid may result in even more pronounced deviations in the interaction pressure when compared to osmotic stress or hydrostatic experiments. We conclude that great care is required when comparing SFA and multilayer experiments on a quantitative level. We further note that, in SFA experiments, the setup is composed of only two interacting lipid layers, usually on mica support, whereas in our simulations they are arranged on a periodic lattice. The periodicity gives rise to additional interactions between head-group dipoles. However, since the dipoles interact across a bilayer with a low dielectric constant, the corresponding dipolar interactions are reduced when compared to the case of dipoles interacting through a water slab.37 Therefore, these interactions are negligible compared to the dominant mechanisms at play between the surfaces of adjacent membranes. Additional van der Waals contribution due to the mica supports was also considered, but it has significant effect only at larger separations and hence influences the bound-state distances, but not the interaction at small distances.38 Solid supports in the SFA experiments also suppress undulation modes and therefore inhibit undulation interactions.39 As we show in this work, this effect does not change the hydration interactions at the studied separations range. Another surprising observation from comparing SFA with the other simulation scenarios is that the larger area per lipid leads to slightly larger interaction pressure, Figure 5b. This means that the total force between the membranes cannot be considered as the sum of independent head-group contributions. Instead, the lateral distance between neighboring head groups decisively influences the individual head group’s contribution. Upon dehydration hydrogen bonds between lipids and water are broken, see the Appendix. The number of hydrogen bonds per surface area tend to decrease slightly more in the SFA scenario than in the osmotic-stress and hydrostatic scenarios. This result is consistent with the slightly stronger hydration repulsion observed in the SFA simulations, see Figure 5b.

simulations do not capture the dehydration-induced phase transition into the L′β phase (as explained in the Methods section). All simulated curves are relatively similar but tend to exhibit stronger repulsion at large distances than do the experimental curves. Fitting the simulated curves to the singleexponential function p = p0 e−Dw/λ yields decay lengths in the range λ = 0.42−0.47 nm; see Appendix for more details. The systematic discrepancies between experimental and simulated data can be attributed to force field limitations, as discussed further below. In the following, we focus on analyzing differences between our simulations under different sets of boundary conditions. Figure 5a shows the difference between interaction pressures in osmotic stress and hydrostatic simulations. The difference is not significant within the error. This indicates that the membrane thickness, which differs significantly between osmotic stress and hydrostatic simulations, has only little influence on the interaction pressure, as long as the area per lipid does not deviate much (see Figure 3a). It should be noted that, from a methodological viewpoint, the pressure−distance curves in the hydrostatic scenario are obtained by two consecutive simulations: the first at p′ = 1 bar, leading to the extrapolated pressure p via thermodynamic extrapolation, eq 3, and the second simulation at the obtained pressure p in the NpT ensemble. In contrast to the first iteration, represented by the pressure−distance curves obtained in the osmotic stress simulations, they account for the structural response of the simulated system to the exerted pressures. It is remarkable that TE yields accurate pressure−distance curves already in the first iteration despite the liquid-crystalline nature of the simulated system, which makes its structure very responsive to pressure. In other words, the osmotic stress and hydrostatic scenarios are equivalent as far as pressure−distance curves are concerned. Figure 5b shows the difference between the interaction pressures in the osmotic stress and SFA simulations. Below a membrane distance of about Dw = 1 nm, SFA simulations exhibit significantly stronger repulsion (up to 20%) than the osmotic stress simulations. The difference can be attributed to the substantial deviations in the area per lipid below about Dw = 1 nm, when comparing the osmotic stress and SFA simulations (see Figure 3a). In our SFA simulations, the area per lipid is fixed, whereas in real SFA experimentsdepending on the grafting method of lipid to substratethe area per lipid may vary under applied pressure. It has been shown in several studies thatdepending on the particular substrate type, temperature, relative humidity, preparation procedure, and so forththe deposited lipids 9130


Langmuir

Article

The Role of van der Waals and Helfrich Undulation Contributions. We now discuss the contribution of van der Waals (vdW) and undulation forces to the interaction pressure. Dispersion interactions are included implicitly in the simulations via Lennard-Jones interactions between the particles, but they are treated with a certain cutoff. Since the used force field in our simulations is ad hoc designed in order to reproduce the correct area per lipid in the liquid-crystalline phase,40 it does not necessarily reproduce the correct longrange vdW interactions. In the standard continuum approach, the van der Waals pressure acting between two parallel semiinfinite slabs separated by distance D is given as p = −H/6πD3, where H denotes the Hamaker constant.41 We estimate the corresponding Hamaker constant in the simulations to be about Hsim = 30 × 10−21 J, as described in the Appendix, whereas experimentally determined Hamaker constants for phospholipid bilayers interacting across water are in the range Hexp = 5 − 15 × 10−21 J.1,9,41 We note that the systematic error associated with this experimental estimate is difficult to assess due to the indirect nature of the method that is used for determination.9,41 Fluid membranes are not rigid plates and therefore subject to thermal excitations of bending modes, which leads to repulsive steric interactions. These Helfrich undulation forces arise due to the increasing steric confinement of undulation modes as two membranes approach each other.42 The corresponding interacting pressure is given as p = 2c(kBT)2/κbD3 at small separations, where c = 0.116 for a stack of only sterically interacting membranes.43 At large separations, the power-law turns into an exponentially decaying function, which happens at a crossover scale of around Dw ∼ 4−1000 nm,43 meaning that our system is not affected by this exponential crossover. Therefore, the undulation pressure has the same scaling as vdW, and can be rewritten as p = HHelf/6πD3, where HHelf = 12πc(kBT)2/κb. Considering the bending rigidity of DLPC to be κb = 10 kBT,44 we get an effective Helfrich Hamaker constant of HHelf = 2 × 10−21 J, which is much smaller than the van der Waals Hamaker constant. In addition, we note that for membranes that are bound by attractive forces, the amplitude of the Helfrich repulsion is expected to be reduced.45 In essence, the fact that the Helfrich undulation repulsion42 is not well represented in the simulations, since the simulation box is too small to account for the relevant undulation modes, is of no concern because of the negligible undulation pressure amplitude. In order to quantitatively assess the importance of the combined vdW and Helfrich contributions, we first analytically subtract the intrinsic vdW interactions from the obtained pressure−distance curves. For this subtraction, we use eq 4 in ref 46, which accounts for dispersion interactions with a cutoff as used in MD simulations. We then add various analytical van der Waals pressures corresponding to various effective Hamaker constants H, Δp = −

H 6πDw 3

Figure 6. Extrapolated pressure−distance curves from osmotic stress simulations (blue data points) and the corresponding singleexponential fit (blue solid line). In the modified curves, the intrinsic vdW contribution is subtracted and modified vdW interactions proportional to various Hamaker constants H are added.

contribution eq 4 with different effective Hamaker constants. The correction remains small compared to the total interaction even when H is scaled by an order of magnitude (up to 300 × 10−21 J), leading to the conclusion that van der Waals and undulation contributions can be neglected in our studied distance range. First, this means that a possibly inaccurate representation of the vdW contribution (due to inaccurate force field or used cutoff) in the employed computer model has no impact on our conclusions. Second, this finding also explains why pressure−distance curves are vastly independent of the membrane thickness: While the hydration effects are generated at the membrane surfaces, the membrane thickness only influences vdW and Helfrich interactions, but the latter two are negligible. We remark, however, that vdW and Helfrich contributions become important at larger distances, where they significantly influence the equilibrium separation of phospholipid bilayers at full hydration.9 In future studies on more strongly hydrated membranes, the treatment of vdW and Helfrich interactions thus needs to be suitably revised and improved. Based on these results, we further conclude that the significant deviation between experimental pressure−distance curves and our simulation results (see Figure 4) is not due to an improper representation of vdW or Helfrich interactions, but most likely due to force field limitations in the correct representation of head-group hydration effects. Here, a careful revision of the lipid force fields will be crucial for future studies. Influence of Lipid Chain Length. As seen in Figure 4, there is no significant difference between the pressure−distance curves obtained in SFA simulations of DLPC and DPPC. The force field representation of the two molecules is identical, apart from the difference in chain lengths. Each DLPC molecule has two alkyl tails, consisting of 12 carbon atoms, whereas the tails in DPPC consist of 16 carbon atoms. Since the area per lipid is dictated by the boundary condition (Alip = 0.65 nm2), the membrane thickness is the only important structural difference between the two systems. We conclude that the lipid chain length has no direct influence on the interaction pressure between phospholipid membranes unless it has an impact on the properties of the membrane surface (e.g., on the area per lipid, or on the head group configurational freedom). In other words, the influence of the chain length on the interaction pressure is only an indirect one. Response of Phospholipid Membranes to Applied Pressure. In the last part, we discuss the structural response of lipid bilayers to pressure on a thermodynamic level. Our results

(4)

The tunable parameter H accounts for the effects of the sum of attractive van der Waals and repulsive undulation interactions. Figure 6 shows the pressure−distance data obtained in osmotic stress simulations of DLPC bilayers (data points) with the best exponential fit (solid line) prior to the subtraction of the intrinsic vdW contribution. The various dashed lines correspond to the fit with subtracted intrinsic dispersion contribution and added infinite-ranged (without a cutoff) vdW 9131


Langmuir

Article

Figure 7. (a) Normalized volume per lipid as a function of applied pressure on the membrane stack at constant water chemical potential. Solid lines are fitted curves. The normalized volumes are compared to bulk SPC/E water normalized volume vw/v0w (blue circles). (b) The corresponding compressibilities evaluated from the volume−pressure fits via eq 5. Here and only here the pressure-dependent water volume vw(p) is used for calculation of Dw and Vlip via eq 1. (c) Fractions of incremental free energy change going into bilayer shear and longitudinal compression deformations ηshear and ηcomp as defined in eqs 8 and 9.

From the known structural quantities and their distance dependence in Figure 3, we can decompose the free energy change into various contributions. We consider the free energy - per lipid molecule as a function of three independent parameters. One variable describes the bilayer separation Dw, the other two correspond to bilayer molecular packing. As packing parameters we choose the area per lipid Alip and the volume per lipid Vlip. Other quantities such as the bilayer thickness Dlip and hydration level nw can be evaluated from the three chosen variables. The differential of the total free energy per lipid - (Dw,Alip,Vlip) at given temperature and chemical potential can be evaluated via mechanical considerations (see the Appendix for more details) and is in the isotropic pressure ensemble given by

show that the lipid bilayer undergoes significant deformation and volume change in the different simulation ensembles. In order to understand the thermodynamic consequences of this, we derive expressions for different contributions to the bilayer free energy change as the water slab thickness varies. Figure 7a shows the normalized volume per lipid as a function of applied pressure in the constant chemical-potential ensembles, that is, in the hydrostatic and SFA scenarios. Note that the lipid volume Vlip = (1/2)DlipAlip reflects the pressure dependence of the water volume vw via the definition of Dw, eq 1. Fitting the volume−pressure curves with the purely empirical function c1 + c2 tanh(c3 log p + c4) allows one to evaluate the lipid compressibility as χlip = −

1 dVlip Vlip dp

(5)

1 1 d- =− pAlip dDw − pDw dAlip − p dVlip 2 2

In both ensembles, the lipid compressibilities are very similar and decay roughly logarithmically with applied pressure, that is, χlip(p) ∝ −log p, Figure 7b. At small pressures the lipid compressibility is around three times larger than the bulk-water compressibility, which is indicated in Figure 7b by a blue dotted line: In this limit, the system compressibility is therefore dominated by the lipid compressibility alone. With increasing pressure, the lipid compressibility drops down and at around p = 103 bar reaches the water compressibility. Due to this fact, we have used in Figure 7a and b (in contrast to all results presented in the other parts of this paper) the pressure-dependent volume per water molecule, as determined in the bulk simulations shown in Figure 7a, to calculate Dw via eq 1.

≡d-sep + dFshear + dFcomp

(6)

The three terms correspond to the work done by separating the bilayers d-sep, pure shear deformation at constant volume d-shear , and compression d-comp at fixed Alip. Under compression at constant chemical potential μ or hydration level nw, the three geometrical variables are not independent. If we choose to keep Dw as control variable, then the area as well as the volume per lipid depend on the lipid separation, hence Alip = Alip(Dw) and Vlip = Vlip(Dw). As a consequence, the total free energy depends only on separation, that is, -(Dw , Alip(Dw ), Vlip(Dw )) ≡ -(Dw ). 9132


Langmuir

Article

function ∂Vlip/∂Alip and ∂Vlip/∂Dlip, plays a significant role in the bilayer deformation energetics.

To compare the different free energy contributions, we define fractions of incremental free energy change that go into bilayer pure shear and longitudinal compression, namely, ηshear

d-shear(Dw ) = , d-(Dw )

ηcomp =

■

CONCLUSIONS We perform DLPC bilayer MD simulations with three different boundary conditions that mimic common experimental setups for the determination of pressure−distance curves: the osmotic stress method, the hydrostatic method, and the surface force apparatus method. While comparing simulations for different corresponding boundary conditions, we find distinct differences in the structural response to the hydration degree, involving bilayer thickness, area per lipid, and volume per lipid. Pressure−distance curves for all sets of boundary conditions are determined using the thermodynamic extrapolation technique. We observe no significant difference in the interaction pressures for osmotic stress and hydrostatic simulations, and our results confirm the equivalence of the two corresponding experimental methods, which has never been proven but routinely assumed. This result also implies that the TE method is robust and enables accurate pressure predictions in one step, even when the simulated system exhibits some structural response to pressure. At short membrane distances, noticeable differences between interaction pressures are observed when osmotic stress and hydrostatic simulations on one side are compared with SFA simulations on the other side. These differences we ascribe to significant differences in the area per lipid. Interestingly, interaction pressures are higher for larger areas per lipid. Hydrogen-bond analysis indicates that this can be attributed to a subtle dependence of the number of lipid/water hydrogen bonds on the area per lipid. Finally, we show that van der Waals and Helfrich contributions play negligible roles for the studied membrane distance range. Our results imply that, in this distance range, the lipid chain length has an indirect influence on the interaction pressure only if it has an impact on the area per lipid. Our findings are illustrated by a decomposition of the free energy change of the bilayer system upon dehydration into the contributions due to bilayer−bilayer separation, lipid layer shearing, and lipid layer compression. In the hydrostatic and SFA ensembles, the compression work makes a sizable hitherto neglected contribution.

d-comp(Dw ) d-(Dw )

(7)

which can be written as ηshear =

ηcomp =

′ (Dw ) Dw Alip ′ (Dw ) + 2V lip ′ (Dw ) Alip + Dw Alip

(8)

′ (Dw ) 2V lip ′ (Dw ) + 2V lip ′ (Dw ) Alip + Dw Alip

(9)

where Alip ′ = dAlip/dDw and Vlip ′ = dVlip/dDw. In order to extract the above fractions from the simulations, we first fit the results for Alip(Dw) and Vlip(Dw) in Figure 3 with the function c1 + c2 tanh(c3Dw + c4). Figure 7c shows the fractions of incremental free energy changes going into pure shear and volume compression defined in eqs 8 and 9 as a function of Dw. In the osmotic stress simulations, where the volume change is minor, the corresponding work of compression accounts for a few percent at maximum at small separations (blue solid line). The percentage of work going into shear deformation exhibits a peak of around 10% at separations around 0.5−1.5 nm (red solid line). Very similar shearing-deformation fraction is seen in hydrostatic simulations (red dashed line). On the other hand, in the hydrostatic method, the compression of membrane is significant and accounts for up to 60% of free energy change at small separations (blue dashed line). In the case of SFA, the area per lipid is fixed and therefore there is no pure shear deformation (red dotted line). The longitudinal compression deformation in the SFA scenario (blue dotted line) exhibits similar fraction of incremental free energy change as in the case of the hydrostatic method. As a main result, we find the compression work to contribute a significant part of the total free energy change as the bilayer separation is changed. To connect to previous formulations1 we also consider the possibility to eliminate only the volume per lipid and assign it to a function Vlip = Vlip(Dw,Alip), so that the free energy -(Dw , Alip) differential reads

■

Simulation Details

∂Vlip ⎞ 1 ⎛ d-(Dw , Alip)= − p⎜Alip + 2 ⎟dDw 2 ⎝ ∂Dw ⎠ ∂Vlip ⎞ 1 ⎛ ⎟dAlip − p⎜⎜Dw + 2 2 ⎝ ∂Alip ⎟⎠

Molecular dynamics simulations are carried out using the GROMACS simulation package.47 We use the Berendsen thermostat with a time constant of 1 ps in order to control the temperature.48 The simulation pressure is controlled by a Berendsen barostat with a time constant of 1 ps and a compressibility parameter 4.5 × 10−5 bar−1. In osmotic stress and hydrostatic simulations, the pressure is applied isotropically but enforced by anisotropic scaling steps of lateral and perpendicular box dimensions. In SFA simulations, only the perpendicular box dimension is scaled, while the lateral box dimensions are fixed. In the simulations, we use a plain Lennard-Jones cutoff of 0.9 nm. Electrostatic interactions are treated using the Particle Mesh Ewald (PME) method49,50 with a 0.9 nm real-space cutoff. We use the SPC/E water model51 and the united-atom Berger force field40 for DLPC. Prior to the production runs, the systems are equilibrated for at least 5 ns. Production runs for each membrane separation have duration of 75 ns. Simulations for thermodynamic integration (used for

(10)

Defining the ratio between the separation and area deformation work in this case yields the generalized expression X =

=

(∂- /∂Dw )Alip dDw (∂- /∂Alip)Dw dAlip dDw /[Dw + 2(∂Vlip/∂Alip)Dw ] dAlip /[Alip + 2(∂Vlip/∂Dw )Alip ]

APPENDIX

(11)

As we show in the present work, the lipid volume compressibility, expressed by the thermodynamic response 9133


Langmuir

Article

⎛ ∂μ ⎞ ⎛ ∂μ ⎞ ⎛ ∂V ⎞ ⎜ ⎟ =⎜ ⎟ ⎜ ⎟ ⎝ ∂p ⎠ N ⎝ ∂V ⎠ N ⎝ ∂p ⎠ N

measurements of the Coulombic part of the water chemical potential μC) have a total duration of 800 ns per data point. Chemical Potential Evaluation in MD

μ′ = kBT ln[ρ(z 0)/ρ0 ] + μC (z 0) + μLJ (z 0)

(17)

The first factor on the right-hand side can be evaluated via the Maxwell relation

In order to perform simulations at prescribed chemical potential, its value has to be determined from the simulations. The chemical potential of water between the bilayers is composed of three contributions,

⎛ ∂ ⎛ ∂F ⎞ ⎞ ⎛ ∂μ ⎞ ⎜ ⎟ =⎜ ⎜ ⎟ ⎟ ⎝ ∂V ⎠ N ⎝ ∂V ⎝ ∂N ⎠V ⎠ N

(12)

⎛ ∂ ⎛ ∂F ⎞ ⎞ ⎛ ∂p ⎞ =⎜ ⎜ ⎟ ⎟ = −⎜ ⎟ ⎝ ∂N ⎠V ⎝ ∂N ⎝ ∂V ⎠ N ⎠V

The first term is the ideal part with ρ(z0) being the water density at position z0 and ρ0 is an arbitrary reference, which yields only an offset. The other two terms correspond to excess Coulomb and Lennard-Jones (LJ) contributions, respectively. In thermal equilibrium, the total chemical potential is homogeneous and thus independent of the position z0, so that we can evaluate it at an arbitrary position. Since in the SPC/E water model the Lennard-Jones interaction is defined by only one site, it is convenient to evaluate μLJ(z0) by the Widom Test Particle Insertion method (TPI).52 The Coulomb part μC(z0) corresponds to a three-site interaction, so that it is more convenient to use the Thermodynamic Integration method (TI).53,54 In order to achieve a pressure accuracy of δp ≃ ± 10 bar at large membrane distances, the chemical potential of water has to be determined with an accuracy of δμ = δp/vw ≃ ±0.02 kJ/ mol, according to eq 3. The bottleneck in this procedure is the TI part, as it converges slowest. The reference chemical potential μ corresponds to the bulk SPC/E water chemical potential at 1 bar and 300 K. We determine its value to be μ = −29.55 ± 0.01 kJ/mol in independent bulk-water simulations, where we used the reference density ρ0 = 1 g/cm3 in eq 12. Note that the choice of ρ0 is irrelevant as it only shifts the chemical potential.

(18)

and the second one via the differential eq 16 for dV = 0, namely, ⎛ ∂V ⎞ ⎛ ∂V ⎞ ⎛ ∂N ⎞ ⎜ ⎟ = −⎜ ⎟ ⎜ ⎟ ⎝ ∂N ⎠ p⎝ ∂p ⎠ ⎝ ∂p ⎠ N V

(19)

This gives the expression ⎛ ∂μ ⎞ ⎛ ∂V ⎞ ⎜ ⎟ = ⎜ ⎟ ≡ vw′ (p) ⎝ ∂p ⎠ N ⎝ ∂N ⎠ p

(20)

where vw′ (p) denotes the partial molecular volume of water between the surfaces at constant pressure. Figure 8 shows the total system size of DLPC simulations at 1 bar (osmotic stress method) as a function of the number of

Thermodynamic Extrapolation

As the system is simulated in the Np′T ensemble, the water chemical potential is a function of the corresponding independent variables, that is, μ′ = μ(N,p′,T). In general the measured chemical potential μ′ deviates from the reference bulk-water chemical potential μ. However, μ′ and the pressure p′ can be used to determine the desired interaction pressure p in the NμT ensemble, p(N , μ , T ) = p(N , μ′, T ) + Δp(N , Δμ , T )

Figure 8. Volume of the simulation box of a system consisting of 72 DLPC lipid molecules at 1 bar as a function of the number of water molecules between the membranes N.

(13)

where we have abbreviated Δμ = μ − μ′. As we deal with constant temperature, we will omit T from further notation. The pressure correction Δp = p − p′ can be expressed as a Taylor series around μ′, 2 ⎛ ∂p ⎞ 1⎛∂ p⎞ Δp = ⎜ ⎟ Δμ + ⎜ 2 ⎟ Δμ2 2 ⎝ ∂μ ⎠ ⎝ ∂μ ⎠ N N

water molecules between the membranes. Due to the low compressibility of water, the trend is almost perfectly linear and yields a partial volume v′w(p) = 0.0305 ± 0.0001 nm3 at 1 bar, in close agreement with the molecular water volume as obtained in independent bulk simulations. Taking the simulated compressibility of water into account, χ = 5 × 10−5 bar−1 (Figure 7b), enables one to predict the water volume at higher pressures, namely, (∂v′w/∂p)N = −χv′w = −1.5 × 10−6 nm3/bar, which yields a volume of vw′ (p) = 0.029 nm3 at 1000 bar, which is a negligible correction for our purposes. Finally, inserting eq 20 into eq 14, the extrapolated pressure at prescribed chemical potential μ reads

(14)

In order to determine (∂p/∂μ)N, we use the following two total differentials, dμ(N , V ) =

⎛ ∂μ ⎞ ⎛ ∂μ ⎞ ⎜ ⎟ dN + ⎜ ⎟ dV ⎝ ∂N ⎠V ⎝ ∂V ⎠ N

⎛ ∂V ⎞ ⎛ ∂V ⎞ d V (N , p) = ⎜ ⎟ d N + ⎜ ⎟ d p ⎝ ∂N ⎠ p ⎝ ∂p ⎠ N

(15)

p(μ , N ) = p(μ′, N ) + (16)

Δμ 1⎛ ∂ 1 ⎞ + ⎜ ⎟ vw′ 2 ⎝ ∂μ vw′ ⎠

Δμ2 N μ′

(21)

The quadratic term in Δμ can be evaluated by considering the water compressibility as discussed above, and it yields

At constant N, the above relations yield 9134


Langmuir

Article

where we consider a fixed water volume vw = 0.030 nm .

between P and N atoms of DLPC, and θ is the angle between this vector and the z-axis. Figure 10a shows cos θ as a function of the membrane distance Dw. Within the given accuracy, all three scenarios exhibit similar behavior with a saturation value of about cos θ = 0.23 at full hydration and a drop towards cos θ = 0 at small distances. This means that head groups align towards the lateral plane at low hydration levels and form an electrostatically favorable two-dimensional lattice. Figure 10b shows the average number of H-bonds per membrane surface area as a function of Dw for all three ensembles. We use the angle-distance criterion as defined in ref 56. In our system, each lipid molecule has 9 acceptor atoms, that is, eight oxygen atoms and one nitrogen atom. The number of H-bonds between lipids and water molecules decreases with decreasing surface separation. In the SFA scenario, the number of H-bonds tends to be slightly smaller at short separations compared to the other two scenarios, which implies that more H-bonds are broken upon dehydration. This appears to be consistent with the slightly stronger repulsion observed in the SFA scenario.

Fits to the Pressure Curves

Hamaker Constant Estimation

⎛ ∂ 1⎞ 1 ⎛ ∂v′ ⎞ ⎛ ∂p ⎞ ⎜ ⎟ =− 2 ⎜ w ⎟ ⎜ ⎟ vw′ ⎝ ∂p ⎠ N ⎝ ∂μ ⎠ N ⎝ ∂μ vw′ ⎠ N =

χ vw′ 2

(22)

The relative pressure correction stemming from the quadratic term Δp(2) = (1/2)(χ/v′w2)Δμ2 is thus Δp(2) 1 = χ Δp 2 Δp

(23)

where Δp = p(μ) − p(μ′). For Δp = 1000 bar, the quadratic correction accounts for only 0.025 of relative pressure correction, and therefore we safely discard it from our calculations. The TE at fixed number of water molecules and temperature is therefore approximated to linear order as p(μ) ≃ p(μ′) +

μ − μ′ vw

(24) 3

In order to quantify the simulation pressure results, we fit data {pi ± δpi} in a log-plot to a single-exponential decaying function

p(Dw ) = p0 e−Dw / λ

The Hamaker constant in MD simulations can be estimated by pairwise summation over LJ interactions. We approximate the lipid membrane by a homogeneous slab that contains particles with a number density ρlip = 40/nm3, which interact via a −Clip−lip/r6 potential, where the prefactor averaged over all lipid atoms in the used force field is around Clip−lip = 5.6 × 10−3 kJ/ mol nm6. Similarly, we approximate the water between the membranes by a slab of thickness Dw with number density ρw = 33/nm3. The average interaction LJ parameter of water molecules with lipid atoms is Clip−w = 3.8 × 10−3 kJ/mol nm6, and between water molecules themselves Cww = 2.6 × 10−3 kJ/mol nm6. The pairwise summation of two semi-infinite solid slabs across water yields a Hamaker constant in the form,41,55

(25)

where p0 and λ are the fitting parameters. We take the error bars into account by minimizing the functional Σi[(log p(D(i) w) − log pi)/(δpi/pi)]2, which is the sum of squared relative deviations. The fitting parameters are shown in Table 2 and the corresponding functions are plotted in Figure 9. Table 2. Fitting Parameters p0 and λ MD type DLPC DLPC DLPC DPPC

(osmotic stress method) (hydrostatic method) (SFA method) (SFA method)

p0 [kbar] 4.3 5.8 5.7 3.1

± ± ± ±

0.7 0.7 0.5 0.3

λ [nm] 0.42 0.41 0.41 0.47

± ± ± ±

Hsim = π 2(C lip − lipρlip2 − 2C lip − wρlip ρw + Cwwρww 2 )

0.02 0.02 0.02 0.02

(26)

This gives an estimate of the Hamaker constant in the simulations of Hsim ≈ 30 × 10−21 J. Pressure Relations

We can describe the free energy per lipid molecule by two independent variables describing the bilayer structure, for example, area per lipid Alip and bilayer thickness Dlip, and one variable describing the water content in between, such as the water-slab thickness Dw. As we are interested only in isothermic changes at constant chemical potential, we omit the temperature and μ dependence from the notation. The free energy per lipid molecule can then be written as a function -(Dw , Alip , D lip). The partial derivatives of geometrical variables can be obtained via mechanical considerations, namely, Figure 9. MD pressure results (symbols) fitted to the singleexponential function p(Dw) = p0 e−Dw/λ (solid lines).

Head-Group Orientation and Hydrogen Bonds

We now analyze the distance dependence of lipid head-group orientation and number of hydrogen bonds (H-bonds) between water and lipid molecules for all sets of boundary conditions. The head-group orientation is defined via the vector 9135

∂1 = − pzw Alip ∂Dw 2

(27)

∂1 = − pz Alip ∂D lip 2

(28)

∂1 = − pxy (Dw + D lip) ∂Alip 2

(29)


Langmuir

Article

Figure 10. (a) Head-group orientation and (b) average number of H-bonds between water and lipid molecules per surface area as a function of separation.

■

The three independent structural parameters Dw, Alip, and Vlip give rise to three coupled pressure terms, pwz , pz, and pxy. Here pz and pwz are perpendicular pressures in the lipid and water phases, and pxy is the pressure acting in lateral directions. We can express the free energy also in the (Dw,Alip,Vlip) ensemble, where the volume per lipid is expressed as Vlip = (1/2)DlipAlip, which yields ∂1 = − pzw Alip ∂Dw 2

(30)

2Vlip ⎞ V 1 ⎛ ∂⎟⎟ + p lip = − pxy ⎜⎜Dw + z 2 ⎝ Alip Alip ⎠ ∂Alip

(31)

∂= −pz ∂Vlip

(32)

The free energy total differential in the (Dw,Alip,Vlip) ensemble is then 1 d-= − pzw Alip dDw − pz dVlip 2 ⎡ 2Vlip ⎞ V ⎤ 1 ⎛ ⎟⎟ − p lip ⎥ dAlip −⎢ pxy ⎜⎜Dw + z ⎢⎣ 2 ⎝ Alip ⎥⎦ Alip ⎠

(33)

For the case of isotropic pressure pz = pwz = pxy ≡ p, that is, in mechanical equilibrium, we obtain eq 6 in the main text

■

1 1 d- = − pAlip dDw − p dVlip − pDw dAlip 2 2

REFERENCES

(1) Parsegian, V. A.; Fuller, N.; Rand, R. P. Measured work of deformation and repulsion of lecithin bilayers. Proc. Natl. Acad. Sci. U.S.A. 1979, 76, 2750−2754. (2) Parsegian, V. A.; Zemb, T. Hydration forces: Observations, explanations, expectations, questions. Curr. Opin. Colloid Interface Sci. 2011, 16, 618−624. (3) Rand, R. P.; Parsegian, A. V. Hydration forces between phospholipid-bilayers. Biochim. Biophys. Acta 1989, 988, 351−376. (4) Marsh, D. Water-adsorption isotherms and hydration forces for lysolipids and diacyl phospholipids. Biophys. J. 1989, 55, 1093−1100. (5) Besseling, N. A. M. Theory of hydration forces between surfaces. Langmuir 1997, 13, 2113−2122. (6) Norde, W.; Lyklema, J. Protein adsorption and bacterial adhesion to solid-surfaces - a colloidchemical approach. Colloids Surf. 1989, 38, 1−13. (7) Lu, D. R.; Lee, S. J.; Park, K. Calculation of solvation interaction energies for protein adsorption on polymer surfaces. J. Biomater. Sci., Polym. Ed. 1991, 3, 127−147. (8) Afshar-Rad, T.; Bailey, A.; Luckham, P.; MacNaughtan, W.; Chapman, D. Direct measurement of forces between lipid bilayers. Faraday Discuss. Chem. Soc. 1986, 81, 239−248. (9) Lis, L. J.; McAlister, M.; Fuller, N.; Rand, R. P.; Parsegian, V. A. Interactions between neutral phospholipid-bilayer membranes. Biophys. J. 1982, 37, 657−665. (10) Marra, J.; Israelachvili, J. Direct measurements of forces between phosphatidylcholine and phosphatidylethanolamine bilayers in aqueous-electrolyte solutions. Biochemistry 1985, 24, 4608−4618. (11) McIntosh, T. J.; Magid, A. D.; Simon, S. A. Steric repulsion between phosphatidylcholine bilayers. Biochemistry 1987, 26, 7325− 7332. (12) Petrache, H. I.; Gouliaev, I.; Tristram-Nagle, C.; Zhang, R.; Suter, R. M.; Nagle, J. F. Interbilayer interactions from high-resolution x-ray scattering. Phys. Rev. E 1998, 57, 7014−7024. (13) Nagle, J. F.; Tristram-Nagle, S. Structure of lipid bilayers. Biochim. Biophys. Acta 2000, 1469, 159−195. (14) Israelachvili, J. N.; Adams, G. E. Measurement of forces between 2 mica surfaces in aqueouselectrolyte solutions in range 0−100 nm. J. Chem. Soc., Faraday Trans. I 1978, 74, 975−1001. (15) Israelachvili, J. N. Solvation forces and liquid structure, as probed by direct force measurements. Acc. Chem. Res. 1987, 20, 415− 421. (16) Israelachvili, J.; Min, Y.; Akbulut, M.; Alig, A.; Carver, G.; Greene, W.; Kristiansen, K.; Meyer, E.; Pesika, N.; Rosenberg, K.; Zeng, H. Recent advances in the surface forces apparatus (SFA) technique. Rep. Prog. Phys. 2010, 73, 036601. (17) Israelachvili, J. N. Intermolecular and Surface Forces, 2nd ed.; Academic Press Inc.: London, 1991. (18) Schneider, M. J. T.; Schneider, A. S. Water in biological membranes: Adsorption isotherms and circular dichroism as a function of hydration. J. Membr. Biol. 1972, 9, 127−140.

(34)

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS M.K. acknowledges the support of the Alexander von Humboldt Foundation. E.S. acknowledges support from a Marie Curie Intra-European Fellowship within the European Commission 7th Framework Program. Financial support by the DFG via SFB 765 is acknowledged. 9136


Langmuir

Article

full hydration, constant pressure, and constant temperature. Biophys. J. 1997, 72, 2002−2013. (41) Parsegian, V. A. Van der Waals Forces: A Handbook for Biologists, Chemists, Engineers, and Physicists; Cambridge University Press: Cambridge, 2006. (42) Helfrich, W. Steric interaction of fluid membranes in multilayer systems. Z. Naturforsch. 1978, 33, 305−315. (43) Lipowsky, R.; Netz, R. R. Stacks of fluid membranes under pressure and tension. Europhys. Lett. 1995, 29, 345−350. (44) Kummrow, M.; Helfrich, W. Deformation of giant lipid vesicles by electric-fields. Phys. Rev. A 1991, 44, 8356−8360. (45) Netz, R. R. Complete unbinding of fluid membranes in the presence of short-ranged forces. Phys. Rev. E 1995, 51, 2286−2294. (46) Sendner, C.; Horinek, D.; Bocquet, L.; Netz, R. R. Interfacial Water at Hydrophobic and Hydrophilic Surfaces: Slip, Viscosity, and Diffusion. Langmuir 2009, 25, 10768−10781. (47) van der Spoel, D.; Lindahl, E.; Hess, B.; Groenhof, G.; Mark, A. E.; Berendsen, H. J. GROMACS: Fast, flexible, and free. J. Comput. Chem. 2005, 26, 1701−1718. (48) Berendsen, H. J. C.; Postma, J. P. M.; van Gunsteren, W. F.; DiNola, A.; Haak, J. R. Molecular-dynamics with coupling to an external bath. J. Chem. Phys. 1984, 81, 3684−3690. (49) Darden, T.; York, D.; Pedersen, L. Particle mesh Ewald - an n.log(n) method for Ewald sums in large systems. J. Chem. Phys. 1993, 98, 10089−10092. (50) Essmann, U.; Perera, L.; Berkowitz, M. L.; Darden, T.; Lee, H.; Pedersen, L. G. A smooth particle mesh Ewald method. J. Chem. Phys. 1995, 103, 8577−8593. (51) Berendsen, H. J. C.; Grigera, J. R.; Straatsma, T. P. The missing term in effective pair potentials. J. Phys. Chem. 1987, 91, 6269−6271. (52) Widom, B. Some topics in theory of fluids. J. Chem. Phys. 1963, 39, 2808−2812. (53) Frenkel, D.; Smit, B. Understanding Molecular Simulation: From Algorithms to Applications; Academic Press Inc.: New York, 2002. (54) Allen, M. P.; Tildesley, D. J. . Computer Simulation of Liquids; Oxford University Press: New York, 1991. (55) Mahanty, J.; Ninham, B. W. Dispersion Forces; Academic Press Inc.: London, 1976. (56) Luzar, A.; Chandler, D. Structure and hydrogen bond dynamics of water-dimethyl sulfoxide mixtures by computer simulations. J. Chem. Phys. 1993, 98, 8160−8173.

(19) LeNeveu, D. M.; Rand, R. P.; Parsegian, V. A. Measurement of forces between lecithin bilayers. Nature 1976, 259, 601−603. (20) Marčelja, S.; Radić, N. Repulsion of interfaces due to boundary water. Chem. Phys. Lett. 1976, 42, 129−130. (21) Israelachvili, J.; Wennerström, H. Role of hydration and water structure in biological and colloidal interactions. Nature 1996, 379, 219−225. (22) Hammer, M.; Anderson, T.; Chaimovich, A.; Shell, M.; Israelachvili, J. The search for the hydrophobic force law. Faraday Discuss. 2010, 146, 299−308. (23) Schneck, E.; Sedlmeier, F.; Netz, R. R. Hydration repulsion between biomembranes results from an interplay of dehydration and depolarization. Proc. Natl. Acad. Sci. U.S.A. 2012, 36, 14405−14409. (24) Tieleman, D. P.; Marrink, S. J.; Berendsen, H. J. C. A computer perspective of membranes: molecular dynamics studies of lipid bilayer systems. Biochim. Biophys. Acta 1997, 1331, 235−270. (25) Schneck, E.; Netz, R. R. From simple surface models to lipid membranes: Universal aspects of the hydration interaction from solvent-explicit simulations. Curr. Opin. Colloid Interface Sci. 2011, 16, 607−611. (26) Eun, C.; Berkowitz, M. Origin of the Hydration Force: WaterMediated Interaction between Two Hydrophilic Plates. J. Phys. Chem. B 2009, 113, 13222−13228. (27) Hayashi, T.; Pertsin, A.; Grunze, M. Grand canonical Monte Carlo simulation of hydration forces between nonorienting and orienting structureless walls. J. Chem. Phys. 2002, 117, 6271−6280. (28) Pertsin, A.; Platonov, D.; Grunze, M. Origin of short-range repulsion between hydrated phospholipid bilayers: A computer simulation study. Langmuir 2007, 23, 1388−1393. (29) Markova, N.; Sparr, E.; Wadsö, L.; Wennerström, H. A calorimetric study of phospholipid hydration. simultaneous monitoring of enthalpy and free energy. J. Phys. Chem. B 2000, 104, 8053−8060. (30) Schubert, T.; Schneck, E.; Tanaka, M. First order melting transitions of highly ordered dipalmitoyl phosphatidylcholine gel phase membranes in molecular dynamics simulations with atomistic detail. J. Chem. Phys. 2011, 135, 055105-1−055105-11. (31) Patrick S. Coppock, P. S; Kindt, J. T. Determination of Phase Transition Temperatures for Atomistic Models of Lipids from Temperature-Dependent Stripe Domain Growth Kinetics. J. Phys. Chem. B 2010, 114, 11468−11473. (32) McIntosh, T. J.; Simon, S. A. Contributions of hydration and steric (entropic) pressures to the interactions between phosphatidylcholine bilayers - experiments with the subgel phase. Biochemistry 1993, 32, 8374−8384. (33) Chen, Y. L. E.; Gee, M. L.; Helm, C. A.; Israelachvili, J. N.; McGuiggan, P. M. Effects of Humidity on the Structure and Adhesion of Amphiphilic Monolayers on Mica. J. Phys. Chem. 1989, 93, 1051− 1059. (34) Chen, Y. L. E.; Helm, C. A.; Israelachvili, J. N. Measurements of the Elastic Properties of Surfactant and Lipid Monolayers. Langmuir 1991, 7, 2694−2699. (35) Neuman, R. D.; Park, S.; Shah, P. Lateral Diffusion of Surfactant Monolayer Molecules Confined Between Two Solid Surfaces. J. Phys. Chem. 1994, 98, 12474−12477. (36) Evans, D.; Wennerström: The Colloidal Domain: Where Physics, Chemistry, Biology, and Technology Meet (Advances in Interfacial Engineering) (37) Netz, R. R. Debye-Hückel theory for slab geometries. Eur. Phys. J. E 2000, 3, 131−141. (38) Parsegian, V. A. Reconciliation of van derWaals Force Measurements between Phosphatidylcholine Bilayers in Water and between Bilayer-Coated Mica Surfaces. Langmuir 1993, 9, 3625−3628. (39) Tsao, Y.; Evans, D. F.; Rand, R. P.; Parsegian, V. A. Osmotic Stress Measurements of Dihexadecyldimethylammonium Acetate Bilayers as a Function of Temperature and Added Salt. Langmuir 1993, 9, 233−241. (40) Berger, O.; Edholm, O.; Jähnig, F. Molecular dynamics simulations of a fluid bilayer of dipalmitoylphosphatidylcholine at 9137


Hydration Interaction between Phospholipid Membranes: Insight into

Recommend Documents