7086
J. Phys. Chem. B 2008, 112, 7086–7094
Interactions of Lipid Bilayers with Supports: A Coarse-Grained Molecular Simulation Study Chenyue Xing and Roland Faller* Department of Chemical Engineering and Materials Science, UniVersity of CaliforniasDaVis, DaVis, California 95616 ReceiVed: September 11, 2007; ReVised Manuscript ReceiVed: January 15, 2008
The study of lipid structure and phase behavior at the nanoscale is of utmost importance due to implications in understanding the role of the lipids in biochemical membrane processes. Supported lipid bilayers play a key role in understanding real biological systems, but they are vastly underrepresented in computational studies. In this paper, we discuss molecular dynamics simulations of supported lipid bilayers using a coarse-grained model. We first focus on the technical implications of modeling solid supports for biomembrane simulations. We then describe noticeable influences of the support on the systems. We are able to demonstrate that the bilayer system behavior changes when supported by a hydrophilic surface. We find that the thickness of the water layer between the support and the bilayer (the inner-water region in the latter part of this paper) adapts through water permeation on the microsecond time scale. Additionally, we discuss how different surface topologies affect the bilayer. Finally, we point out the differences between the two leaflets induced by the support. Introduction Phospholipid bilayers are abundant and critical components to all biological membranes, including cellular and intracellular membranes. They play an important role in biological systems. Cell membranes consist of phospholipid molecules with inclusions like proteins and sterols.1–4 Supported lipid bilayers provide a simplified yet useful model for investigating the architecture of membranes and for understanding the interaction of cellular membranes with the extracellular media.5,6 Supported bilayers can be formed by vesicle fusion on solid substrates including silica-based surfaces (e.g., glass, aerogels, xerogels) and mica.5,7–9 They maintain the structural and dynamic properties of the biomembranes and provide a bridge between natural systems and artificial materials. Thus, supported lipid bilayers have attracted extensive interest experimentally.10,11 Moreover, they have attracted attention because of their potential applications in biotechnology. Molecular computer simulation is a powerful tool to reveal the detailed structures and physical mechanisms of biomembranes.12–16 It can provide direct insight into highly complex biological systems, elucidate principles, and test ideas for developing biotechnical products.17 Simulations of lipid bilayers have involved a large number of different lipids, including DPPC, DLPC, DSPC, and so forth, and many key phenomena in membranes have been investigated, such as lateral mobility, self-assembly, phase transitions, or phase separation.2,17–22 The interactions between lipids and other molecules, for example, proteins, sugars, or alcohols, have also been studied computationally.23–27 Almost all studies were performed by modeling free-standing lipid bilayers in solution. Although supported lipid bilayers have been extensively studied in experiments, we are only aware of one recent simulation of a bilayer system on a substrate.28 The understanding of the changes of structural and dynamic characteristics of the bilayer upon depositing on a solid support requires large system size and long simulation time. * To whom correspondence should be addressed.
Therefore, there is an urgent need for supported bilayer models, especially on the mesoscopic scale. In this contribution, we perform coarse-grained molecular simulations to study supported lipid bilayer systems and various properties of the lipids supported on a hydrophilic surface. Coarse-graining (CG) is a very powerful approach that balances the specificities that researchers are interested in and the computing capacity.17,20,29,30 Mapping 4-8 heavy atoms into one interaction site makes the models feasible to describe large time and length scale properties, that is, the collective phenomena of membranes.17,18,20,24,31 Collective phenomena involve a large number of molecules and occur on length scales from several tens of nanometers to micrometers and time scales of microseconds or even milliseconds. Most of the collective phenomena, for instance, the phase behaviors, are not dictated by the detailed molecular structure but by the more generic nature of the molecule, for example, the amphiphilicity. By eliminating atomistic details, CG raises the computational efficiency and thus can be used to study large-scale behaviors more intensively. Many different computational techniques and models have been developed and dedicated to different questions of interest. Atomistic models are accurate in describing not only the molecular structure but also the intermolecular interactions such as chemical bonds and electrostatic and van der Waals interactions. Because of the nature of the interactions, atomistic models employ a time step as short as the a tenth of the period of the fastest mode in the system,32 often the stretching of a covalent bond or angle. Atomistic models include at least every nonhydrogen atom into the system, and with the short time step, they are able to study membranes of a few tens of nanometers over tens of nanoseconds using currently available computing facilities. They are widely used in the study of the local structure and dynamics of a membrane16,33,34 as well as the effect of inclusions such as sterols on the structure of the membrane.26,27 As processes like self-assembly, phase transitions, and phase separations occur on length and time scales beyond atomistic
10.1021/jp077305l CCC: $40.75 2008 American Chemical Society Published on Web 05/08/2008
Interactions of Lipid Bilayers with Supports
J. Phys. Chem. B, Vol. 112, No. 23, 2008 7087
Figure 1. Coarse-grained representation of supported dipalmitoylphosphatidylcholine (DPPC) bilayers.
capabilities, CG models have to be applied. Reducing the number of degrees of freedom by grouping several atoms into one effective particle and eliminating short-range dynamics, CG models speed up simulations and access the length and time scales of collective phenomena.17,20,29–31,35,36 The CG model that we use in this paper was developed by Marrink et al.31 and has been proven to be semiquantitatively accurate in many lipid simulations.18,20,31 It has not yet been used to study supported systems. For our particular supported bilayer system, the way in which water is modeled has to be adjusted. In simulations of supported lipid systems, we have demonstrated that the water model is crucial.37 Water is essential as an environment for biomembranes. Its particular, chemical and physical characteristics lead to a complex phase behavior. This results in the well-known anomalous density behavior in the liquid phase, the possible existence of a liquid-liquid (LL) phase transition, and the existence of many solid phases.38 A LL transition of water has been seen in simulations of a system with a solid support39 and hence is closely related to our study of supported biomembranes. In this paper, we use DPPC, dipalmitoylphosphatidylcholine, to study the influence of a hydrophilic substrate on the lipid bilayer. As the zwitterionic headgroups will be attracted by the hydrophilic substrate, we would expect low mobility of the lipids near the substrate. There is experimental evidence showing that a free-standing bilayer has faster dynamics than a supported one.40 The rest of this paper is organized as follows. The Simulation Details are presented in the next section. Structural and dynamical properties of the simulated systems are examined in the following. Comparisons between free-standing and supported lipid systems are also discussed, and we end with Conclusions. Simulation Details We performed molecular dynamics simulations using the GROMACS simulation suite.41 Lipids and water were simulated using the coarse-grained (CG) model proposed by Marrink, and the initial configuration was derived from a study of a pure unsupported dipalmitoylphosphatidylcholine (DPPC) bilayer.31 Model parameters were directly transferred from that work for the free-standing DPPC bilayer simulation, while some modifications were made to water molecules in the supported bilayer simulation. In one of our earlier works, we have shown that the water in Marrink’s CG model easily freezes when a solid phase is present.37 The same observation was reported in one of Marrink’s recent papers.42 Defining the original water-water nonbonded interaction strength in the CG model31 as 100%, we
TABLE 1: Nonbonded Interaction Strength between the Support Particles (Su) and All Other Particles31a N Type Su
subtype
Q
Su
P
0
d
a
da
C
0
d
a
da
N/A
In
In
In
In
In
R
In
In
In
In
a
Level of interaction: N/A (no interaction), In (intermediate), or R (repulsive).
Figure 2. Schematic representation of simulation geometries in the study of a supported lipid bilayer. Left: view from the X-Z plane; Right: view from Y-Z plane. All numbers are in nanometers (nm).
have tested several lower interaction strengths (90, 82, 76, and 68%) to study the relationship between the water freezing and the water-water interaction strength. Details can be found in our earlier publication.37 The water-water nonbonded interaction has been weakened to 76% of the original value in the following simulations of supported bilayers. This was determined so that water exhibits bulk properties in the supported bilayer system, except for the region very close to the support where clear water layering was observed. A pure CG water system of 4129 water particles on a support finds, in spite of the ordered water structure, the water particles in the fluid phase with a diffusion coefficient of (4.71 ( 0.03) × 10-5 cm2/s at 323 K. The original CG model yields a diffusivity of 2 × 10-5 cm2/s for water at 300 K and a density of 990 kg/m3.31 The bulk density is about 900 kg/m3 with the scaled interaction strength. After testing a few scaling factors for the interaction strength,37 we found the best compromise between fluidity and density. Similar water behavior in the presence of a support has been reported earlier.39,44–46 Let us now discuss our specific system. Figure 1 shows the coarse-grained representation of DPPC. The phospholipid DPPC consists of 12 CG interaction sites, with 2 charged sites (type Q) representing the PC headgroup, positive for choline and negative for phosphate. The glycerol linkage is modeled by two intermediate hydrophilic sites (type N). Each tail contains four hydrophobic sites (type C). The solvent, that is, water, is
7088 J. Phys. Chem. B, Vol. 112, No. 23, 2008
Xing and Faller
Figure 3. Different surface topologies of the solid support. From left: (a) flat surface, system A, (b) four dips on surface corresponding to 1/4 of the total surface area, system B, (c) random distributed rough surface, system C.
Figure 4. Density profiles of parts in DPPC bilayers at the end of 12 µs simulations (black, water; blue dotted, outer leaflet; red, inner leaflet); (a) supported DPPC system A, (b) free-standing DPPC, (c) supported DPPC system B, (d) supported DPPC system C.
modeled by single hydrophilic sites representing four real water molecules. The CG model above was parametrized to match atomistic simulations on the structural, dynamic, elastic, and thermodynamic properties. The CG sites interact pairwise via a Lennard-Jones (LJ) potential. Five different LJ interaction strengths were defined to achieve the interactions ranging from strongly repulsive to strongly attractive. A screened Coulomb interaction models the electrostatic interaction between the zwitterionic headgroups. The choline group bears a +1 charge, and the phosphate group bears a -1. The screening of the Coulomb potential was performed by a cutoff at rc ) 1.2 nm, and the force was smoothly shifted to zero at the cutoff from r ) 0. We used a neighbor-list cutoff of rNL ) 1.2 nm, and the interaction cutoff was rc ) 1.2 nm, shifting the Lennard-Jones potential smoothly to zero from rs ) 0.9 nm.41 As the neighborlist cutoff is generally required to be 0.1-0.3 nm larger than the interaction cutoff, we also did a short simulation (40 ns) with rNL ) 1.4 nm. The results obtained did not differ significantly from those of the simulation with rNL ) 1.2 nm. We therefore used rNL ) 1.2 nm for all production runs in order to directly follow the original implementation.31 For computa-
tional efficiency, all CG sites have the same mass of 72 amu. The periodic boundary condition (PBC) was applied in all simulations. To better understand a supported lipid bilayer system, we first performed simulations of a free-standing bilayer. All data discussed below were obtained from bilayer systems containing 512 DPPC lipids (256 lipids for each leaflet) unless otherwise stated. The free-standing bilayer was simulated under NPT condition. The system was controlled by a Berendsen thermostat with a time constant of T ) 1.0 ps and a Berendsen barostat with T ) 1.0 ps47 using a compressibility of 5 × 106 bar-1, so that the system was tension-free. The reference temperature was 323 K, above the main phase transition temperature of DPPC, which is Tc) 314K.48 The pressure coupling was applied to the three Cartesian directions independently at 1 bar. The CG model allows a larger time step than that of traditional atomistic simulations. We used 40 fs. On the basis of comparison of diffusion constants in the CG model and in atomistic models, a conversion factor of 4 is needed due to the effective speed up in the dynamics of CG simulations.31 All times reported in the remaining of the paper are effective times (e.g., 20 ns of CG
Interactions of Lipid Bilayers with Supports
J. Phys. Chem. B, Vol. 112, No. 23, 2008 7089
Figure 5. Left (a): mass of inner-water (system A). Right (b): number of water particles in the inner-water region versus time (system A).
Figure 6. Number of water particles in the inner-water region versus time (black: simulation D, starting with 938 inner-water particles; red: simulation E, starting with 1068 inner-water particles).
simulation time corresponds to 80 ns of effective time). After 80 ns of equilibration, the system reached an average area per headgroup of 62.3 Å2. A free bilayer with our changed water interaction led to about 76 Å2 (see below). The system contained 10372 CG waters. As we later insert a support in close proximity to the bottom of the bilayer once the free-standing bilayer system reaches equilibrium, we need a thicker bulk water (>3 nm) above, such that under PBC, the nearest-replicated support above the bilayer will be distant enough. The equilibrated free-standing bilayer system was transferred onto a support of 1760 completely fixed particles of type N in the Marrink model31 on a grid with a spacing of 0.3 nm. The support particle has the same 0.47 nm radius as all other CG. The lattice is smaller than the particle radius so that adjacent particles overlap, leaving no gap in between. In the x and y directions, the box size was carefully chosen to be commensurate with the lattice spacing under PBC, leading to a highly smooth supporting surface, infinite in the x-y plane. Table 1 lists the nonbonded interaction strength between the support and all other particles. Support particles do not interact among themselves. The initial system size was 12 × 13.2 × 12.7 nm3 (cf. Figure 2). The area was chosen as a best compromise between the
equilibrium area per headgroup of a free-standing DPPC bilayer (62.3 Å2) and the commensurability with the smooth surface. Constant pressure simulations in the bilayer plane were not possible because the interaction density of the support would change with rescaling. We therefore kept the area of the system constant throughout the simulations while allowing the normal pressure to fluctuate around the reference value (1 bar). We call this an NAPzT ensemble. Note that in most simulations, we only use a single support wall, which, through PBC, interacts as two effective walls. We have also performed simulations with two identical solid supports, one above the other, both fixed spatially with 1 nm of vacuum in between. After 40 ns, the vacuum region stayed the same as not a single particle managed to pass the supports. In addition, the density profile of this system did not differ significantly from that of a single-solid-support system, which implies that the possible interaction of water particles with waters separated by the solid support does not affect the whole system. All other simulation parameters were kept the same as those in the free-standing bilayer simulation. Analyses of both systems were carried out after an equilibration period of 500 ns. We additionally performed simulations with topologically rough surfaces (cf. Figure 3) as real surfaces are not perfectly flat. Figure 3a is a perfectly flat surface, as described above (system A). Figure 3b shows a surface with four big dips which, altogether, equal to 1/4 of the total surface area; the depth of the dips is 0.3 nm (system B). For a third system C, we randomly picked half of the support particles and increased their z coordinates by 0.3 nm, which produced a rough surface, shown in Figure 3c. The choice of the dip size is not random but rather guarantees that any pair of adjacent particles on different levels is only (0.32 + 0.32)1/2 ) 0.42 nm apart from each other, allowing no solvent molecules to pass through. Results and Discussion Density Profiles. One of the most obvious influences of a support to the lipid bilayer is the change in density profiles. Figure 4 shows the comparison of density profiles of the supported and free-standing systems. At the beginning of the production run, there were 895 waters below (i.e., the innerwater region) and 9477 waters above (i.e., the outer-water region) the bilayer. The number of waters in the simulations
7090 J. Phys. Chem. B, Vol. 112, No. 23, 2008
Xing and Faller
Figure 7. Snapshots of DPPC systems at the end of the 12 µs simulation. Left (a): supported DPPC system A. Right (b): free-standing DPPC system.
Figure 8. Double logarithmic plot of mean-squared displacement (MSD) of supported DPPC system A and free-standing DPPC system.
was chosen such that the thickness of the water between the support and the bilayer corresponded to about two layers of water molecules. There is no direct experimental evidence on the exact amount of the inner-water, but two layers is a commonly hypothesized value.49–51 Figure 5a shows the mass
of the inner-water in system A at the end of the production run (12 µs). By doing a rough weighted average of the z positions of choline groups in the inner leafet, the bilayer is 0.8 nm away from the support. Before approaching the bilayer, there are two plateaus corresponding to two water layers due to the influence of the support at z ) 0. In Figure 4a, there is a region of about 2 nm of water at the vicinity of the surface exhibiting structured behavior. Sharp peaks in the density profile water layers are found. Because of the periodic boundary condition, a corresponding region can be seen at the top part of the box as well. The ordering of water induced by the support was described as a liquid-liquid (LL) transition in ref 39. They have studied thermodynamic and structural properties of water confined between two smooth hydrophobic plates. Although their results are not fully comparable, as we have a semihydrophilic surface while they had hydrophobic ones, there still is evidence that the water behaves liquid-like in its dynamics while arranging in a solid-like local structure.44–46 We have measured the meansquared displacement (MSD) of the ordered water in system A, and its lateral diffusivity was found to be Dwater ) (1.18 ( 0.12) × 10-5 cm2/s. The lateral diffusivity of nonordered bulk water in the same system was Dwater ) (4.23 ( 0.09) × 10-5 cm2/s. The ordered and nonordered water are both fluid, given that their mobility is on the same order of magnitude, though
Figure 9. Pressure tensor components in system A over the course of several microseconds.
Interactions of Lipid Bilayers with Supports
Figure 10. Distribution of NC3 groups (choline) along the z direction in supported DPPC systems. From top: (a) supported DPPC system A, (b) system B, (c) system C.
the ordered water is slower. The lateral diffusivity of the original-strength interacting water in the vicinity of a free bilayer reported in ref 31 is 2.0 × 10-5 cm2/s at 300 K. Therefore, we believe that such a LL transition of water was induced by the insertion of a solid support. The waters are exposed to the attractive force from the support, and the support breaks the symmetry of the chemical environment of the bilayer in the normal direction. As a comparison, in Figure 4b, all water shows bulk density. We found very slow transport of water through the membrane. The inner-water region slowly thickened. The number of the inner-water particles changed as time increased, as plotted in Figure 5b. Along the 12 µs trajectory, the number of the innerwater particles increased roughly linearly until 6 µs when it reached equilibrium. The permeability coefficient of water through the bilayer can be calculated based on the equation P ) KD/d,52,53 where D (cm2/s) is the diffusion coefficient of
J. Phys. Chem. B, Vol. 112, No. 23, 2008 7091 water, K is the partition coefficient, and d is the thickness of the hydrophobic region of the bilayer. Taking K ) 4.2 × 10-5, D ) 2.8 × 10-5 cm2/s, and d ) 3.8 nm, we got P ) 3.1 × 10-3 cm/s. This result drops in the range reported by both computational and experimental studies.52–57 According to Figure 5b, we conclude that the equilibrium amount of the innerwater particles is 990 ( 10. To verify that 990 ( 10 is a convincing value for the equilibration of the inner-water, we have done two more simulations (i.e., simulations D and E). Starting from the final configuration from the 12 µs production run, we moved 60 water particles from the inner-water region to the outer-water region to make the initial configuration for simulation D. For simulation E, we moved the water particles in the opposite way, moving 60 waters from the outer-water region to the inner-water. After a few hundreds of ns of relaxation, we ran both simulations for 5 µs. By tracking the number of inner-water particles, we found that both systems head toward the proposed equilibrium inner-water amount, 990 ( 10 (cf. Figure 6). This particular discovery provides strong evidence for approximately two layers as a typical amount for the inner-water region. It also implies that it takes at least several microseconds before the supported bilayer reaches equilibrium. In a supported system, the two leaflets no longer have the same environment. Thus, in Figure 4a, we clearly see the heterogeneity in the density profiles of the lipids. The inner (proximal to the surface) leaflet has similar properties as those of the ordered water layers. The headgroups of the inner leaflet are attracted to the support, leading to a rearrangement into several vertical levels, that is, we find locally almost solid-like order in the vicinity of the support. The interaction between the support particles and the moving particles is not long-ranged. Therefore, the outer leaflet is weakly affected by the support, and its density profile is similar to a leaflet in the free-standing bilayer system (Figure 4b). We have to keep in mind here that both leaflets have the same number of lipids, and flip-flop is not observed in the time scales reported here. The change in density profiles can also be validated by the structure visualization. Figure 7a and b shows snapshots of the supported and free-standing bilayers. Although the lateral size of the systems is only around 10 nm on the side, we can still observe weak undulations in the free-standing bilayer (Figure 7b). The bilayer in Figure 7a is flatter because the headgroups of the inner leaflet are arranged at two z positions, suppressing the undulation. Diffusion. Experiments have shown that a supported lipid bilayer is at least three times less mobile compared to a freestanding bilayer.40 To validate this observation computationally and to examine the dynamic behavior of the lipids, we have calculated the MSD of each leaflet in both supported and freestanding systems, respectively. From the MSD, one can calculate the lateral diffusion coefficient using the Einstein relation
Dlateral ) lim
1
tf∞ 4t
〈[rb(t)]2 〉
where 〈[r b(t)]2〉 is the mean-squared displacement (MSD) averaged over all lipids in the selected group. In Figure 8, the MSDs for the supported bilayer system A and the free-standing bilayer system have been plotted up to10 µs. We have found that in the free-standing bilayer system, D(free) ) (2.25 ( 0.14) × 10-7 cm2/s is in good agreement with experimental findings that range from 1.0 × 10-7 to 1.5 × 10-7 cm2/s.58,59 The free-standing bilayer with scaled water has a faster lateral diffusion (see below). The results are also in accord with the lateral diffusivity reported in ref 31. In the supported bilayer system, the inner leaflet was severely slowed
7092 J. Phys. Chem. B, Vol. 112, No. 23, 2008
Xing and Faller
Figure 11. Supported pure water (W) system with n ) 0.1 antifreeze particles (WF). Left (a): snapshot of the final configuration. Right (b): density profiles of the system (black), normal CG water W (red), and antifreeze particles WF (green).
down by the attractive surface, while the outer leaflet was only weakly affected. Dinner(sup) ) 1.95 × 10-8 cm2/s and Douter(sup) ) (2.78 ( 0.01) × 10-7 cm2/s. The inner leaflet is thus slower than the outer leaflet by about an order of magnitude. This is because the lipids in the inner leaflet have been restrained by the attractive support, which therefore restrained the motion of the lipids in the bilayer plane. The extent of the decreasing mobility is also in reasonable agreement with experiments. In addition, the lateral diffusion coefficient reflects the change in supporting surface topologies. The lateral diffusion coefficient that we have found for the inner leaflets in systems B and C were 6.41 × 10-8 and (1.79 ( 0.01) × 10-7 cm2/s, respectively. This is faster than that of the inner leaflet in system A but slower than that in a free-standing bilayer. This trend in changing mobility of lipids correlates with decreasing ordering of the lipids induced by the supports from system A to C. The outer leaflet in the supported bilayer is slightly faster than that in the control simulation. This may be due to the weaker coupling to the bottom leaflet due to the suppressed undulations or due to the different water interactions. Pressure and Surface Tension. Free-standing bilayers are tension-free in a NPT ensemble. However, the NAPzT ensemble, which was used in the study of the supported system, does not ensure this. The system was built from a tension-free state from the free-standing bilayer simulation, and the lateral area was then kept constant. Under such conditions, the lateral pressures obtained negative average values with the normal pressure still coupled to 1 bar. Pressure components are plotted for the last 2 µs of a 12 µs simulation of system A in Figure 9. In the simulations presented in this paper, the introduction of the support induced the inhomogeneity of the systems and the anisotropy in the pressure. Similar trends have been reported by Kaznessis et al. for the simulation of lipid monolayers exposed to a wall potential.60 A surface tension arises from the difference between the normal pressure and the lateral pressure.
γsys )
∫-∞∞ (PN(z) - PT(z))dz
where the normal pressure is PN ) Pzz and the lateral pressure is PT ) (Pxx+Pyy)/2. The negative lateral pressure in the supported bilayer system therefore suggests a positive surface
tension and furthermore implies a smaller equilibrium area per headgroup (at least for the proximal leaflet) than that of a freestanding system. Surface Topology and Bilayer Structure. We now discuss the influence of surface topology on bilayer structure (cf. Figure 3). The peak values of the density of water layers and lipid particles decrease as the roughness of the surface increases (Figure 4). Histograms for the number of choline groups (NC3 particles) located at different z coordinates are shown in Figure 10. Two clear layers of NC3 can be seen in the flat surface system A (Figure 10a). Layers can still be observed in system B, but they are farther away from the surface than those in A (Figure 10b). Seventy-five percent of the support particles are at z ) 0.3, and the rest are at z ) 0.0. The equilibrium distance between the support and the headgroups determined the higher position of the lipids. An even clearer decrease in the ordering can be seen in the system C (Figure 10c), in which the distribution of NC3 particles in the inner leaflet becomes similar to that in the outer leaflet. In short, all density-related data imply that the lipid bilayer structure follows the surface topology. In addition, the mobility of lipids is also influenced by topology. With the increase in roughness from system A to C, the ordering of the lipids induced by the supports decreased; thus, the lipids became more mobile. The results of three systems with different surface topologies are particularly important in that the supported bilayer model in this paper will be extended to porous supports in future so that the simulations can be helpful in reconciling the experimental study of supported bilayers on different substrate materials. Discussion of CG Water Models. As discussed above, in order to accommodate our solid support, the CG water was adjusted. The original CG water in Marrink’s model31 has a known problem of its freezing temperature, which is too high compared to that of real water. Instead of 273 K, it generally freezes at 290 ( 5 K.31 It freezes very rapidly once a nucleation site is formed in the system, even at a high temperature, and the freezing process is irreversible.42 The solid support in our lipid system acts as a large nucleation site, triggering freezing. Before scaling down the water-water interaction strength, we tried a series of other operations to reduce the nucleation effect,
Interactions of Lipid Bilayers with Supports for example, decreasing the hydrophilicity of the support or increasing the interaction size of the water beads. None has resolved the freezing problem. Very recently, Marrink et al. released an extended version42 of the original model. In the extended version, a new particle type was introduced as an antifreeze agent. With merely a mole fraction of n ) 0.1 of such particles in the pure water system, the freezing temperature was lowered to below 250 K. It lowered both the density and self-diffusion constant by 10% but did not have a severe impact on bilayer properties.42 Hence, we tested the model and replaced 10% of the original CG waters (W) with the antifreeze particles (WF) and simulated that pure water system with a solid support. Figure 11a and b shows a snapshot and the density profile of supported pure water with n ) 0.1 antifreeze particles. The middle part of the system stays in a fluid state with most of the antifreeze particles. The solid support still affects more than a 5 nm thick water layer. Furthermore, the ordered water region contains very few of the antifreeze particles This implies that the antifreeze particles demix from normal water and do not overcome the driving force of nucleation at the region close to the support. We tried to increase the amount of antifreeze particles to 15%. That did not reduce the nucleation but brought down the water density to as low as 800 kg/m3. Therefore, we have to come to the conclusion that the antifreeze particles only help when there is no large nucleation site or surface present. The best alternative method is to scale down the interaction strength between water molecules so that the nucleation impact of the solid support could be shorter-ranged. Of course, the scaled water model has a problem in reproducing some properties of the free-standing bilayer. For example, the scaled water model yields a larger area per headgroup (76 Å2) and a slightly thinner bilayer (4 nm) for the free-standing DPPC bilayer. The lateral diffusion constant of lipids is also increased to (3.58 ( 0.11) × 10-7 cm2/s. The density of bulk water in such a system is roughly 900 kg/m3. As in the supported bilayer system, the heterogeneity in pressure components tends to shrink the lateral area. A simulation in a NPT ensemble is not currently possible for supported systems as we want to keep the interaction strength of the support constant. We would expect the real equilibrium area per headgroup to be a product of the two competing forces, the surface tension due to the support and the loose packing due to weaker water. Conclusion We have presented results of coarse-grained MD simulations of a DPPC bilayer on solid supports. We have found significant differences between supported DPPC and free-standing DPPC systems. First, the density profiles of the bilayer were affected by the insertion of the support. The spatial arrangement of the lipid particles in the vicinity of the support is strongly ordered. We also found that water permeates through on the microsecond time scale. The resulting water cushion is about 2 coarse-grained water layers thick. This corresponds to about 1 nm. The dynamic properties of the supported bilayer are different from those of a free-standing bilayer. The inner leaflet is restrained by the attractive interaction of the support. Thus, the lateral diffusion of the inner leaflet is slower than that of the outer leaflet by roughly 1 order of magnitude. The mobility of the outer leaflet is similar to that in the free-standing system. Therefore, the supported system experienced heterogeneity in the dynamics. Due to the NAPzT ensemble in the simulation of supported bilayer systems, the lateral pressures obtained negative average values when keeping the lateral area at the equilibrium area of a free-standing system with, however, a different water interac-
J. Phys. Chem. B, Vol. 112, No. 23, 2008 7093 tion. The system was no longer tension-free. Given that, one expects a smaller equilibrium area for a supported DPPC bilayer. The influence of the surface topology of the support to the bilayer was also studied. Three difference surface topologies were implemented. The densities and positions of the headgroups demonstrate that a flat surface has the strongest effect on confining the lipids. Though still slower than that in a freestanding bilayer, the mobility of lipids increases with the roughness of the surface. The current state of the computational study of the supported bilayer systems is very encouraging, despite the fact that the current scaled down water model needs further improvement. The key properties are at least qualitatively in agreement with experimental and atomistic results. The coarse-grained model applied in this particular study is sufficient to semiquantitatively study the structural and dynamical properties of a supported bilayer system. It is also possible to extend it to larger systems and to understand the system behavior on longer time scales. On the basis of this first attempt of simulating a supported bilayer, work on simulating such a system on an atomistic level is underway in order to provide deeper insight into the nature of the atomic-scale phenomena, which the current experimental techniques are not able to examine. Acknowledgment. This project was supported by the NSF NIRT Program (Grant CBET 0506662). References and Notes (1) Kusumi, A.; Tsuda, M.; Akino, T.; Ohnishi, S.; Terayama, Y. Biochemistry 1983, 22, 1165. (2) Falck, E.; Patra, M.; Karttunen, M.; Hyvonen, M. T.; Vattulainen, I. Biophys. J. 2004, 87, 1076. (3) Pebay-Peyroula, E.; Rosenbusch, J. P. Curr. Opin. Struct. Biol. 2001, 11, 427. (4) London, E. Biochim. Biophys. Acta 2005, 1746, 203. (5) Sackmann, E. Science 1996, 271, 43. (6) Hianik, T. Acta Phys. SloVaca 2006, 56, 687. (7) Groves, J. T.; Ulman, N.; Boxer, S. G. Science 1997, 275, 651. (8) Richter, R. P.; Brisson, A. R. Biophys. J. 2005, 88, 3422. (9) Weng, K. C.; Stalgren, J. J. R.; Risbud, S. H.; Frank, C. W. J. Non-Cryst. Solids 2004, 350, 46. (10) Groves, J. T.; Dustin, M. L. J. Immunolog. Methods 2003, 278, 19. (11) Cornell, B. A.; BraachMaksvytis, V. L. B.; King, L. G.; Osman, P. D. J.; Raguse, B.; Wieczorek, L.; Pace, R. J. Nature 1997, 387, 580. (12) Hogberg, C. J.; Lyubartsev, A. P. J. Phys. Chem. B 2006, 110, 14326. (13) Berger, O.; Edholm, O.; Jahnig, F. Biophys. J. 1997, 72, 2002. (14) Feller, S. E.; Yin, D. X.; Pastor, R. W.; MacKerell, A. D. Biophys. J. 1997, 73, 2269. (15) Marrink, S. J.; Mark, A. E. J. Phys. Chem. B 2001, 105, 6122. (16) Tieleman, D. P.; Berendsen, H. J. C. J. Chem. Phys. 1996, 105, 4871. (17) Muller, M.; Katsov, K.; Schick, M. Phys. Rep. 2006, 434, 113. (18) Faller, R.; Marrink, S. J. Langmuir 2004, 20, 7686. (19) Goetz, R.; Lipowsky, R. J. Chem. Phys. 1998, 108, 7397. (20) Marrink, S. J.; Risselada, J.; Mark, A. E. Chem. Phys. Lipids 2005, 135, 223. (21) Mouritsen, O. G. Chem. Phys. Lipids 1991, 57, 179. (22) Kandt, C.; Ash, W. L.; Tieleman, D. P. Methods 2007, 41, 475. (23) Jakobsson, E. Trends Biochem. Sci. 1997, 22, 339. (24) Sum, A. K.; Faller, R.; de Pablo, J. J. Biophys. J. 2003, 85, 2830. (25) Dickey, A. N.; Faller, R. J. Polym. Sci., Part B: Polym. Phys. 2005, 43, 1025. (26) Hofsass, C.; Lindahl, E.; Edholm, O. Biophys. J. 2003, 84, 2192. (27) Tu, K. C.; Klein, M. L.; Tobias, D. J. Biophys. J. 1998, 75, 2147. (28) Heine, D. R.; Rammohan, A. R.; Balakrishnan, J. Mol. Simul. 2007, 33, 391. (29) Lopez, C. F.; Moore, P. B.; Shelley, J. C.; Shelley, M. Y.; Klein, M. L. Comput. Phys. Commun. 2002, 147, 1. (30) Shelley, J. C.; Shelley, M. Y.; Reeder, R. C.; Bandyopadhyay, S.; Moore, P. B.; Klein, M. L. J. Phys. Chem. B 2001, 105, 9785. (31) Marrink, S. J.; de Vries, A. H.; Mark, A. E. J. Phys. Chem. B 2004, 108, 750. (32) MacKerell, A. D. J. Comput. Chem. 2004, 25, 1584.
7094 J. Phys. Chem. B, Vol. 112, No. 23, 2008 (33) Aman, K.; Lindahl, E.; Edholm, O.; Hakansson, P.; Westlund, P. O. Biophys. J. 2003, 84, 102. (34) PasenkiewiczGierula, M.; Takaoka, Y.; Miyagawa, H.; Kitamura, K.; Kusumi, A. Biophys. J. 1999, 76, 1228. (35) Shelley, J. C.; Shelley, M. Y.; Reeder, R. C.; Bandyopadhyay, S.; Klein, M. L. J. Phys. Chem. B 2001, 105, 4464. (36) Nielsen, S. O.; Lopez, C. F.; Srinivas, G.; Klein, M. L. J. Phys.: Condens. Matter 2004, 16, R481. (37) Bennun, S. V.; Dickey, A. N.; Xing, C. Y.; Faller, R. Fluid Phase Equilib. 2007, 261, 18. (38) Vega, C.; Abascal, J. L. F.; Sanz, E.; MacDowell, L. G.; McBride, C. J. Phys.: Condens. Matter 2005, 17, S3283. (39) Kumar, P.; Buldyrev, S. V.; Starr, F. W.; Giovambattista, N.; Stanley, H. E. Phys. ReV. E 2005, 72. (40) Przybylo, M.; Sykora, J.; Humpolickova, J.; Benda, A.; Zan, A.; Hof, M. Langmuir 2006, 22, 9096. (41) Lindahl, E.; Hess, B.; van der Spoel, D. J. Mol. Model. 2001, 7, 306. (42) Marrink, S. J.; Risselada, H. J.; Yefimov, S.; Tieleman, D. P.; de Vries, A. H. J. Phys. Chem. B 2007, 111, 7812. (43) Deleted in proof. (44) Ruan, C. Y.; Lobastov, V. A.; Vigliotti, F.; Chen, S. Y.; Zewail, A. H. Science 2004, 304, 80. (45) Giovambattista, N.; Rossky, P. J.; Debenedetti, P. G. Phys. ReV. E 2006, 73. (46) Israelachvili, J. N. Surf. Sci. Rep. 1992, 14, 109.
Xing and Faller (47) Berendsen, H. J. C.; Postma, J. P. M.; van Gunsteren, W. F.; Dinola, A.; Haak, J. R. J. Chem. Phys. 1984, 81, 3684. (48) Silver, B. Physical Chemistry of Membranes: An Introduction to the Structure and Dynamics of Biological Membranes; Kluwer Academic Publishers: Dordrecht, The Netherlands, 1985. (49) Charitat, T.; Bellet-Amalric, E.; Fragneto, G.; Graner, F. Eur. Phys. J. B 1999, 8, 583. (50) Johnson, S. J.; Bayerl, T. M.; McDermott, D. C.; Adam, G. W.; Rennie, A. R.; Thomas, R. K.; Sackmann, E. Biophys. J. 1991, 59, 289. (51) Bayerl, T. M.; Bloom, M. Biophys. J. 1990, 58, 357. (52) Finkelstein, A. J. Gen. Physiol. 1976, 68, 127. (53) Paula, S.; Volkov, A. G.; VanHoek, A. N.; Haines, T. H.; Deamer, D. W. Biophys. J. 1996, 70, 339. (54) Jansen, M.; Blume, A. Biophys. J. 1995, 68, 997. (55) Marrink, S. J.; Berendsen, H. J. C. J. Phys. Chem. 1994, 98, 4155. (56) Sugii, T.; Takagi, S.; Matsumoto, Y. J. Chem. Phys. 2005, 123. (57) Bemporad, D.; Essex, J. W.; Luttmann, C. J. Phys. Chem. B 2004, 108, 4875. (58) Vaz, W. L. C.; Clegg, R. M.; Hallmann, D. Biochemistry 1985, 24, 781. (59) Konig, S.; Pfeiffer, W.; Bayerl, T.; Richter, D.; Sackmann, E. J. Phys. II 1992, 2, 1589. (60) Kaznessis, Y. N.; Kim, S. T.; Larson, R. G. Biophys. J. 2002, 82, 1731.
JP077305L