INTRODUCTION

TOP ABSTRACT INTRODUCTION SIMULATION METHODS AND... RESULTS DISCUSSION REFERENCES

The role of cholesterol in membranes has been extensively investigated for many years (Yeagle, 1993; Presti, 1985; Finegold, 1993), but many unanswered questions remain. For example, whereas the manner by which cholesterol affects the lipid chain phase transition in phosphatidylcholines has been thoroughly documented (McMullen et al., 1994) experimentally, the underlying molecular interactions are not well understood. At very low (less than ~10%) cholesterol concentrations, fluid phase lipid bilayers retain fluid-like properties. With higher cholesterol concentrations, the bilayer becomes a mixture of fluid-like and "liquid-ordered" phases. The liquid ordered phase is characterized by well-ordered hydrocarbon chains and increased bilayer thickness and increased area compressibility modulus (Needham et al., 1988; Touard et al., 1999). Despite general agreement on the qualitative aspects of lipid-cholesterol phase behavior and the properties of the phases, recently proposed phase diagrams (Vist and Davis, 1990; Almedia et al., 1992; Thewalt and Bloom, 1992; McMullen and McElhaney, 1995) disagree in the locations of boundaries between phases, the location of critical points, and the nature of coexisting solid-ordered phases. At a molecular level, the complex interplay of hydrophobic and electrostatic interactions, coupled with entropic changes due to chain ordering at the cholesterol surface, make it very difficult to use unambiguous theoretical modeling to aid in the interpretation of experimental data or to discern the molecular bases for the observed thermodynamic properties of these systems.

The effect of cholesterol on the lateral organization of membranes and on the domain structure of simple lipid bilayers are topics of high interest in membrane biophysics. Investigations of the lateral organization of lipid-cholesterol bilayers using fluorescence microscopy to study phospholipid (Slotte, 1995; Keller and McConnell, 1999; Radhakrishnan and McConnell, 1999a,b) and sphingolipid (Radhakrishnan et al., 2000, 2001) monolayer films at an air-water interface reveal large (up to 10 µm) domains of cholesterol-poor regions surrounded by cholesterol-rich boundary regions ~1 µm in width. The size of the regions and the overall lateral organization within the monolayers varies sensitively with pressure and is suggestive of critical mixing at certain pressures. McConnell and co-workers have used a free energy model to describe the domains seen in monolayer mixtures. The model is based upon the formation of complexes of cholesterol and ordered lipid molecules (condensed complexes). Complexes of p lipid molecules and q cholesterol molecules form single molecular clusters of cholesterol plus ordered lipid. Under lateral pressure, complexes aggregate to form the observed domains. The rich and complex properties of lipid-cholesterol mixtures revealed by experiment arise from energetic and entropic interactions at the submolecular level. However, current experimental methodologies have insufficient resolution to offer insight into interactions at this level, so that the atomic level mechanisms, which drive the observed nano- and larger-scale properties, are not well understood.

It is well documented that cholesterol plays a major role in the regulation of the fluidity and mechanical stiffness of membranes (Needham et al., 1988). This ability of cholesterol to modulate the fluidity of membranes is thought to be important in many biological processes and medical conditions, including cell fusion (Nakanishi et al., 2001), development of Alzheimers disease (Yip et al., 2001; Chochina et al., 2001; Kawahara and Kuroda, 2001), activity of the sodium pump (Cornelius, 2001), phase separation leading to raft formation (Brown and London, 1998; Brown, 1998), modulation of oxygen transport across red blood cell membranes (Buchwald et al., 2000)), neurotransmitter potency (Sooksawate and Simmonds, 2001), and modulation of uptake of mycobacteria by macrophages (Gatfield and Pieters, 2000), to name just a few. In fact, the biological significance of cholesterol-induced changes in membrane fluidity is so vast and varied that we could find no single review article that encompasses this entire topic. Recent observations that HIV-1 particles interact with cell surfaces at the location of cholesterol-rich rafts (Nguyen and Hildreth, 2000; Liao et al., 2001) further underscore the importance of understanding lipid-cholesterol interactions at the microscopic level.

Using large scale computer simulations, it is possible to gain insights into complex membrane interactions at the atomic level. The field of lipid bilayer simulation has grown rapidly in the past several years (Chiu et al., 1995, 1999a,b; Pastor, 1994; Stouch, 1993; Damodaran and Merz, 1994; Huang et al., 1994; Feller et al., 1994, 1997; Tu et al., 1995; Venable et al., 1993; Egberts et al., 1994; Armen et al., 1998; Tieleman and Berendsen, 1996; Berger et al., 1997; Husslein et al., 1998; Smondryev and Berkowitz, 1999b; for reviews, see Feller, 2000; Forrest and Sansom, 2000; Pastor and Feller, 1996; Merz and Roux, 1996; Tieleman et al., 1997; Jakobsson, 1997). More recent simulations of bilayers now examine the effects of helices (Tieleman et al., 1997), membrane proteins (Shrivastava et al., 2000; Shrivastava and Sansom, 2000), and cholesterol.

Early atomic level simulations of relatively small systems of lipids and cholesterol (Chol) were carried out using Monte Carlo (MC) (Scott, 1991) and molecular dynamics (MD) (Edholm and Nyberg, 1992; Robinson et al., 1995) methods. Improvements in models and in computing power have produced more details and longer and larger simulations in recent years. The simulations of Tu et al. (1998) examined the perturbing effect of 12.5% Chol on a dipalmitoyl phosphatidyl choline (DPPC) bilayer. They find relatively little change in bilayer structure and order produced by the addition of cholesterol. This finding differs from our results, reported in this paper, although we will argue that for cholesterol concentrations between ~8% to 12% the bilayers may be near a phase boundary between a fluid phase similar to a pure fluid DPPC bilayer, and a "liquid ordered" phase. This means that the two states may be nearly equally likely to appear in simulations, and slight differences in parameters, system size, or boundary conditions could force the appearance of one phase over the other. Smondryev and Berkowitz (1999a) have studied DPPC bilayers at 11% Chol and 50% Chol concentration. They also observe the condensing effect of cholesterol, which increases with cholesterol concentration. They find decreases in the area per molecule in their simulations and increases in hydrocarbon chain order parameters. The area per molecule (counting both DPPC and Chol) of their DPPC-Chol bilayer at 12% Chol is ~55 Ų. They also find, in 2-ns simulations, that cholesterol preferentially forms hydrogen bonds with DPPC carbonyl oxygens and one phosphate oxygen. Smondryev and Berkowitz (2000) have also recently carried out simulations of bilayers of dimyristoyl phosphatidylcholine (DMPC) bilayers containing cholesterol and two cholesterol analogs, ergosterol and lanosterol. They find that the three sterols have about the same effect on the chain order parameters, but they observe differences in the average depth of the three sterols, as well as in their hydrogen-bonding patterns. Róg and Pasenkiewiscz-Gierula (Pasenkiewiscz-Gierula et al., 1997; Róg and Pasenkiewicz-Gierula, 2001) have recently published a 15-ns simulation of a DMPC bilayer containing 22% Chol. They find that the condensing effect of cholesterol is to reduce the area per membrane molecule (counting both DMPC and Chol) from ~60 Ų to 53 Å² upon addition of this concentration of cholesterol. In two papers, they also examine the effect of cholesterol on polar and hydrocarbon regions, respectively. They find considerable hydrogen bonding between the cholesterol hydroxyl and DMPC choline methyls, with lesser bonding to DMPC carbonyl and phosphate oxygens. Recently, our group has carried out simulations of DPPC bilayers containing cholesterol and palmitoyl-oleyl phosphatidylcholine bilayers with cholesterol (Chiu et al., 2001a,b). These studies suggested that 1) cholesterol prefers to interact with lipid molecules, avoiding direct contact with other cholesterol molecules even at 1:1 concentrations, and forming hydrogen bonds with lipid carbonyl or (less often) phosphate oxygen atoms; 2) a single cholesterol molecule is able to partially or fully (depending on concentration) immobilize trans-gauche isomerism in the plateau region of the chains (carbons 3-12 for DPPC); 3) cholesterol molecules are well covered by close contact with lipid chains for both DPPC and palmitoyl-oleyl phosphatidylcholine.

Whereas the above studies considered a variety of lipid-Chol bilayers of varying composition under varying types of boundary conditions and with different force fields, a systematic investigation of lipid-Chol bilayers of varying composition, using one force field and one set of boundary conditions, has not yet been done. The purpose of this paper is to fill in this gap by presenting results of additional simulations of DPPC-Chol bilayers at different concentrations. These simulations, coupled with extended simulations of systems previously published by our group, will allow us to complete a full data set of lipid-cholesterol atomic level interactions, over the entire relevant range of cholesterol concentrations. By also extending the simulations to at least 5 ns, preceded by an extensive equilibration period including MD and configurational bias Monte Carlo (CBMC), we have observed new details of lipid-Chol interactions that were not apparent in our earlier runs.

	SIMULATION METHODS AND PROCEDURES

TOP ABSTRACT INTRODUCTION SIMULATION METHODS AND... RESULTS DISCUSSION REFERENCES

We constructed bilayer patches of DPPC-Chol for the 1:1, 2:1, 3:1, and 4:1 systems from previously constructed DPPC bilayers by replacing selected lipid molecules with cholesterol molecules. The replacement was done using DPPC bilayers that were in an intermediate state between fluid and gel, to ease the perturbation effects of the replacement of lipids by Chol molecules. Bilayer patches with DPPC:Chol ratios of 7:1, 8:1, and 11.5:1 were built successively from the equilibrated 4:1 system by replacing Chol molecules with DPPC molecules. The 24:1 system was constructed from a fluid DPPC bilayer by replacing 4 of 100 DPPC molecules with Chol molecules. Table 1 lists the bilayers that were constructed. The 1:1, 2:1, and 47:3 DPPC:Chol systems are extensions of earlier simulations, described in previous publications (Chiu et al., 2001a,b).

These simulations have all been equilibrated and continued to give trajectories of longer duration. In all cases the area per molecule evolved in time from an initial value, over a period of approximately a nanosecond, to the final values reported here. As Fig. 1 shows, the final values of the areas per molecule then remained stable for the 5-ns production phase of the simulations. Our general strategy for placing cholesterols in all of the systems was to spread them roughly uniformly over the bilayer patch with the same number in each leaflet. In this approach we avoid a priori assumptions about the nature of phase separations or cholesterol associations in bilayers. However, by studying a range of DPPC:Chol ratios from 24:1 up to 1:1, we in fact sample a comprehensive range of possible lipid-cholesterol configurations. The simulations are not long enough to directly observe any phase separations, but as we report below, by analyzing DPPC-Chol and Chol-Chol radial distribution functions we gain information about the interactions that can lead eventually to large scale phase transformations.

View larger version (36K): [in this window] [in a new window]   FIGURE 1   Plot of area per molecule versus time for all simulations, showing the 5 ns of data after each system had attained full equilibration. In all cases the final equilibrated areas per molecule differed from the starting areas per molecule before equilibration. Trajectories were interrupted at each 1000 ps for CBMC moves. Labels are (a) 1-1 DPPC:Chol; (b) 2-1 DPPC:Chol; (c) 3-1 DPPC:Chol; (d) 4-1 DPPC:Chol; (e) 7-1 DPPC:Chol; (f) 8-1 DPPC:Chol; (g) 11.5:1 DPPC:Chol; (h) 16-1 DPPC:Chol; (i) 24-1 DPPC:Chol; (j) pure DPPC. Dotted line is experimental value for area per molecule for pure DPPC from Nagle and Tristram-Nagle (2000).

All of the bilayers were initially energy minimized to remove bad steric contacts. Then over 30 cycles consisting of 100 ps of MD followed by 20,000 CBMC steps (Seipmann and Frenkel, 1992) were run on each system using a temperature of 323 K. The MD was run in all cases under constraints of constant pressure (1 atm) normal to the membrane, constant surface tension of  = 80 dyn/cm (40 dyn/cm per leaflet) in the membrane plane, and constant temperature of 323 K. There is at present no consensus on what is the appropriate surface tension to apply to a simulated membrane patch of a given size. In this study, we determined that a surface tension of 80 dyn/cm with our force field parameters gives a surface area of 63 Å²/phospholipid for a pure DPPC bilayer, in agreement with experiment (Nagle and Tristram-Nagle, 2000). Therefore, we elected to use 80 dyn/cm as the surface tension boundary condition in this study. Temperature coupling was applied separately to lipid and water components.

For the CBMC trial moves, one of the two separate DPPC chains (sn-1, sn-2) was selected at random on a randomly selected molecule for regrowth. A bond was selected at random on the chosen chain, and regrowth was attempted from that point toward the end of the chain. No attempts were made to regrow entire lipid molecules, and the CBMC procedure was confined to DPPC hydrocarbon chains only. The CBMC procedure consisted of 120 sampling trial moves at each chain position, starting from the randomly chosen carbon. During CBMC, no lateral diffusive moves were done, as the procedure was intended for lipid chain conformational equilibration only. In earlier work we have described our CBMC procedure fully (Scott et al., 1998). MD boundary constraints used the weak coupling method with a coupling constant of 0.2 ps for temperature and 4.0 ps for pressure and surface tension scaling. Discarding the first few picoseconds to allow for equilibration after CBMC, this procedure was repeated to give a total of three separate continuous trajectories, which were then used for averaging for each system.

Details of force fields used are given in our earlier papers (Chiu et al., 1999a,b,c). Noteworthy features include: partial charges on lipid head groups were determined by a Gaussian calculation (Chiu et al., 1995), united atoms were used for hydrocarbon groups on linear chains with van der Waals parameters determined to fit specific volume data for liquid hydrocarbons (Chiu et al., 1999c), and explicit hydrogens were used for hydrocarbon groups on cholesterol rings. For all MD runs reported here we used the GROMACS (Berendsen et al., 1995; Lindahl et al., 2001) molecular simulation code package. For computing torsion angle motions around lipid and cholesterol chain saturated bonds, third neighbor 6 to 12 interactions were replaced with the dihedral potential function due to Ryckaert and Bellmans (1978). For the hydrocarbon chain interactions between atoms on different chains, we used a modified set of parameters developed by our group by fitting to densities of a series of alkanes, and, for the C==C double bond, to 5-decene (Chiu et al., 1999c). Group-based spherical cutoffs of 20 Å for both electrostatics and 6 to 12 forces were used. Neutral charge groups were used. We have discussed the rationale for these choices in an earlier paper (Chiu et al., 1999a).

	RESULTS

TOP ABSTRACT INTRODUCTION SIMULATION METHODS AND... RESULTS DISCUSSION REFERENCES

The evolution of the area per membrane molecule over the last 5 ns of the simulations is shown in Fig. 1. Not shown on this plot is the initial evolution of each system from its starting area per molecule. The area plots for each system in Fig. 1 indicate the level of stability in lateral area reached by each simulation after equilibration, and the areas reached by each system represent predictions of the simulations. The plot also shows the condensing effect of cholesterol, as the area per molecule decreases with increasing cholesterol concentration. The least stable systems are the 7:1 and 8:1 DPPC:Chol bilayers (Fig. 1, e and f). Interestingly, the 11.5:1 system (Fig. 1 g) has become more highly ordered than other simulations at nearby Chol concentrations. This system is, like all others, stable in area over the 5-ns calculated trajectory. It is likely that the appearance of an highly ordered state in this concentration region is indicative of possible coexistence on a larger scale of ordered and fluid states, in that both fluid and ordered states are nearly equally likely to evolve in a simulation. Whereas current simulation time and length scales are inadequate to fully test this conjecture, analysis of trends in radial distribution functions reported below lends support to the concept. In a region of the lipid-cholesterol phase diagram near a phase boundary, one expects that fluctuations in area per molecule will be greater and equilibration to be slower.

Fig. 2 is a plot of the steady-state area per molecule (from Fig. 1) versus cholesterol concentration. Data points for cholesterol concentrations between 50% and ~12% are very nearly co-linear. At lower cholesterol concentration the area per molecule rises above the line drawn between the higher concentration points. The linear portion of the plot suggests a simple relation between the area per cholesterol, the area per DPPC, and the overall area per molecule over the concentration range between 0.1 and 0.5:

(1)

in which A_Mol is the overall area per molecule for given cholesterol and DPPC concentrations, X_Chol and X_DPPC, respectively. At the highest cholesterol concentrations, the area per cholesterol, A_Chol and the area per DPPC, A_DPPC, are independent of the cholesterol concentrations. If we suppose that this independence extends over the entire linear region of Fig. 2, then, because X_DPPC + X_Chol = 1, we can extract the areas per DPPC and Chol from the slope and intercept of the line. We find that A_DPPC = 50.7 Å² and A_Chol = 22.3 Å². This gives an area per DPPC-Chol pair of ~73 Ų, which is consistent with the value of 78 Å² found by Smondryev and Berkowitz (2000). The fact that these values are similar for two simulations performed in different ensembles using different force fields is indicative of the growing reliability and consistency of simulations of lipid bilayers.

View larger version (11K): [in this window] [in a new window]   FIGURE 2   Plot of area per molecule versus cholesterol concentration with linear regression line using areas from Fig. 1.

We note that the area per cholesterol deduced from Fig. 2 is substantially smaller than the crystalline cholesterol molecular area of ~37 Ų (Shieh et al., 1981). The cholesterol cross-sectional area in the bilayer membrane represents a steric minimal molecular cross-section, and is expected to be lower than the area per molecule in a cholesterol crystal, where rigid cholesterol molecules cannot fill space as efficiently as the flexible lipid chains allow to happen in a bilayer. The lower value is similar to that measured by Rothman and Engleman (1972) by dipping plasticine-filled Corey-Pauling-Koltun models in water. Although low-tech, this method is quite effective in determining a theoretical minimum cross-section. The measured cross-section area profile varied between ~12 Ų near the hydroxyl end and the tail chain to between 20 to 28 Å² in the fused ring region. We hypothesize that our low effective cross-sectional area, which is near the theoretical minimum shown by Rothman and Engelman, arises from the ability of flexible acyl chains to pack very closely around individual cholesterol molecules, enabling tighter packing than in a cholesterol crystal. Fig. 3 illustrates typical lipid-cholesterol packing configurations. As Fig. 3, A and B illustrates, for both low and high cholesterol concentrations, lipid chains are able to assume conformations that allow for maximal coverage of the surface of the cholesterol molecule. For low cholesterol concentrations chains that are not in direct contact with cholesterols assume more fluid-like conformations. At the highest (1:1) concentration, we also observe some direct contact between cholesterols.

View larger version (34K): [in this window] [in a new window]   FIGURE 3   Snapshots of lipid-cholesterol molecular clusters: (A) from 1:1 simulation; (B) from 8:1 simulation. Cholesterols are shown as spacefilled molecules, whereas DPPC molecules are shown as stick molecules

Fig. 4 is a plot of calculated order parameter profiles for all of the systems simulated. The plots represent averages over both sn-1 and sn-2 chains. With the exception of the 11.5:1 system, the profiles follow the same progression observed for perdeuterated palmitoyl-oleyl phosphatidylcholine:Chol bilayers observed by Lafleur et al. (1990). The profiles for the 1:1 and 2:1 systems agree well with data of Sankaram and Thompson (1990), as we pointed out in our earlier paper describing these systems (Chiu et al., 2001a). To further quantify the nature of lipid-cholesterol associations in bilayers we have calculated radial distribution functions (rdf) for various atom pairs on different molecules. The rdf between atom x and atom y is defined as the average over all x atoms in the system of the distance from an x atom to each other y atom. The distances between designated atoms on different molecules were binned, and the resulting rdf was normalized by dividing by 4r²dr in which r is the midbin distance variable and dr is the bin width, set at 0.5 Å. Fig. 5 is a plot of the distribution of nearest cholesterol distance to each carbonyl oxygen. The figure shows a sharp peak at ~2.5 Å for all concentrations with the peak height diminishing with cholesterol concentration for concentrations below 1:1. The primary peak is a direct measure of hydrogen bonding between cholesterol hydroxyl groups and DPPC carbonyl oxygens, and shows that this bonding is strong from 50% cholesterol down to 12.5% cholesterol. The overall shape of the curves in Fig. 5 is a measure of the range of the ordering induced by the carbonyl O-Chol OH affinity. The curves are liquid-like in structure with weaker secondary peaks and almost no structure beyond 5 Å in all of the systems. The fact that the first peak in the rdf plots is at the same position at all concentrations of cholesterol suggests that lipid-cholesterol membranes at all concentration ratios have the same structure at the nanoscale (the scale of lipid-cholesterol pairs). Fig. 6 is a plot of an rdf between DPPC phosphate oxygens and cholesterol O atoms (taking the nearest phosphate O on each DPPC molecule). In this plot the strong peak at 2.5 Å for the 1:1 system is due to competition between phosphate and carbonyl oxygens for bonding with the cholesterol hydroxyl group. This competition is only seen at the 50% DPPC:Chol system among the simulations we ran. For cholesterol concentrations below 50% comparison of Figs. 5 and 6 shows that cholesterol exhibits a clear preference for hydrogen bonding with the DPPC carbonyl oxygens over the phosphate oxygens.

View larger version (24K): [in this window] [in a new window]   FIGURE 4   Plot of order parameter profiles for the systems simulated.

View larger version (18K): [in this window] [in a new window]   FIGURE 5   Plot of rdf between DPPC carbonyl O - cholesterol O for all simulations. (A) Higher cholesterol concentrations; (B) lower cholesterol concentrations.

View larger version (20K): [in this window] [in a new window]   FIGURE 6   Plot of rdf between DPPC phosphate O - cholesterol O. (A) Higher cholesterol concentrations; (B) lower cholesterol concentrations.

As a check that the calculated rdfs between specific atoms reveals important details not otherwise available, we have also calculated rdfs between the center of mass of different molecules. Fig. 7 shows plots of rdfs between the centers of mass of DPPC and Chol molecules. The center-of-mass rdf for pure lipid in the fluid phase shows essentially no structure. At the various cholesterol concentrations much less clearly defined structure is seen in the center-of-mass rdfs in Fig. 7 than is seen in the earlier rdfs involving potential hydrogen-bonding partners shown in Figs. 5 and 6. Comparison of Fig. 7 with Figs. 5 and 6 underscores the fact that correlations between large and flexible molecules depend strongly on the atoms within the molecules used for the calculations. By choosing atoms likely to be involved in strong electrostatic (hydrogen bonding) interactions we are able to observe structural properties of the systems that might otherwise go unnoticed. To address the issue of direct cholesterol-cholesterol interactions we have calculated the rdf between pairs of cholesterol oxygens, and Fig. 8 shows the resulting plots of this rdf versus r. The plots in Fig. 8 show that even at low cholesterol concentrations there is a tendency for pairs of cholesterols to be located at separations of 5, 10, and 15 Å. To further investigate this surprising result we have examined the time evolution of the Chol O-Chol O rdfs. These are plotted in Fig. 9 (A-C) for the 1:1, 3:1, and 7:1 DPPC:Chol runs, respectively. The 7:1 and 3:1 data clearly show that the peaks at multiples of 5 Å begin to appear during the latter stages of the simulations. This implies strongly that the patterns we see in the radial distribution function emerge from the properties of the system rather than being a function of the relative distribution of the cholesterol and lipid at the beginning of the simulation. For the 1:1 system the DPPC-Chol packing is already quite close so there is little change over time for this simulation. For the 3:1 and 7:1 systems the data suggest the emergence of a lateral organizational structure consisting of dynamical "building blocks" of cholesterol plus lipid, which consist of a single cholesterol plus a single lipid closely associated by a hydrogen bond primarily between the DPPC carbonyl O and the Chol hydroxyl. For two Chol molecules to achieve a 5-Å separation between their O atoms, two building blocks can pack this close but not directly hydrogen bond to each other. The number of molecules available for participating in hydrogen bonding can be directly calculated from the area under the initial peak in each of the rdfs. Fig. 10 is a plot of these data, for water and cholesterol oxygen, water and DPPC carbonyl oxygen, and water and DPPC phosphate oxygen atoms, respectively. The plots show that nearly two waters participate in the first hydration shell for each cholesterol, whereas slightly less than one DPPC molecule is in the first Chol coordination shell. By comparison, in their simulation of DMPC-Chol at 22% cholesterol, Pasenkiewiscz-Gierula et al. (2000) find a similar association between DMPC and cholesterol, and 1.1 ± 0.1 bonds between cholesterol and water in their system. However Pasenkiewiscz-Gierula et al. (2000) use both a distance bond angle criterion for hydrogen bond determination, whereas we use only a distance criterion.

View larger version (23K): [in this window] [in a new window]   FIGURE 7   Plots of rdf between the centers of mass of DPPC and cholesterol for a subset of the simulations with a similar plot for pure DPPC (between centers of mass of DPPC molecules) shown for reference.

View larger version (21K): [in this window] [in a new window]   FIGURE 8   Plot of rdf between cholesterol O - cholesterol O. (A) Higher cholesterol concentrations; (B) lower cholesterol concentrations averaged over the entire 5 ns of data taking.

View larger version (17K): [in this window] [in a new window]   FIGURE 9   Plots of rdf between cholesterol O - cholesterol O decomposed into separate plots for different periods during the simulation. (A) 1:1 DPPC:Chol simulation; (B) 3:1 DPPC:Chol simulation; (C) 7:1 DPPC:Chol simulation.

View larger version (16K): [in this window] [in a new window]   FIGURE 10   Plots of average numbers of molecules in the first coordination shell of cholesterol oxygens. Contributions from water oxygens, DPPC carbonyl oxygens, and DPPC-phosphate oxygens are shown. All data are calculated from integrals of the radial distribution functions.

Fig. 11 shows a snapshot of two Chol molecules in this category from the 7:1 simulation. The Chol hydroxyl groups are hydrogen bonded to different lipids but are also quite close to each other. For a larger scale picture of lateral organization, which is emerging in the simulations, Fig. 12 shows a snapshot of the Chol molecules as they have evolved at the end of the 3:1 simulation. For most Chol oxygens in this simulation the nearest neighbor Chol O atom is close to a multiple of 5 Å. Whereas the distribution of nearest neighbor distances is regular, the lateral arrangement of the Chol molecules in Fig. 12 is not lattice-like on the 5-ns timescales of our simulations. Rather, Chol molecules tend to associate in a linear fashion, forming curvilinear microdomains across the simulation cell. This pattern appears in all of our simulations at 4:1 and higher DPPC:Chol ratios.

View larger version (73K): [in this window] [in a new window]   FIGURE 11   Snapshot from 7:1 DPPC:Chol system showing two closely associated cholesterol molecules (space-filled molecules) with surrounding lipids (stick molecules).

View larger version (70K): [in this window] [in a new window]   FIGURE 12   Snapshot of one leaflet from 3:1 DPPC:Chol system with water removed, illustrating lateral arrangements of cholesterols. Cholesterols are shown as space-filled molecules, whereas DPPC molecules are shown as stick molecules. Note the three sets of pairs of cholesterols with tightly associated lipid molecules.

	DISCUSSION

TOP ABSTRACT INTRODUCTION SIMULATION METHODS AND... RESULTS DISCUSSION REFERENCES

The plot in Fig. 2 shows that a DPPC-cholesterol bilayer is in a liquid-ordered state with essentially the same lateral area per lipid and cholesterol molecules from the highest (50%) cholesterol concentration down to approximately 10% cholesterol concentration. This suggests that a single cholesterol can on average partially or fully order, by direct and indirect contact, as many as eight or nine lipid molecules. Our area data also show the surprisingly low cross-sectional area profile presented by cholesterol in a lipid bilayer. The reduced cross-section is seen to be a consequence of the "wetting" of the surface of cholesterol molecules by lipid chains, allowing much closer molecular contact than in a cholesterol crystal. Theoretical models based on a larger cross-section for cholesterol in a bilayer are not consistent with this result (Kessel et al. 2001). Simulations of DMPC bilayers containing 22% Chol by Róg and Pasenkiewicz-Gierula (RP) (2001) show that this system reaches a steady-state area of 53 Å²/mol, compared with a value of ~60 Ų for pure DMPC (Pasenkiewiscz-Gierula et al., 1997). This reduction in area is not as great as the reduction we obtain. Possible reasons for the difference include both simulation artifacts and real systemic differences. The simulations of RP involved insertion of Chol into a fluid bilayer by replacement of DMPC molecules, whereas our procedure was to carry out these replacements in a bilayer in an intermediate state between fluid and ordered states. We believe our procedure will lead to better equilibration because of the relative ease with which chains can disorder, compared with an expected longer timescales for spontaneous ordering. Other differences between the simulations we present here and those of RP involve simulation size (RP use 72 lipids with 23 waters/lipid, compared with at least 100 lipids with 32 waters/lipid for our work), interaction cutoff distance (RP use 12 Å whereas we use a cutoff of 20 Å for both 6 to 12 and electrostatic interactions), and boundary conditions (we use constant normal pressure and lateral surface tension whereas RP use constant isotropic pressure). An additional factor is our use of configurational bias Monte Carlo both for more complete equilibration and for enhanced sampling. There are also likely to be real systemic differences due to the different chain lengths of DMPC and DPPC. The hydrophobic length of the DMPC chains is similar to the length of a cholesterol molecule, whereas that of DPPC chains is longer. DPPC chains are for this reason better able to fold under the cholesterol tail chain, which may facilitate closer packing of DPPC-Chol compared with DMPC-Chol.

Interactions between lipid and cholesterol molecules in bilayers involve direct electrostatic forces, induced dipole forces (6-12 potentials), solvent interactions, and entropic contributions from the lipid acyl chains. The free energy of a lipid bilayer containing cholesterol is reduced by hydrogen bonding interactions between lipids and cholesterol and by the increased penetration of water into the bilayer near cholesterol hydroxyls. It is also reduced by strong 6-12 interactions between lipid chains and cholesterol molecules. Competing with these interactions is the increase in free energy due to acyl chain entropy loss caused by the steric interactions with cholesterols. The complex interplay between these mechanisms determines the lateral organization of the membrane and its thermophysical properties. The appearance of regular peaks in the Chol O-Chol O rdf plots in the latter stages of the simulations suggests a complex ordering is in progress involving both lipids and cholesterols. The ordering involves both single-molecule states, i.e., chain ordering, and lateral positional ordering of the molecules within the leaflets, i.e., the formation of molecular clusters, or complexes. The ordering is probably driven by a number of factors including hydrogen bonding between Chol hydroxyl and DPPC carbonyl O, and to a lesser extent DPPC phosphate O, atoms, and also including entropically driven attractive interactions that can force complexes of Chol plus hydrogen-bonded DPPC molecules together. From the current simulations it is unclear whether the observed structural effects are transient in time or are dependent on the size of the system. We are currently carrying out extended simulations on systems four times larger to examine this issue.

Our results provide a possible molecular basis for the condensed complex model of Radhakrishnan and McConnell (Keller and McConnell, 1999; Radhakrishnan and McConnell, 1999a,b) and possibly for suggestions of "superlattice" formation in lipid:Chol bilayers (Chong, 1994; Chong et al., 1996). In the simulations, each cholesterol preferentially is hydrogen bonded to at least one lipid molecule. Evidently this two-molecule complex can associate with other similar complexes via the mechanisms mentioned above, giving rise to nanoscale "building blocks." The ultimate size of the blocks should depend on the Chol concentration, but it is possible that these blocks can, at some concentrations, "tile" the membrane. Although it is not possible to determine the cooperativity of interactions between complexes (Anderson and McConnell, 2001) on the current molecular dynamics simulation timescale, we note that the cholesterol mole fraction below which condensed complexes do not form (10%) is precisely the mole fraction below which our simulated phase change occurs, as shown in Fig. 2. Propagation of these simulations to larger length and time scales is clearly warranted and is currently in progress.