We show, via molecular simulations, that not only does cholesterol induce a lipid order, but the lipid order also enhances cholesterol localization within the lipid leaflets. Therefore, there is a strong interdependence between these two phenomena. In the ordered phase, cholesterol molecules are predominantly present in the bilayer leaflets and orient themselves parallel to the bilayer normal. In the disordered phase, cholesterol molecules are mainly present near the center of the bilayer at the midplane region and are oriented orthogonal to the bilayer normal. At the melting temperature of the lipid bilayers, cholesterol concentration in the leaflets and the bilayer midplane is equal. This result suggests that the localization of cholesterol in the lipid bilayers is mainly dictated by the degree of ordering of the lipid bilayer. We validate our findings on 18 different lipid bilayer systems, obtained from three different phospholipid bilayers with varying concentrations of cholesterol. To cover a large temperature range in simulations, we employ the Dry Martini force field. We demonstrate that the Dry and the Wet Martini (with polarizable water) force fields produce comparable results.

Cell plasma membranes, the semipermeable interfaces that allow cells to interact with the extracellular environment, are made of proteins, lipids, sterols, and other molecules.1,2 The lateral mobility and function of membrane proteins are understood to be dictated by the structure and dynamics of the membrane.3,4 Cholesterol, the main sterol in the mammalian cell membrane, affects physical properties, such as bending modulus,5 fluidity,6 permeability,7 and phase behavior of the membrane,8 and it is understood to play a role in the formation of ordered domains.9–11 It is well-established that cholesterol solubility in the bilayers is a function of the molecular characteristics of the lipids. For instance, lipid head groups shield cholesterol molecules from water, and thus, the ones with large head groups, like sphingomyelins, display higher solubility of cholesterol.12,13 Saturated lipids are reported to enhance cholesterol solubility as opposed to mono/un-saturated lipids.14–18 

Cholesterol is well known to affect the phase behavior of lipid bilayers. To a first approximation, pure phospholipid bilayers predominantly display two main phases: ordered and disordered, as distinguished by their spatial order and chain configurational order in the bilayers.10,12 In the ordered phase, both the spatial and chain configurational order are high. While in the disordered domain, both these orders are low. The pure phospholipid bilayers undergo a first order phase transition from the ordered to the disordered phase as the temperature is increased. Cholesterol molecules induce conformational ordering of the lipids in the disordered phase. Even though at high concentrations, cholesterol molecules perturb lipid packing in the ordered phase, their overall effect is to favor the formation of the ordered phase by enhancing the chain configurational order.19–23 

In our previous work, we showed via coarse-grained simulations that cholesterol distributes between the leaflets of equimolar asymmetric lipid bilayers in a manner so as to relieve the compressive and tensile mechanical stresses that arise due to differences in the conformational ordering of the leaflets.24 Thus, the conformational ordering of phospholipids likely affects cholesterol distribution. This is also evidenced from neutron scattering studies that have shown that in phospholipid bilayers comprising polyunsaturated fatty acid (PUFA), cholesterol is sequestered near the bilayer center rather than in the leaflets, due to the highly disordered structure of PUFA.10,14,25 This interdependence of conformational order and cholesterol distribution within the bilayer has not been established in a quantitative manner. In this work, we investigate the relationship between the distribution of cholesterol molecules within the leaflets and the midplane region of the bilayers and the conformational ordering of the lipids in various symmetrical bilayers. To systematically study this relationship, we simulate the lipid bilayers at small temperature increments to impart a conformational disorder. While cell membranes include membrane proteins, carbohydrates, sterols as well as many types of phospholipids, it is complicated to study such a system, and, therefore, we restrict ourselves to studying lipid bilayers of only one or two kinds of phospholipids with cholesterols.

To study the large temperature range in simulations, in some cases, we need to access unphysical temperatures. Therefore, we employ an implicit solvent coarse-grained lipid model, termed as Dry Martini force field.26 However, an important concern is whether the Dry Martini force field is amenable to capturing the structural properties of lipid bilayers and the distribution of cholesterol that are observed in more detailed simulations of lipid bilayers with explicit water. Therefore, as a first step, we have performed a systematic comparison of the results obtained from Dry Martini with an improved coarse-grained Wet Martini force field with a polarizable water model as well as known experimental results.27 We show that the results obtained from the two force fields are comparable for the purpose of this study. For all the studied lipid bilayer systems, we find that at the melting temperature, the concentration of cholesterol in the leaflets and the midplane is the same. This result, in support of our previous work,24 shows that while experiments have identified numerous factors that affect cholesterol localization in the bilayers, the main effect boils down to how these factors change the ordering of the lipids in the bilayers, which ultimately dictates cholesterol localization.

For comparing the results obtained from Dry Martini version 2.1 and Wet Martini version 2.2 (study-I), we have studied three symmetric lipid bilayers, namely, DOPC (1, 2-dioleoyl-sn-glycero-3-phosphocholine)/CHOL (Cholesterol), 16:0SM (d18:1/16:0 sphingomyelin)/CHOL, and DOPC/16:0SM/CHOL for a temperature range of 290 to 350 K and cholesterol concentrations varying from 0% to 50% mole-ratio. DOPC in an unsaturated phospholipid that is in the disordered phase at room temperature. SM, on the other hand, is a saturated phospholipid that has a higher conformational order than DOPC. The system of a saturated lipid, cholesterol, and an unsaturated lipid is a simple mimic of the mammalian membrane and has been used in many biophysical studies to mimic the mammalian cell membrane.28–32 

After establishing the suitability of the Dry Martini force field for studying the structural properties of lipid bilayers at different cholesterol concentrations, we examine the spatial distribution of cholesterol in DOPC/CHOL, 20:0SM/CHOL, and DOPC/20:0SM/CHOL lipid bilayers over a wider range of temperatures and cholesterol concentrations using the Dry Martini force field (study-II). The temperature of the system and the concentration of cholesterol in study-II are varied from 210 to 400 K and 0% to 50% mole-ratio, respectively. We have chosen temperatures beyond the range that is biologically relevant so as to clearly identify phase transitions in all the studied systems. For instance, the melting temperature of a pure DOPC bilayer is quite low (∼253 K)33,34 while that of a pure 20:0SM bilayer is relatively high (∼319 K).35 

In both studies-I and -II, the simulation box size is 125 Å in all directions and periodic boundary conditions are applied. The total number of lipids (two-chain lipids + cholesterol) for all the systems is kept at 510 ± 2. In Wet Martini simulations, there are 10 850 ± 40 water molecules. Figure S1(a) in the supplementary material54 shows a snapshot of the equilibrium configuration of the DOPC/16:0 SM/CHOL lipid bilayer at 290 K. Figure S1(b) in the supplementary material54 shows topologies of the lipids studied along with the name and type of each bead. Initial configurations of the bilayer systems are constructed using the insane.py script developed by Marrink and co-workers.36 The simulations are performed using the GROMACS/5.1.2 molecular simulation package.37 To check for the potential finite size effects, in our previous work, we performed simulations with a larger box size of 250 × 250 × 125 Å3 and confirmed that the smaller box size is sufficient to capture the physical characteristics of the studied lipid bilayers.24 

The initial configuration is energy minimized using the steepest descent with the maximum force tolerance of 1000.0 KJ/mole/nm. To pre-equilibrate the bilayer structure, a 25 ns canonical ensemble (constant number of particles N, constant volume V, and constant temperature T) simulation followed by a 25 ns isothermal-isobaric ensemble (constant N, pressure P, and T) simulation is performed at each temperature. The time step for the pre-equilibration simulations is 10 fs. After pre-equilibration, a production run of 1 μs is performed at each temperature in the NPT ensemble, and the data are collected for the last 500 ns. A 500 ns of equilibration time ensures that the system is in equilibrium even at the lowest temperatures studied by us.

For the systems where the Dry Martini force field is used, simulation parameters are taken from the Arnarez et al. study that parametrized this force field based on the original Martini force field.26 The Wet Martini force field with the polarizable water model parameters are taken from the work of Michalowsky et al.27 Table S1 in the supplementary material54 summarizes the simulation parameters used in all the simulations of lipid bilayers studied in this work. Tieleman and co-workers have shown that the free energy barrier associated with the cholesterol flip-flop is lower than that found in fully atomistic simulations.20 While a lower free energy barrier increases the kinetics, it does not affect the equilibrium distribution of cholesterol in the leaflets and midplane regions. We have employed a coarse-grained model to efficiently sample the long time- and length-scales needed in these simulations. However, a coarse-grained model may not capture some specific chemical interactions between lipids and water. Furthermore, we acknowledge that apart from the Martini force field, there are other coarse-grained models available as well to study the role of cholesterol in phospholipid bilayers.38 

Area per lipid (apl) is the ensemble-averaged area of the plane of the simulation box in which the bilayer resides divided by the total number of two-chain lipids + cholesterol molecules in each leaflet.39 Bilayer thickness (DP-P) is calculated as the ensemble-averaged distance between the phosphate (PO4) beads of lipids in each leaflet.26 Orientation of a cholesterol molecule within a bilayer is defined as the angle between vector, a ¯ [the vector joining the ROH bead to the C1 bead of cholesterol, see Fig. S1(b) in the supplementary material54 and vector z ¯ (the bilayer normal). To quantify the chain configurational order of lipids, we calculate the ensemble-averaged order parameter of the acyl chain of lipids (Schain) as1,24,40
S c h a i n = 1 2 ( 3 co s 2 θ 1 ) ,
(1)
where θ is the angle between each bond vector of the lipid beads and the bilayer normal (z axis). A large value of the Schain (∼1) indicates a large chain configurational order characteristic of the ordered phase,12,41 while a small value is a signature of the disordered phase. Translational ordering of lipids is quantified by calculating two-dimensional (2D) radial distribution functions [g2D(r)] of the C1B beads [Fig. S1(b) in the supplementary material54] that represent the center of mass of lipids.1 

To ensure that the Dry Martini force field is adequate for studying the structural properties and cholesterol distribution in lipid bilayers, we compare its performance with the Wet Martini force field as well as the known experimental results. While the Dry Martini force field has been parameterized on the Wet Martini force field only, it is important to ensure that their performance is comparable at different temperatures and compositions of lipid bilayers. Table I lists the different systems and conditions for which the two force fields are compared.

TABLE I.

Temperature and composition of the lipid bilayers that are compared using the Dry and Wet Martini force fields. In the ternary mixtures, DOPC and 16:0SM are in the 1:1 molar ratio.

Lipid bilayerNo. of simulationsTemperature (K)CHOL (mol. %)
DOPC/CHOL 7 T × 6% CHOL = 42 290–350 0–50 
16:0SM/CHOL 42 290–350 0–50 
DOPC/16:0SM/CHOL 42 290–350 0–50 
Lipid bilayerNo. of simulationsTemperature (K)CHOL (mol. %)
DOPC/CHOL 7 T × 6% CHOL = 42 290–350 0–50 
16:0SM/CHOL 42 290–350 0–50 
DOPC/16:0SM/CHOL 42 290–350 0–50 

Figure 1 compares the results obtained from the Dry and Wet Martini force fields for apl, DP-P, Schain, and g2D(r) of DOPC/16:0SM/CHOL bilayers with 30% mole of cholesterol at different temperatures. These results for DOPC/CHOL and 16:0SM/CHOL bilayers with 30% mole-ratio of cholesterol are shown in (Figs. S2 and S3 in the supplementary material54), respectively. The results obtained for other %mol. of CHOL using the Dry and Wet Martini force fields are also found to be similar. From this analysis, it is established that the results obtained from the Dry and the Wet Martini simulations are in good agreement with the studied lipid bilayers, the temperature range, and CHOL concentrations. The largest difference in the results is observed in the case of apl of the DOPC/CHOL bilayer [Fig. S2(a) in the supplementary material]54 where apl from the Dry Martini simulations is ∼0.04 nm2 lower than the Wet Martini results, but the rate of change of apl with the temperature is the same (Fig. S3 in the supplementary material54). The apl and DP-P values obtained for pure DOPC and pure 16:0 SM at different temperatures from the Dry and Wet Martini force fields are in excellent agreement with each other (Table S2, supplementary material54). These values are also in good agreement with the reported experimental,35,42,43 coarse-grained,24,44 and atomistic simulation45–47 results (see the text accompanying Table S2 in the supplementary material54). From these results, we conclude that the Dry Martini force field is adequate for studying the structural properties of these lipid bilayers at different temperatures.

FIG. 1.

Comparison of results obtained from the Dry and Wet Martini force fields for DOPC/16:0SM/CHOL lipid bilayers for temperatures ranging from 290 to 350 K, and a cholesterol concentration of 30% mole. (a) Area per lipid, apl; (b) bilayer thickness, DP-P; (c) chain order parameter, Schain; 2D RDFs of C1B-C1B beads of the two-chain lipids, g2D(r) obtained from the (d) Dry Martini; and (e) Wet Martini simulations. Successive g2D(r) have been shifted vertically by 0.4 units for the sake of clarity. Error bars in the figures are obtained by performing block averages over three equally sized blocks.

FIG. 1.

Comparison of results obtained from the Dry and Wet Martini force fields for DOPC/16:0SM/CHOL lipid bilayers for temperatures ranging from 290 to 350 K, and a cholesterol concentration of 30% mole. (a) Area per lipid, apl; (b) bilayer thickness, DP-P; (c) chain order parameter, Schain; 2D RDFs of C1B-C1B beads of the two-chain lipids, g2D(r) obtained from the (d) Dry Martini; and (e) Wet Martini simulations. Successive g2D(r) have been shifted vertically by 0.4 units for the sake of clarity. Error bars in the figures are obtained by performing block averages over three equally sized blocks.

Close modal

We employ the Dry Martini force field to perform a systematic investigation of how the domain/phase characteristics of lipid bilayers affect the distribution of cholesterol molecules within the bilayer. For this purpose, we have studied three types of lipid bilayers: DOPC/CHOL, 20:0 SM/CHOL, and DOPC/20:0 SM/CHOL with the cholesterol molar concentration ranging from 0 to 50%. At 50% mol. cholesterol concentration, we are already close to the cholesterol solubility limit in the lipid bilayers, and therefore, do not study higher concentrations. Each of these systems have been studied for a temperature range of 210–400 K. Overall, 360 different simulations are performed. First, we discuss how the phase transition temperature, Tm, of the lipid bilayers is calculated. Then, we illustrate the effect of cholesterol concentration on Tm, and finally, we examine the relationship between Tm and cholesterol localization and orientation within the bilayers.

Figure 2 shows apl, DP-P, Schain and g2D(r) of pure DOPC, pure 20:0SM, and 1:1 mixture of DOPC/20:0SM bilayers in the absence of any cholesterol, termed as pure phospholipid bilayers. It is observed that the 20:0SM and the 1:1 mixture of DOPC/20:0SM display the ordering of lipids at low temperatures, indicated by a large value of Schain. The ordering of lipids is accompanied by a reduction in apl and an increase in DP-P as the lipids get tightly packed. As the temperature increases, these bilayers undergo an ordered → disordered transition. The pure 20:0 SM undergoes a sharp transition, while the 1:1 DOPC/20:0 SM shows a gradual transition to the disordered phase. In contrast, a pure DOPC bilayer does not show any significant ordering even at the lowest temperatures studied, and therefore remains in the disordered phase. The transition temperature, Tm, of the bilayers is determined as the point of inflection of apl, DR-R, and Schain as a function of temperature.

FIG. 2.

(a) Area per lipid, apl; (b) bilayer thickness, DP-P; (c) chain configurational order, Schain; 2D RDFs of C1B-C1B beads, g2D(r) in (d) pure DOPC bilayer; (e) g2D(r) in pure 20:0SM bilayer; and (f) g2D(r) in 1:1 DOPC/20:0SM bilayer. Successive g2D(r) have been shifted vertically by 0.4 units for the sake of clarity.

FIG. 2.

(a) Area per lipid, apl; (b) bilayer thickness, DP-P; (c) chain configurational order, Schain; 2D RDFs of C1B-C1B beads, g2D(r) in (d) pure DOPC bilayer; (e) g2D(r) in pure 20:0SM bilayer; and (f) g2D(r) in 1:1 DOPC/20:0SM bilayer. Successive g2D(r) have been shifted vertically by 0.4 units for the sake of clarity.

Close modal

Figure S4 in the supplementary material54 shows a plot of temperature derivatives of these structural properties. Tm of pure 20:0SM is estimated to be 310 K, which agrees well with the experimentally reported value of 310–320 K.35  Tm of 1:1 mixture of the DOPC/20:0SM system is estimated to be 250 K. For pure DOPC, there is no point of inflection in the studied temperature range. For the 20:0SM and DOPC/20:0SM systems, g2D(r) are shown in Figs. 2(e) and 2(f), respectively. g2D(r) clearly shows a loss of long-range structure as the temperature increases beyond Tm. g2D(r) of DOPC does not show any long-range order. The behavior of g2D(r) and the absence of Tm in DOPC in the studied temperature range implies that this bilayer does not attain the ordered phase even at the lowest temperatures studied. It is observed in Figs. 2(a)2(c) that the properties of the binary mixture of the DOPC/20:0SM system are in between the two single-lipid bilayers.

Next, we study the effect of cholesterol concentration on the structural properties of the lipid bilayers and the phase transition temperature, Tm. Figure 3 shows the variation in the structural properties of the three studied bilayers with the cholesterol concentration at different temperatures. Figure 3(a) shows that apl decreases with the cholesterol concentration, indicating the well-known condensation effect of cholesterol.24,32,48 Figure 3(b) shows that at high temperatures, that is in the disordered phase, bilayers with higher cholesterol concentrations are thicker.16,24,49 However, in the case of 20:0SM/CHOL and DOPC/20:0SM/CHOL bilayers, in the ordered phase, the trend is reversed with the bilayer thickness decreasing with the cholesterol concentration. Similarly, for these two bilayers, it is observed that in the disordered phase (high temperatures), Schain increases with the cholesterol concentration. In the ordered phase, this trend is reversed with Schain decreasing with the cholesterol concentration [Fig. 3(c)]. The decrease in Schain in the ordered phase is due to the perturbation of lipid packing by cholesterol, as reported before.10,41

FIG. 3.

Structural properties of 20:0SM/CHOL bilayers as a function of temperature and cholesterol concentration. (a) Area per lipid, apl; (b) bilayer thickness, DP-P; and (c) chain order parameter, Schain. Error bars in all figures are estimated from the block averages of three equally sized blocks.

FIG. 3.

Structural properties of 20:0SM/CHOL bilayers as a function of temperature and cholesterol concentration. (a) Area per lipid, apl; (b) bilayer thickness, DP-P; and (c) chain order parameter, Schain. Error bars in all figures are estimated from the block averages of three equally sized blocks.

Close modal

Tm, estimated as the point of inflection temperature of Schain for 20:0SM/CHOL and DOPC/20:0SM/CHOL bilayers with different cholesterol concentrations are shown in Fig. S5 in the supplementary material54. Tm increases with cholesterol concentration, as has been reported before for the DPPC/CHOL bilayer.3 The DOPC/CHOL bilayers do not show any point of inflection temperature for any cholesterol concentration, and so do not manifest any ordered → disordered transition. The Tm values of all the studied systems are shown in Fig. 4.

FIG. 4.

Phase transition temperature, Tm, of the studied bilayers. The DOPC/CHOL bilayers do not have Tm within the temperature range studied. Lines are a guide to the eye.

FIG. 4.

Phase transition temperature, Tm, of the studied bilayers. The DOPC/CHOL bilayers do not have Tm within the temperature range studied. Lines are a guide to the eye.

Close modal

1. Cholesterol distribution within the bilayer

We have shown previously that in equimolar asymmetric lipid bilayers, cholesterol distributes between the leaflets so as to relieve the compressive and tensile mechanical stresses.24 These mechanical stresses arise due to differences in the conformational order of the lipids in the two leaflets. The interdependence of cholesterol distribution and lipid order should also manifest as a function of the phase behavior of the lipid bilayers, which is the focus of this work. In the current work, we relate the distribution of cholesterol in the bilayer to the phase of bilayers. For this purpose, we divide the bilayers into the leaflet and the midplane regions. As suggested for the DPPC/CHOL bilayer,50 the midplane region is defined as ∼14 Å from the center of the bilayers, and the remaining top and bottom parts form the leaflet region (Fig. S6, supplementary material54).

Figure 5 shows the fraction of cholesterol molecules in the leaflet and midplane regions of the three studied bilayers as a function of temperature and the molar concentration of cholesterol. It is observed that as the temperature increases, the cholesterol concentration decreases in the leaflets and increases in the midplane. The temperature when the concentration of cholesterol is equal in the leaflets and the midplane is termed as Tcross. It is observed that for a bilayer, Tcross increases with cholesterol concentration. It should be noted that the presence of cholesterol in the midplane has been reported in previous simulation studies as well.10,14,51–53

FIG. 5.

Fraction of cholesterol molecules in the leaflets and midplane regions for (a) (blue) DOPC/CHOL, (b) (red) 20:0SM/CHOL, and (c) (black) DOPC/20:0SM/CHOL lipid bilayers for different cholesterol concentrations. Error bars are estimated via block averages over three independent simulations. The temperature when the fraction of cholesterol is equal in the leaflets and the midplane is termed as Tcross.

FIG. 5.

Fraction of cholesterol molecules in the leaflets and midplane regions for (a) (blue) DOPC/CHOL, (b) (red) 20:0SM/CHOL, and (c) (black) DOPC/20:0SM/CHOL lipid bilayers for different cholesterol concentrations. Error bars are estimated via block averages over three independent simulations. The temperature when the fraction of cholesterol is equal in the leaflets and the midplane is termed as Tcross.

Close modal

Figure 4 reveals that Tm increases with cholesterol concentration, which highlights the role of cholesterol in imparting conformational order in the lipid bilayers. Figure 5, on the other hand, reveals that the change in the conformational order, induced via temperature, affects the cholesterol distribution. Thus, these two results illustrate the interdependence of cholesterol concentration and lipid order. In fact, we note that Tm = Tcross (Fig. 6). This implies that when the bilayer is in the disordered phase, that is T > Tm, the cholesterol molecules are predominantly in the midplane of the bilayer. In the ordered phase (T < Tm), the cholesterol molecules are mainly in the leaflets of the bilayer. The crossover in the region that dominates in the cholesterol concentration happens exactly at Tm. How does one rationalize this observation? In the ordered phase, the phospholipid chains are arranged parallel to each other, and so the cholesterol molecules prefer to align parallel to the phospholipids in the leaflets due to the favorable hydrophobic interactions. However, in the disordered phase, when the conformational order of the phospholipid chains is lost, cholesterol is unable to retain hydrophobic interactions by aligning with the chains. Thus, it is more favorable for cholesterol to get sequestered in the midplane region, which is less crowded and so entropically favorable. This result aligns with the previous reports that cholesterol molecules prefer to lie parallel to the plane of the lipid bilayer in the midplane region.41,50 We have also analyzed the orientation of cholesterol molecules as a function of their location in the lipid bilayer.

FIG. 6.

Relationship between the melting temperature of the lipid bilayers, Tm, and the temperature when the distribution of cholesterol molecules is equal in the leaflets and the midplane, Tcross. It is found that Tm = Tcross.

FIG. 6.

Relationship between the melting temperature of the lipid bilayers, Tm, and the temperature when the distribution of cholesterol molecules is equal in the leaflets and the midplane, Tcross. It is found that Tm = Tcross.

Close modal

Figure 7 shows the average orientation of a cholesterol molecule as a function of its location in the three studied lipid bilayers, namely, DOPC/CHOL, 20:0SM/CHOL, and DOPC/20:0SM/CHOL with 10–50% mole of cholesterol at three different temperatures, one above the phase transition temperature (Tm), one close to Tm, and one below Tm. For this calculation, we divide the bilayers into 50 slabs parallel to the XY plane. In each slab, the ensemble-averaged angle of the cholesterol molecules with the z axis (bilayer normal) is calculated. The reported angles correspond to the average over configurations of three equally sized blocks of simulations. Figure 7 shows that irrespective of the conditions (that is, temperature, type of the lipid bilayer, cholesterol concentration), the cholesterol molecules that are located in the leaflets orient with bilayer normal at an angle θ < 30° or θ > 150°, while, at the midplane, they orient with the angle 30° < θ < 50°. An important observation from Fig. 7 is that the cholesterol orientation changes almost linearly with the distance from the center of the bilayer. While we have studied three different lipid bilayers, one would surmise that our reported results are generally applicable to other lipid bilayer systems as well. The rationale for this assertion is that in the three systems studied by us, cholesterol distribution seems to be strongly dependent only on the conformational order of the lipids rather than on other chemical details of the lipids. Our previous work, wherein we studied 17 different asymmetric bilayers, also showed that irrespective of the bilayer type, the distribution of cholesterol molecules depended strongly on the lipid order.24 An important result of our work is that at Tm the concentration of cholesterol in the leaflets and the midplane is comparable, and in the disordered phase, there is a higher concentration of cholesterol in the midplane region. Poly-unsaturated lipid bilayers are disordered at biological temperatures, and these bilayers indeed exhibit cholesterol in the midplane region.14,53 Our work explains this experimental observation by demonstrating that cholesterol localization is dictated by the ordering of lipids in the bilayers. Figure 8 shows snapshots of the three kinds of lipid bilayers studied by us with 30% mole cholesterol and at 290 K. These snapshots illustrate the differences in bilayer thickness, and the orientation and location of cholesterol molecules in these three bilayers.

FIG. 7.

Ensemble-averaged angle of a cholesterol molecule with the bilayer normal as a function of the location of the hydroxyl headgroup of cholesterol in different molar concentrations of cholesterol in the lipid bilayers: (a) (blue) DOPC/CHOL, (b) (red) 20:0SM/CHOL, and (c) (black) DOPC/20:0SM/CHOL. The shaded area represents the midplane of the bilayers. Error bars are obtained from block averages of three equally sized blocks.

FIG. 7.

Ensemble-averaged angle of a cholesterol molecule with the bilayer normal as a function of the location of the hydroxyl headgroup of cholesterol in different molar concentrations of cholesterol in the lipid bilayers: (a) (blue) DOPC/CHOL, (b) (red) 20:0SM/CHOL, and (c) (black) DOPC/20:0SM/CHOL. The shaded area represents the midplane of the bilayers. Error bars are obtained from block averages of three equally sized blocks.

Close modal
FIG. 8.

Simulation snapshots of the three lipid bilayers studied by us with a molar cholesterol concentration of 30% and temperature of 290 K: (a) 20:0SM/CHOL, (b) DOPC/20:0SM/CHOL, and (c) DOPC/CHOL. In these snapshots, the blue beads represent the choline and phosphate groups of 20:0SM and DOPC. The green and gray beads represent the amide beads of 20:0SM and glycol beads of DOPC, respectively. The white and yellow colors represent cyclic rings and an alkyl tail of cholesterol, respectively. For the 20:0SM/CHOL system (a), the thickness of the bilayer is the largest of the other two systems, and most of the cholesterol molecules are parallel to the bilayer normal. For the DOPC/20:0SM/CHOL system (b), some cholesterol molecules are at the midplane. The cholesterol molecules at the midplane are perpendicular to the bilayer normal while the ones in the leaflets are parallel to the bilayer normal. For the DOPC/CHOL system (c), the bilayer thickness is the smallest, most cholesterol molecules are at the midplane, and they are randomly oriented with respect to the bilayer normal. In these snapshots, we have not shown the alkyl tails of the phospholipids for the sake of clarity.

FIG. 8.

Simulation snapshots of the three lipid bilayers studied by us with a molar cholesterol concentration of 30% and temperature of 290 K: (a) 20:0SM/CHOL, (b) DOPC/20:0SM/CHOL, and (c) DOPC/CHOL. In these snapshots, the blue beads represent the choline and phosphate groups of 20:0SM and DOPC. The green and gray beads represent the amide beads of 20:0SM and glycol beads of DOPC, respectively. The white and yellow colors represent cyclic rings and an alkyl tail of cholesterol, respectively. For the 20:0SM/CHOL system (a), the thickness of the bilayer is the largest of the other two systems, and most of the cholesterol molecules are parallel to the bilayer normal. For the DOPC/20:0SM/CHOL system (b), some cholesterol molecules are at the midplane. The cholesterol molecules at the midplane are perpendicular to the bilayer normal while the ones in the leaflets are parallel to the bilayer normal. For the DOPC/CHOL system (c), the bilayer thickness is the smallest, most cholesterol molecules are at the midplane, and they are randomly oriented with respect to the bilayer normal. In these snapshots, we have not shown the alkyl tails of the phospholipids for the sake of clarity.

Close modal

2. Diffusivity of phospholipids and cholesterol as a function of temperature and cholesterol concentration

We have plotted the in-plane diffusivity of phospholipids and cholesterol molecules as a function of temperature and cholesterol concentration [Fig. 9]. In the disordered phase, the diffusivity shows Arrhenius behavior (see Fig. S7 in the supplementary material54) showing the plots of ln (Diffusivity) versus 1000/T]. As expected, diffusivity decreases with cholesterol concentration, as an increase in the cholesterol concentration increases the lipid order.

FIG. 9.

Diffusivity of phospholipids (A1, B1, C1, C2) and cholesterol (A2, B2, C3) as a function of temperature and different cholesterol concentrations. In the disordered phase, the diffusivity follows Arrhenius behavior. The diffusivity at each temperature decreases with the increase in cholesterol concentration.

FIG. 9.

Diffusivity of phospholipids (A1, B1, C1, C2) and cholesterol (A2, B2, C3) as a function of temperature and different cholesterol concentrations. In the disordered phase, the diffusivity follows Arrhenius behavior. The diffusivity at each temperature decreases with the increase in cholesterol concentration.

Close modal

Along with the in-plane diffusivity, flip-flop rates of cholesterol across the leaflets are another measure of dynamics in the bilayers. In accordance with our previous work, we find that the flip-flop rates are slow in the ordered bilayers and are higher in the disordered bilayers. The flip-flop rate increases with temperature, as expected.24 

In this work, we show that the structural properties of lipid bilayers (area per lipid, bilayer thickness, and chain order parameter) obtained from the Dry Martini force field match well with those obtained from the Wet Martini force field. The Dry Martini force field is computationally more efficient, and thus, is useful for assessing longer time and length-scales than those afforded by the Wet Martini force field. By simulating three different lipid bilayers, namely, DOPC/CHOL, 20:0 SM/CHOL, and DOPC/20:0 SM/CHOL as a function of cholesterol concentration and temperature, we show that the spatial distribution of cholesterol molecules is closely related to the ordered and disordered phases of the bilayers. In the ordered phase, cholesterol molecules are preferentially distributed in the leaflets and align themselves parallel to the bilayer normal. In the disordered phase, cholesterol molecules are preferentially located in the midplane region of the bilayer. Interestingly, the melting temperature of the lipid bilayers, Tm is found to be the same as Tcross, defined as the temperature when the amount of cholesterol in the leaflets and the midplanes is equal. Thus, there is an intimate relationship between cholesterol distribution and lipid order. Not only does cholesterol impart a lipid order, but an increase in the lipid order also facilitates the localization of cholesterol in the leaflets rather than in the midplane region. These results can have important biological applications, as the distribution of cholesterol within the bilayer may impact the transport of material and the activity of membrane proteins. An interesting future line of inquiry would be if our findings will hold for other sterols, such as ergosterols.

This work was supported by the NSF CDS&E grant (No. 1953311). Computational resources for this work were provided by the National Science Foundation XSEDE grant (No. DMR190005) and the Ohio Supercomputer Center.

The authors have no conflicts to disclose.

No ethics approval is required for this work.

Mohammadreza Aghaaminiha: Data curation (lead); Formal analysis (lead); Investigation (lead); Methodology (equal); Software (lead); Validation (lead); Visualization (lead); Writing – original draft (lead). Amir M. Farnoud: Conceptualization (supporting); Methodology (supporting); Supervision (supporting); Writing – review & editing (equal). Sumit Sharma: Conceptualization (lead); Funding acquisition (lead); Project administration (lead); Resources (lead); Supervision (lead); Writing – review & editing (lead).

We have shared simulation files for the system of 16:0SM/CHOL with 30% cholesterol at the temperature of 290 K with the Wet Martini force field. The files are located here: https://github.com/ssumit12/lipid_bilayers. The files include the system topology files, the configuration files of lipids and water, the GROMACS simulation input files and the scripts for submitting simulation jobs to a supercomputer. These files are sufficient to generate the entire simulation trajectories. Other data are available on request from the authors.

1.
Y.
Zhang
,
A.
Lervik
,
J.
Seddon
, and
F.
Bresme
,
Chem. Phys. Lipids
185
,
88
(
2015
).
2.
N.
Kučerka
,
M.-P.
Nieh
, and
J.
Katsaras
,
Biochim. Biophys. Acta, Biomembr.
1808
,
2761
(
2011
).
3.
Y.
Wang
,
P.
Gkeka
,
J. E.
Fuchs
,
K. R.
Liedl
, and
Z.
Cournia
,
Biochim. Biophys. Acta, Biomembr.
1858
,
2846
(
2016
).
4.
A. G.
Lee
,
Biochim. Biophys. Acta, Biomembr.
1666
,
62
(
2004
).
5.
J.
Pan
,
T. T.
Mills
,
S.
Tristram-Nagle
, and
J. F.
Nagle
,
Phys. Rev. Lett.
100
,
198103
(
2008
).
6.
M. J.
Ueda
,
T.
Ito
,
T. S.
Okada
, and
S. I.
Ohnishi
,
J. Cell Biol.
71
,
670
(
1976
).
7.
L.
Bagatolli
and
P. B.
Sunil Kumar
,
Soft Matter
5
,
3234
(
2009
).
8.
G. W.
Feigenson
,
Biochim. Biophys. Acta, Biomembr.
1788
,
47
(
2009
).
9.
X.
Zhuang
,
J. R.
Makover
,
W.
Im
, and
J. B.
Klauda
,
Biochim. Biophys. Acta, Biomembr.
1838
,
2520
(
2014
).
10.
D.
Marquardt
,
N.
Kučerka
,
S. R.
Wassall
,
T. A.
Harroun
, and
J.
Katsaras
,
Chem. Phys. Lipids
199
,
17
(
2016
).
11.
P.
Khakbaz
and
J. B.
Klauda
,
Biochim. Biophys. Acta, Biomembr.
1860
,
1489
(
2018
).
12.
D.
Marquardt
,
F. A.
Heberle
,
J. D.
Nickels
,
G.
Pabst
, and
J.
Katsaras
,
Soft Matter
11
,
9055
(
2015
).
13.
J. T.
Buboltz
and
G. W.
Feigenson
,
Biochim. Biophys. Acta, Biomembr.
1417
,
232
(
1999
).
14.
T. A.
Harroun
,
J.
Katsaras
, and
S. R.
Wassall
,
Biochemistry
47
,
7090
(
2008
).
15.
A. B.
García-Arribas
,
A.
Alonso
, and
F. M.
Goñi
,
Chem. Phys. Lipids
199
,
26
(
2016
).
16.
N.
Kucerka
,
M.
Nieh
,
J.
Pencer
,
J.
Sachs
, and
J.
Katsaras
,
Gen. Physiol. Biophys.
28
,
117
(
2009
).
17.
S. O.
Yesylevskyy
and
A. P.
Demchenko
,
Eur. Biophys. J.
41
,
1043
(
2012
).
18.
M. R.
Ali
,
K. H.
Cheng
, and
J.
Huang
,
Proc. Natl. Acad. Sci. U.S.A.
104
,
5372
(
2007
).
19.
K.
Simons
and
D.
Toomre
,
Nat. Rev. Mol. Cell Biol.
1
,
31
(
2000
).
20.
W. F. D.
Bennett
,
J. L.
MacCallum
,
M. J.
Hinner
,
S. J.
Marrink
, and
D. P.
Tieleman
,
J. Am. Chem. Soc.
131
,
12714
(
2009
).
21.
F. J.
Alvarez
,
L. M.
Douglas
, and
J. B.
Konopka
,
Eukaryot. Cell
6
,
755
(
2007
).
22.
R. F. M.
de Almeida
,
A.
Fedorov
, and
M.
Prieto
,
Biophys. J.
85
,
2406
(
2003
).
23.
M.
Aghaaminiha
,
S. A.
Ghanadian
,
E.
Ahmadi
, and
A. M.
Farnoud
,
Biochim. Biophys. Acta, Biomembr.
1862
,
183350
(
2020
).
24.
M.
Aghaaminiha
,
A. M.
Farnoud
, and
S.
Sharma
,
Soft Matter
17
,
2742
(
2021
).
25.
S. R.
Wassall
and
W.
Stillwell
,
Chem. Phys. Lipids
153
,
57
(
2008
).
26.
C.
Arnarez
,
J. J.
Uusitalo
,
M. F.
Masman
,
H. I.
Ingólfsson
,
D. H.
de Jong
,
M. N.
Melo
,
X.
Periole
,
A. H.
de Vries
, and
S. J.
Marrink
,
J. Chem. Theory Comput.
11
,
260
(
2015
).
27.
J.
Michalowsky
,
L. V.
Schäfer
,
C.
Holm
, and
J.
Smiatek
,
J. Chem. Phys.
146
,
054501
(
2017
).
28.
K.
Pluhackova
,
S. A.
Kirsch
,
J.
Han
,
L.
Sun
,
Z.
Jiang
,
T.
Unruh
, and
R. A.
Böckmann
,
J. Phys. Chem. B
120
,
3888
(
2016
).
29.
H. I.
Ingólfsson
et al,
J. Am. Chem. Soc.
136
,
14554
(
2014
).
30.
S.
Nazemidashtarjandi
,
A.
Vahedi
, and
A. M.
Farnoud
,
Langmuir
36
,
4923
(
2020
).
31.
S.
Nazemidashtarjandi
,
V. M.
Sharma
,
V.
Puri
,
A. M.
Farnoud
, and
M. M.
Burdick
,
ACS Nano
16
,
2233
(
2022
).
32.
B.
Ramstedt
and
J. P.
Slotte
,
Biochim. Biophys. Acta, Biomembr.
1758
,
1945
(
2006
).
33.
M.
Caffrey
,
Lipid Thermotropic Phase Transition Database (LIPIDAT). User’s Guide Version 1.0
(
U.S. Department of Commerce, National Institute of Standards and Technology, Standard Reference Data Program
,
Gaithersburg, MD
,
1993
).
34.
J. R.
Silvius
,
Thermotropic Phase Transitions of Pure Lipids in Model Membranes and Their Modifications by Membrane Proteins
(
Wiley
,
New York
,
1982
).
35.
X.-M.
Li
,
J. M.
Smaby
,
M. M.
Momsen
,
H. L.
Brockman
, and
R. E.
Brown
,
Biophys. J.
78
,
1921
(
2000
).
36.
T. A.
Wassenaar
,
H. I.
Ingólfsson
,
R. A.
Böckmann
,
D. P.
Tieleman
, and
S. J.
Marrink
,
J. Chem. Theory Comput.
11
,
2144
(
2015
).
37.
M. J.
Abraham
,
T.
Murtola
,
R.
Schulz
,
S.
Páll
,
J. C.
Smith
,
B.
Hess
, and
E.
Lindahl
,
SoftwareX
1–2
,
19
(
2015
).
38.
C. M.
MacDermaid
,
H. K.
Kashyap
,
R. H.
DeVane
,
W.
Shinoda
,
J. B.
Klauda
,
M. L.
Klein
, and
G.
Fiorin
,
J. Chem. Phys.
143
,
243144
(
2015
).
39.
S. Y.
Bhide
,
Z.
Zhang
, and
M. L.
Berkowitz
,
Biophys. J.
92
,
1284
(
2007
).
40.
S. J.
Marrink
,
A. H.
de Vries
, and
A. E.
Mark
,
J. Phys. Chem. B
108
,
750
(
2004
).
41.
C. L.
Armstrong
,
D.
Marquardt
,
H.
Dies
,
N.
Kučerka
,
Z.
Yamani
,
T. A.
Harroun
,
J.
Katsaras
,
A.-C.
Shi
, and
M. C.
Rheinstädter
,
PLoS One
8
,
e66162
(
2013
).
42.
H. I.
Petrache
,
S.
Tristram-Nagle
,
K.
Gawrisch
,
D.
Harries
,
V. A.
Parsegian
, and
J. F.
Nagle
,
Biophys. J.
86
,
1574
(
2004
).
43.
P. R.
Maulik
and
G. G.
Shipley
,
Biochemistry
35
,
8025
(
1996
).
44.
A.
Liu
and
X.
Qi
,
Comput. Mol. Biosci.
02
,
78
(
2012
).
45.
S. Y.
Bhide
and
M. L.
Berkowitz
,
J. Chem. Phys.
123
,
224702
(
2005
).
46.
A. A.
Polyansky
,
P. E.
Volynsky
,
D. E.
Nolde
,
A. S.
Arseniev
, and
R. G.
Efremov
,
J. Phys. Chem. B
109
,
15052
(
2005
).
47.
P.
Niemelä
,
M. T.
Hyvönen
, and
I.
Vattulainen
,
Biophys. J.
87
,
2976
(
2004
).
48.
A.
Tsamaloukas
,
H.
Szadkowska
, and
H.
Heerklotz
,
Biophys. J.
90
,
4479
(
2006
).
49.
D.
Marquardt
et al,
Soft Matter
12
,
9417
(
2016
).
50.
Y.
Oh
and
B. J.
Sung
,
J. Phys. Chem. Lett.
9
,
6529
(
2018
).
51.
M. D.
Weiner
and
G. W.
Feigenson
,
J. Phys. Chem. B
122
,
8193
(
2018
).
52.
G.
Khelashvili
,
G.
Pabst
, and
D.
Harries
,
J. Phys. Chem. B
114
,
7524
(
2010
).
53.
S. J.
Marrink
,
A. H.
de Vries
,
Thad. A.
Harroun
,
J.
Katsaras
, and
S. R.
Wassall
,
J. Am. Chem. Soc.
130
,
10
(
2008
).
54.
See supplementary material at https://www.scitation.org/doi/suppl/10.1116/6.0002489 for simulation parameters for Dry and Wet Martini; a comparison of physical properties obtained for the Dry and Wet Martini simulations; a snapshot of equilibrium configuration of a studied lipid bilayer; estimation of melting temperature of lipid bilayers; and the relationship between log(diffusivity) of cholesterol and phospholipids and temperature.

Supplementary Material