Simulations of Asymmetric Membranes Illustrate Cooperative Leaflet Coupling and Lipid Adaptability

Biological membranes are composed of lipid bilayers that are often asymmetric with regards to the lipid composition and/or aqueous solvent they separate. Studying lipid asymmetry both experimentally and computationally is challenging. Molecular dynamics simulations of lipid bilayers with asymmetry are difficult due to finite system sizes and time scales accessible to simulations. Due to the very slow flip-flop rate for phospholipids, one must first choose how many lipids are on each side of the bilayer, but the resulting bilayer may be unstable (or metastable) due to differing tensile and compressive forces between leaflets. Here we use molecular dynamics simulations to investigate a number of different asymmetric membrane systems, both with atomistic and coarse-grained models. Asymmetries studied include differences in number of lipids, lipid composition (unsaturated and saturated tails and different headgroups), and chemical gradients between the aqueous phases. Extensive analysis of the bilayers’ properties such as area per lipid, density, and lateral pressure profiles are used to characterize bilayer asymmetry. We also address how cholesterol (which flip-flops relatively quickly) influences membrane asymmetries. Our results show how each leaflet is influenced by the other and can mitigate the structural changes to the bilayer overall structure. Cholesterol can respond to changes in bilayer asymmetry to alleviate some of the effect on the bilayer structure, but that will alter its leaflet distribution, which in turn affects its chemical potential. Ionic imbalances are shown to have a modest change in bilayer structure, despite large changes in the electrostatic potential. Bilayer asymmetry can also induce a modest electrostatic potential across the membrane. Our results highlight the importance of membrane asymmetry on bilayer properties, the influence of lipid headgroups, tails and cholesterol on asymmetry, and the ability of lipids to adapt to different environments.

Introduction

Lipid membranes separate organelles and cells from their exterior environments. Lipids can spontaneously assemble into a bilayer, which is the basic structure of cell membranes. These lipid bilayers function as a selective barrier to molecular entry/exit for the cell and serve as planar platforms for colocalizing various cell machinery. Because of their small scale, complexity and fluid phase behavior, it is difficult to study lipid bilayers directly, particularly in vivo. It has been long established that most cellular membranes have an asymmetric distribution of lipids between the inner and outer bilayer leaflets. Quantifying lipid compositions in either leaflet, transleaflet exchange (flip-flop), how cells create and maintain asymmetry, and the biological and physical effects of membrane asymmetry remain open research areas.

There are a number of ways that a membrane can exhibit asymmetry: different numbers of lipids, different types of lipids, asymmetric proteins, and by separating chemically different aqueous compartments. Removing lipids from one leaflet creates a stress in the membrane and can result in membrane curvature. Extreme examples are small unilamellar vesicles, which have high curvature and different numbers of lipids in either leaflet. In cells, regions of high curvature, such as membrane tubules, exosomes, and during vesicle fission/fusion, can results in lipid asymmetry (Callan-Jones et al., 2011; Skotland et al., 2017; Gruenberg, 2020). Cells tightly control the transbilayer distribution of lipids. Mammalian plasma membranes have phosphatidylcholine and sphingomyelin enriched in the outer leaflet and phosphatidylethanolamine and phosphatidylserine enriched in the inner leaflet (Van Meer et al., 2008; Lorent et al., 2020). The localization of cholesterol and rate of flip-flop remains debated (Lange and Steck, 2008; Mondal et al., 2009; Steck and Lange, 2018). Disruption of the membrane asymmetry with phosphatidylserine exposure on the outer leaflet is an important cell signal, for example during apoptosis (Van Meer et al., 2008). In vivo this is even more complex, because domains or lipid rafts can cause a heterogenous lipid distribution within the same leaflet. Considerable effort has been spent studying the interplay between lipid phase separation in one leaflet and asymmetric distribution of lipids across the membrane (Collins and Keller, 2008; Perlmutter and Sachs, 2011). A main question is if and how ordered domains enriched in cholesterol and sphingomyelin in the external leaflet can induce domains in the inner leaflet. Complexity is furthered by the active transport of lipids across membranes by transmembrane proteins. Additionally, cells tightly control the ion and osmotic gradients across the membrane with transmembrane protein channels and transporters. Key examples of how cells actively use the electrostatic potential difference are nerve cell propagation and energy production. The situation is even more complex, because the activity of the channels and transporters that create and maintain these asymmetries are in turn affected by their local lipid and ionic environment. Teasing out the interplay between the physicochemical effects of asymmetry and biological processes is ongoing.

Advances in computer simulations of lipid bilayers have made simulating complex and more biologically relevant membranes possible. This is exemplified by recent simulations on a realistic plasma membrane (Marrink et al., 2019), a neuronal membrane (Ingólfsson et al., 2017), bacterial membranes (Khalid et al., 2015), and entire organelles, such as a chromatophore vesicle (Singharoy et al., 2019). While at the outset it seems trivial for one to set up an asymmetric membrane in silico, as one can simply place lipids on one side or the other, there are problems that can create artifacts (Park et al., 2015; Huber et al., 2019). Finite system sizes means small bilayer patches are necessary, which can influence asymmetry effects, such as membrane bending. Simulation times that are very short compared to slow lipid movements, such as phospholipid flip-flop (McNamee and McConnell, 1973; Tieleman and Marrink, 2006; Bennett et al., 2009b; Marquardt et al., 2017), can influence simulation results. Including transmembrane proteins that have asymmetric shape further complicates this problem. Therefore, one can simulate nearly any asymmetry one chooses, but deciding what is relevant is a more difficult problem. Using the Martini coarse-grained (CG) model (Marrink et al., 2007), it was shown that even completely removing one leaflet (i.e., infinite asymmetry) was possible in relatively short and small simulations due to the strong stabilizing effect of periodic boundary conditions (Esteban-Martín et al., 2009), and slow time scale for phospholipid flip-flop. Another example showed that a CG simulation of a POPC membrane with 15 mol% difference in the number of lipids between leaflets was metastable, and applying a small pressure to break the metastability resulted in a membrane bud forming (Siggel et al., 2019). A transmembrane scramblase protein was then shown to alleviate the asymmetry.

To date most molecular dynamics (MD) simulations of lipid bilayers have been symmetric in the composition of their two leaflets. Many simulations of complex asymmetric bilayers have chosen the composition based on matching the average area per lipid (APL) of symmetric bilayers that represent each of the different desired leaflets (Huber et al., 2019). Recent work by Doktorova and Weinstein (2018) suggest that an asymmetric membrane should be constructed not just so the overall membrane tension is zero but so the tension in each leaflet is zero. Subsequently Hossein and Deserno (2020) showed that differential stress can still exist within membranes with non-zero tension, suggesting a detailed balancing of compositional and lateral stress is needed for balancing the physical effects of asymmetry. There has also been a number of simulations investigating the role of ion and lipid gradients on the transmembrane electrostatic potential (Gurtovenko and Vattulainen, 2007; Gurtovenko and Lyulina, 2014). Large ion imbalances lead to membrane electroporation, where a water filled pore forms across the bilayer, which has also been studied extensively with simulations (Gurtovenko and Vattulainen, 2005; Böckmann et al., 2008; Gurtovenko et al., 2010). MD simulations have also looked at cholesterol in lipid domains and have provided physical insight into domain registration for asymmetric membranes (Perlmutter and Sachs, 2011; Tian et al., 2016). Simulations have shown how liquid-ordered domains on one leaflet can promote more order in the opposite leaflet (Perlmutter and Sachs, 2011; Tian et al., 2016).

We have conducted a systematic investigation of different types of membranes simulated with multiscale simulation protocols. This work provides detailed physical insights in membrane asymmetry with biological implications as well as interesting technical details for the simulation community. MD simulations using atomistic and CG models of lipid bilayer asymmetry are presented and the bilayers’ properties are analyzed and compared. We address a number of asymmetric membrane scenarios: a gradual change in asymmetry, different numbers of lipids, different types of lipid, and ion imbalances, with and without cholesterol. We calculate lipid properties and bilayer properties, including lipid density profiles, electrostatic potentials, lipid order parameters, APL, cholesterol flip-flop and lateral pressure profiles. These results show how changes in lipid asymmetry can have modest to no effect on some membrane properties but rather large changes for other properties, emphasizing when more detailed attention to asymmetry is necessary. By comparing atomistic and CG models, we can assess model and chemical differences and inform future work. Simulations with cholesterol illustrate how lipids that can flip-flop readily can alter the asymmetry and influence the resulting membrane properties.

Materials and Methods

Large CG Lipid Membrane Simulations and Setup

We set up the larger CG lipid bilayers with varying degrees of asymmetry using the Martini tool insane (Wassenaar et al., 2015). These where simulated at the CG resolution only, contained up to 222 lipids in each leaflet and were run with an asymmetric DPSM (d18:1/18:0-sphingomyelin), POPC (16:0/18:1-PC), and PAPE (16:0/20:4-phosphatidylethanolamine) lipid mixture with a 40/40/20 ratio in the upper leaflet and 10/40/50 in the lower leaflet. Simulations were done with and without a high 50 mol% cholesterol content with lipid count asymmetry ranging from 0 to 35 lipids from the lower leaflet. All simulations were setup with 150 mM NaCl (∼75 beads each) and >15 CG water beads per lipid (∼6700 beads), corresponding to >60 waters per lipid.

Small CG and AA Lipid Membrane Simulation Setup

Systems with up to 64 lipids in each leaflet were setup with DPPC (16:0-PC) and DIPC (18:2-PC) in both leaflets, with 0, 5, or 10 lipids removed from one leaflet (e.g., DPPC_5). Asymmetric membranes with DPPC in one leaflet and DIPC in the other were simulated [DIPC(upper)/DPPC(lower)], and with 5 and 10 lipids removed from the DIPC leaflet (e.g., DPPC_0_DIPC_5). Each system was simulated at the atomistic and CG resolution and with and without ∼30 mol% cholesterol (relative to the upper layer) placed in each leaflet at the start of the simulation (replacing the PC lipid to maintain 64 lipids). Simulations were done with an asymmetric bilayer with DPPE (16:0-PE) in the upper leaflet and DPPC in the lower leaflet (headgroup asymmetry), both with 64 lipids. Additionally, at the CG resolution a larger range of asymmetries for the DIPC(upper)/DPPC(lower) were simulated (see Supplementary Figure S3). Each system contained 150 mM NaCl (∼22 beads each) and ∼2000 CG water beads. We set up double bilayer systems by replicating a single bilayer twice in the z-dimension. For asymmetric bilayers we also flipped the bilayer leaflets, so that the upper and lower are on opposite sides for the two bilayers. To create ion imbalances, we simply changed one sodium to a chloride in one water bath and chloride to sodium in the other, leading to a 4e imbalance.

Martini Simulations

We used the Martini CG force field v2.0 (Marrink et al., 2007) and a recent (unpublished) version of cholesterol that updated the shape of the original cholesterol model (Marrink et al., 2007) similar to the updated cholesterol model described in Melo et al. (2015) but without using virtual-sites. A 20 fs time step was used with the leap-frog algorithm and the new-rf Martini simulation parameters (De Jong et al., 2016). Each of the 64 lipids per leaflet simulations was simulated for 2 μs at 333 K and the 222 lipids per leaflet simulations for 10 μs at 320 K, both using the velocity rescale thermostat (Bussi et al., 2007) with a 12 ps time constant. The systems with DPPE were simulated at 340 K to be above the melting temperature for DPPE. Semi-isotropic pressure coupling was used with Berendsen barostat (Berendsen et al., 1984) with 1 bar pressure, a compressibility of 3 × 10^–4 bar^–1, and a relaxation time of 12 ps. Electrostatics were computed using the reaction-field method, with a dielectric constant of 15 within a 1.1 nm cutoff and infinite beyond the cutoff.

Atomistic Simulations

CHARMM36 lipids (Klauda et al., 2010) with CHARMM specific (Durell et al., 1994; Neria et al., 1996) TIP3P water (Jorgensen et al., 1983) were simulated with GROMACS v5.1.4 (Van Der Spoel et al., 2005; Abraham et al., 2015). The starting structure was taken from the 200 ns frame of the Martini CG model and back mapped to CHARMM36 using the Backward method (Wassenaar et al., 2014). Simulations were then continued using a 2 fs time step, with four replicates of 250 ns, with the final 200 ns from each used for analysis. Lennard-Jones interactions were cut-off after 1.2 nm, with a force-switch-function from 1.0 to 1.2 nm. Long-range electrostatic interactions were calculated with the particle mesh Ewald method (Darden et al., 1993; Essmann et al., 1995). Temperature was maintained at 333 K using the Nose-Hoover thermostat (Martyna et al., 1992) with a time constant of 1 ps. Semi-isotropic pressure coupling was used to maintain a constant pressure of 1 bar lateral and normal to the plane of the bilayer. The Parrinello-Rahman barostat (Parrinello and Rahman, 1981) was used with a compressibility of 4.5 × 10^–5 and a coupling constant of 5 ps. Water was constrained with the SETTLE method (Miyamoto and Kollman, 1992), and hydrogens on lipids with P-LINCS (Hess et al., 1997; Hess, 2008).

Simulation Analyses

After the simulations we used a variety of computational methods to analyze the lipid membrane systems. The different analysis methods use a combination of GROMACS tools, custom Python scripts using MDAnalysis (Michaud-Agrawal et al., 2011; Melcr et al., 2018; Antila et al., 2019). as well as the FATSLiM program (Buchoux, 2016) for APL calculations using Voronoi tessellation and membrane thickness calculations. The atomistic lipid tail order parameter was calculated using:

$$S_{\text{CH}}=<3co{s}^2\mathrm{\theta }-1>/2$$

where θ is the angle between the C-H bond and the normal to bilayer. This was calculated using an analysis script in MDAnalysis (Michaud-Agrawal et al., 2011; Melcr et al., 2018; Antila et al., 2019). Density and electrostatic profiles were calculated with an in-house script using MDAnalysis (Michaud-Agrawal et al., 2011). The electrostatic profile was determined using the method of Sachs et al. (2004), where the partial charge distribution across the membrane is integrated twice to give the electrostatic potential. For electrostatic potentials with asymmetry, we used a double bilayer set-up, allowing the two water chambers to have different ion concentrations.

The lateral pressure profile (LPP) was calculated using a modified version of GROMACS (Vanegas et al., 2014; Torres-Sánchez et al., 2015). The local stress is calculated by post-processing the trajectory file which returns a local stress file for the simulation averaged over the user-selected time range. We used the Goetz-Loetsky force decomposition and a grid spacing of 0.1 nm for atomistic and 0.15 nm for CG simulations. The local stress tensor (σ) is extracted from the binary file using the program tensortools (Vanegas et al., 2014; Torres-Sánchez et al., 2015), the LPP estimated as (σ_xx + σ_yy)/2 and the leaflet surface tension (ST) for the upper and lower leaflet as the integral from the bilayer center of mass to the box end on either side. Note, for asymmetric bilayers the center of mass is not the exact center of the integrated LPP (where either leaflet has equal amplitude but opposite sign) but for the mixture simulated in Figure 1 the difference was within error of the LPP calculations see upper/lower symmetry in ST (Figure 1F).

FIGURE 1

Figure 1. Tuning lipid number asymmetry for an asymmetric DPSM/POPC/PAPE 40:40:20 upper and 10:40:50 lower mixture, with and without 50% cholesterol. (A) Snapshot showing a sideview of a simulation with cholesterol and a difference in number asymmetry of 20 lipids. Tails are colored white, gray and black for 0, 1, and 4 double bonds, cholesterol is yellow and the headgroups are gray, blue and cyan for SM, PC and PE, respectively. (B) Cholesterol fraction in the upper and lower leaflet of each system. (C) Upper and lower leaflet APL for increased lipid number asymmetry. APL for symmetric upper (red) and symmetric lower (black) are shown in solid lines. Suggested number asymmetry value for this mixture from ratio of symmetric simulations and minimal leaflet ST are shown in dashed and dotted lines, respectively. (D) Bilayer thickness. (E) LPPs averaged over repeat simulations and outliers removed (Supplementary Figure S1). (F) Upper and lower leaflet ST. All analysis is averaged over the last 2 ms of the 10 ms simulations and there are 4–5 repeats for each lipid number asymmetry value.

Cholesterol exhibits the unique ability among the lipids used in this paper to flip from one leaflet to another over the course of a simulation. To investigate the effect of cholesterol mobility on asymmetry, we calculated the number of cholesterol molecules in each leaflet during the CG simulations. We considered cholesterols belonging to each leaflet if the ROH beads were >1.2 nm from the leaflet as defined by MDAnalysis leaflet finder using the non-flip-flopping lipids PO4 beads.

The combined analysis tools and a pipeline for setting up CG simulations is available, as is, on GitHub¹.

Results

CG Membrane Asymmetry

To explore the desired lipid number asymmetry for complex asymmetrical bilayers we systematically simulated the asymmetric DPSM/POPC/PAPE mixture with 40:40:20 lipid ratio in the upper and 10:40:50 ratio in the lower leaflet, both with and without added 50 mol% cholesterol (Figure 1). In this mixture the upper leaflet has significantly more lipids with saturated tails (DPSM) and less with high level of unsaturation (PAPE) compared to the lower leaflet. The upper leaflet should therefore want to have a smaller APL than the lower leaflet. In a periodic system with a fixed number of lipids this would translate into a lipid number asymmetry with fewer lipids on the lower leaflet than the upper one. To determine what the number asymmetry should be, we systematically reduced the number of lipids in the lower leaflet, simulating 4–5 repeats at each asymmetry value ranging from 0 to 35 random lipids removed from the lower leaflet (reaching a number asymmetry of ∼16%). For the cholesterol-containing systems, both leaflets start with an equal amount of cholesterol, which then flip-flops between the leaflets. Figure 1B shows the cholesterol fraction in the upper and lower leaflets average over the last 2 ms (∼5% of the cholesterol is transitioning or resides in the bilayer middle and is excluded from this analysis). There is a net drive to move cholesterol from the lower leaflet into the upper leaflet, until a leaflet number asymmetry of 31 (14%) where the drive reverses.

Figure 1C shows the APL for systems with (right) and without (left) cholesterol as the number asymmetry increases. The APL increases for the lower leaflet and decreases for the upper leaflet in both systems as lipids are removed from the lower leaflet. The increase in the lower leaflet APL is ∼2% higher in amplitude than the decrease in the upper leaflet APL. For the same systems the bilayer thickness changes very little (∼0.02 nm) with increased asymmetry (Figure 1D, note the y-axes narrow range). With enough simulations the trend in thickness change can be picked up showing a reduced bilayer thickness with fewer lipids in both system types, with an interesting initial increase in thickness in the systems without cholesterol.

Figure 1E shows the LPP along the z-axis with increased number asymmetry, where each profile is the average of the 4–5 repeated simulations with a few outliers removed (Supplementary Figures S1, S2). For the systems without cholesterol in the upper leaflet (on the left of the LPP plots), the LPP starts with a deeper minimum and a shallower maximum which then reverses as the number asymmetry increases. This is the opposite for the lower leaflet (right side of the LPP plots). For the cholesterol containing systems the trend is less clear due to the lower relative change in the LPP profile and more complex shape. With increased number asymmetry the minimum for the upper leaflet (left side of the LPP plots) reduces in amplitude as well as the central maximum and the opposite holds for the maximum of the lower leaflet (right side of the LPP plots). Figure 1F quantifies the ST for each leaflet in all the systems (the integral of the LPP to the bilayer center). Both with and without cholesterol the amplitude and the ST difference between the leaflets reduces as the number asymmetry increases, and then flips sign and increases again. The amplitude in the cholesterol-containing systems is notably smaller in spite of the high range of the LPPs. The zero ST is reached at an asymmetry of 15–16 without cholesterol and 14–15 with cholesterol (Figure 1C dashed lines). For comparison asymmetry estimates from the APL ratio of symmetric upper and symmetric lower simulations give 12.1 (without cholesterol) and 11.4 (with cholesterol) (Figure 1C dotted lines).

To address the local effects of leaflet asymmetry in more detail, we compared atomistic and CG simulations of small bilayer patches composed of pure DPPC and pure DIPC. Figure 2 shows the density plots for the different CG bilayers. For the pure bilayers, a ∼8% (5 lipids) and ∼15% (10 lipids) difference in composition had a negligible effect on the lipid bilayer density. Particularly for the DIPC bilayer, the profiles are nearly identical. We then simulated a mixed bilayer with one leaflet DPPC and the other DIPC. The DIPC(upper)/DPPC(lower) bilayers had substantial changes, with a density in between the pure DPPC and DIPC density profiles (Figure 2C). Both leaflets respond in opposite directions, indicating that the DPPC leaflets became more fluid and the DIPC leaflets more structured. Similar to the pure bilayers, changing the number of lipids in the mixed bilayers (removing 5 DIPC), also had little effect on the density.

FIGURE 2

Figure 2. Partial density profiles for the CG lipid bilayers without cholesterol. (A) Profiles for DPPC. (B) Partial densities for DIPC. (C) DIPC(upper)/ DPPC(lower) bilayers, with thin lines showing the individual DIPC or DPPC densities.

Table 1 lists the APL for each of the lipids in either the upper or lower leaflet. As expected, the leaflets with fewer lipids have a higher APL, and the ones with more lipids have a lower APL, compared to the pure bilayers. The mixed DIPC(upper)/DPPC(lower) bilayer has an APL in between the pure DPPC and DIPC membranes, as expected from the density profile. Supplementary Figure S3 plots the APL for the DIPC(upper)/DPPC(lower) bilayers as a function of the number of lipids removed from the DIPC(upper) leaflet. For reference, we also plot the APL for each type of lipid in pure symmetric bilayers (horizontal lines). The vertical line indicates the suggested asymmetry from symmetrical upper/lower leaflet simulations N(upper) = N(lower) × apl (lower)/apl(upper), which gives 51.8, or ∼12 fewer lipids for the DIPC(upper) leaflet.

TABLE 1

Table 1. Area per phospholipid for the CG lipid bilayers.

We tested the effect of adding 30 mol% cholesterol to each of the bilayers. Cholesterol is able to flip-flop on a time scale accessible to CG simulations. This means that cholesterol is able to equilibrate between leaflets and affect the asymmetric bilayers structure. Supplementary Figure S4 shows how the number of cholesterol molecules changes during the CG simulations. Supplementary Figure S4A shows relatively slow exchange of cholesterol in the saturated DPPC bilayer and very fast exchange in the DIPC bilayer. For the asymmetric DIPC(upper)/DPPC(lower) bilayer, the cholesterol molecules concentrate in the DPPC leaflet, depleting the DIPC leaflet. These results match expectations from numerous previous studies on cholesterol flip-flop and chemical potential (Bennett et al., 2009a; Jo et al., 2010; Javanainen and Martinez-Seara, 2019). We note that the final structure at 200 ns was used to start the atomistic simulations. For the atomistic saturated DPPC bilayers, the cholesterol does not undergo flip-flop over the 1 μs of simulation, and rarely for the DIPC bilayer (Supplementary Figure S5). For the DPPC_0 system, we know that the cholesterol distribution must be equal, so we picked a frame (175 ns) that has an equal distribution of cholesterol, as opposed to the frame at 200 ns, which has an asymmetric distribution of cholesterol. For the asymmetric membranes, it is not clear that we have the correct distribution of cholesterol, because its chemical potential might be different in the two leaflets. It is not clear that the chemical potential difference between the CG and AA model would be the same, so we chose the frame at the end of the simulation. Future work could incorporate free energy calculations and/or simulate an ensemble of CG frames for AA conversion.

Figure 3 shows the partial density profiles for the CG bilayers with 30 mol% cholesterol. These profiles show more substantial differences compared to the pure PC bilayers, but the same trends as bilayers without cholesterol. Notably, the DIPC and the DPPC densities have substantial differences compared to the symmetric bilayers, as well as cholesterol, but the overall bilayer density is similar for the asymmetric bilayers. Figure 3C illustrates how cholesterol is concentrated in the DPPC leaflet in the DIPC(upper)/DPPC(lower) bilayer. The APL for each CG bilayer are listed in Table 1 and show very modest changes with respect to the symmetric bilayers. Overall, the APL is much lower than the bilayers without cholesterol, due to cholesterols condensing effect.

FIGURE 3

Figure 3. Partial density profiles for the CG lipid bilayers with 30 mol% cholesterol. (A) Densities for DPPC bilayers. (B) Partial densities for DIPC. (C) The DIPC(upper)/DPPC(lower) bilayer. Thin lines are for the whole bilayer, medium lines are the PC lipids, and cholesterol is thick lines.

Lateral pressure profiles are an important property of lipid membranes, as they show the stress/pressure at each slice through the lipid membrane. Figure 4 shows the LPPs for the CG model with pure DPPC, pure DIPC, and both with 30 mol% cholesterol. For the pure bilayers, there is a large negative peak at the headgroup region and a shallow peak at the center of the membrane. Removing lipids from one leaflet of the DPPC bilayers has little effect on the height of the lateral pressure peak and trough (Figure 4A). There is a shift in the LPPS in the region corresponding to the middle of the lipid tails. A similar trend is observed for the polyunsaturated DIPC bilayers, although a slight decrease in pressure at the headgroup is also observed (Figure 4B). The magnitude of the pressures for the unsaturated DIPC bilayers are substantially lower than the saturated DPPC bilayers. Compared to the pure DPPC bilayers, LPPs for DPPC with 30 mol% cholesterol have a similar trend when removing lipids, but a much lower trough in pressure at the bilayer center (Figure 4C). The LPPs for the DIPC bilayers with and without cholesterol are very similar (Figure 4D compared to Figure 4B).

FIGURE 4

Figure 4. Lateral pressure profiles for CG bilayers with asymmetric number of lipids between leaflets. (A) DPPC. (B) DPPC:CHOL. (C) DIPC. (D) DIPC:CHOL. (E) DIPC(upper)/DPPC(lower). (F) DIPC(upper)/DPPC(lower):CHOL.

AA Membrane Asymmetry

We tested the effect of membrane asymmetry on the structure of atomistic lipid membranes by simulating bilayers with the CHARMM36 lipids (AA). For each system, after 200 ns of CG simulation, the membrane structure was used to back map to an atomistic bilayer, using the Backward method (Wassenaar et al., 2014), and the AA simulation started. Figure 5 shows the resulting partial density profiles for the phospholipids in the membranes. These profiles show that without cholesterol the bilayer’s partial density is modified slightly by relatively large changes in the number of lipids in each leaflet. As expected, the leaflet with fewer lipids resulted in a slightly lower density, and the one with more resulted in a larger density compared to the symmetric bilayers. The differences in density are more pronounced than the CG density differences (Figure 2).

FIGURE 5

Figure 5. Partial density profiles for the AA bilayers without cholesterol. (A) DPPC density profiles. (B) DIPC density profiles. (C) Asymmetric bilayers with DPPC in one leaflet (lower) and DIPC in the other (upper). Thin lines are for densities for DPPC or DIPC lipids separately.

The APL was also calculated for each AA lipid bilayer (Table 2). The results show a trend expected from the shift in density profiles. By removing lipids in one leaflet, the opposite leaflet compensates, so the APL increases in the leaflet with fewer lipids and decreases in the leaflet with more lipids, compared to the pure bilayer. The DIPC(upper)/DPPC(lower) membrane has an APL in between the pure DPPC and pure DIPC areas per lipid.

TABLE 2

Table 2. Area per phospholipid for the AA bilayers.

Partial density profiles for the AA bilayer with cholesterol are shown in Figure 6. Similar to the CG results, we observe substantial differences in the PC lipid density, by removing lipids in one leaflet, but the cholesterol density compensates, and the overall bilayer density has little change. The DIPC(upper)/DPPC(lower) bilayer with cholesterol, shows a substantial concentration of cholesterol in the lower DPPC bilayer, causing a large difference in the overall bilayers density between leaflets. This matches the results from the CG model. The DIPC(upper)/DPPC(lower) bilayer with five fewer lipids in the DIPC bilayer has a similar trend as the pure PC membranes, with cholesterol compensating for asymmetric PC density. This results in very similar overall bilayer densities for the DIPC(upper)/DPPC(lower) bilayer with and without equal number of lipids in either leaflet. The APL for the AA cholesterol membranes (Table 2) have the same trends as the CG membranes, with little changes in area when lipids are removed from only one leaflet.

FIGURE 6

Figure 6. Partial density profiles for the AA lipid bilayers with 30 mol% cholesterol. (A) Profiles for DPPC. (B) Partial densities for DIPC. (C) Bilayers with asymmetric DPPC in one leaflet and DIPC in the other. Thin lines are for the whole bilayer, medium lines are the PC lipids, and cholesterol is thick lines.

To assess the effect of asymmetry on the individual lipid’s structures, we calculated the lipid tail’s order parameter. The order parameter is a measure of the membranes ordering, based on the lipid tails alignment with the membrane normal. A value of 0 indicates no ordering and 1 completely ordered. Figure 7A shows that the bilayers compensate for the change in number of lipids by the leaflet with more lipids becoming more ordered, and the leaflet with fewer lipids becoming less ordered. Figure 7B shows that the DPPC lipids are more ordered with cholesterol, and the effect of asymmetry is reduced. For the 5 DPPC lipid difference, we observe the opposite trend as the other membranes, with the leaflet with more lipids having a slightly lower order parameter, and the one with fewer lipids more ordered. This could be due to an artifact of having the wrong distribution of cholesterol from the starting CG structure, although the difference is quite small. Future work explicitly studying cholesterol’s chemical potential would be of interest to properly assess this discrepancy.

FIGURE 7

Figure 7. DPPC lipid chain order parameters for the DPPC’s sn-2 chain. (A) Order parameters for the pure DPPC bilayer, comparing the top and bottom leaflets, and (B) with 30 mol% cholesterol. (C) Comparing order parameters for the top leaflet in the different membranes with compositional asymmetry, and (D) with 30 mol% cholesterol.

We calculated the order parameter for DPPC lipids in the mixed composition membranes (Figure 7C). In the DIPC(upper)/DPPC(lower) membrane the DPPC molecules have less order, as they are forced to compensate for the more disordered lipids in the opposite leaflet. We also tested the effect of asymmetric membranes with different headgroups, with DPPC in one leaflet and DPPE in the other. The DPPC in the DPPE(upper)/DPPC(lower) bilayer had high order, with a very similar profile as the DPPC_5 upper leaflet. This illustrates how the leaflets can compensate for the other, with the DPPC becoming more ordered to accommodate the DPPE in the opposite leaflet, with a much lower APL. With 30 mol% cholesterol (Figure 7D), there are interesting differences, with the DPPC lipids in the DIPC(upper)/DPPC(lower) mixed membranes having a higher order parameter compared to the pure DPPC bilayer. This is due to the increased concentration of cholesterol in the DPPC leaflet.

The LPPs for the AA bilayers are shown in Figure 8. These pressure profiles have large local pressures of hundreds of bars. For the pure DPPC bilayers, there is a very large trough in the pressure at the headgroup region, and a large positive peak at the bilayer center. Removing lipids from one leaflet has an effect on the depth of the trough at the DPPC headgroup region and no effect on the height of the peak at the bilayer center. Similar to the CG model, the biggest change in the pressure profile is observed in the lipid tail region of the DPPC bilayers. Similar behavior is observed for the pure DIPC bilayers (Figure 8B), but with less of an effect compared to DPPC. The bilayers with 30 mol% cholesterol have qualitatively different pressure profiles, with larger pressures, and a second headgroup trough and headgroup peak, compared to the bilayers without cholesterol. This behavior has been observed in other simulations of cholesterol (Shahane et al., 2019). The cholesterol bilayers show a substantial shift in the headgroup region for the asymmetric bilayers. There is also no effect on the magnitude of the pressure peak at the membrane center for asymmetric bilayers.

FIGURE 8

Figure 8. Lateral pressure profiles for AA bilayers with asymmetric number of lipids between leaflets. (A) DPPC. (B) DPPC:CHOL. (C) DIPC. (D) DIPC:CHOL. (E) DIPC(upper)/DPPC(lower). (F) DIPC(upper)/DPPC(lower):CHOL.

Ion and Electrostatic Imbalances

The membranes electrostatic potential can have important consequences on many biological functions of membranes and membrane proteins. Figure 9A shows the system set-up and the electrostatic potential (Figure 9B) across lipid membranes that do not have a charge imbalance. The bilayers have a relatively large positive potential at the bilayer center. The strength of the potential depends on the bilayer structure and composition, with a lower potential for the more disordered, unsaturated bilayers. For symmetric bilayers, the potential is equal to zero on both sides of the membrane. The poly unsaturated DIPC bilayer has a lower potential at the bilayer center compared to the DPPC bilayer. Asymmetric membranes cause the potential to shift between the two aqueous environments, even though both contain 150 mM NaCl, and there is no charge imbalance. Reducing the number of lipids on one side (DPPC_10) causes an increase in potential of ∼100 mV, while having DIPC or DPPE on the other leaflet causes a potential of ∼-70 mV. We also show the electrostatic potential across a Martini DIPC bilayer (Figure 9B). The Martini model does not reproduce the positive potential at the membrane center, which is a well-known short coming of the model.

FIGURE 9

Figure 9. Electrostatic potential for the AA lipid membranes. (A) Snapshot of the double bilayer set-up. The DPPC lipids are shown as gray lines, ions are red and blue balls, and water is a translucent surface. (B) Bilayers with no charge imbalance between aqueous compartments. The potentials were set to zero on the left side of the bilayer. (C) Double bilayer systems with a 4e charge imbalance. Water regions have a flat potential due to conducting counter ions. (D) Lateral pressure profile for the DIPC double bilayer with and without a 4e charge imbalance.

Another source of membrane asymmetry is imbalances in the composition of the aqueous phases that the bilayer separates. Many biological membranes have an electrostatic potential across the bilayer. To investigate this effect, we simulated lipid bilayers with a net charge imbalance (4e) across the membrane. Figure 9C shows the electrostatic potential across pure DPPC and pure DIPC lipid bilayers with a 4e charge imbalance. The charge imbalance causes a potential difference of ∼1.2 V across the DPPC bilayer and ∼0.75 V across the DIPC bilayer. This illustrates the importance of the bilayer structure on the electrostatic potential gradients. The relatively large charge imbalance has a small effect on the bilayer structure, as shown by the lateral pressure profile in Figure 9D. The APL for DPPC increases to 0.63 from 0.62 nm² without a charge imbalance. Supplementary Figure S6 shows the order parameters for the 4e DPPC bilayer are slightly reduced compared to the pure DPPC bilayer. This modest change in structure is somewhat surprising, given that if we simulated long enough, a pore would likely form across the membrane, to allow the ions to equilibrate between water compartments.

Discussion

Biological membranes contain incredible complexity and diversity. While it has long been known that biological membranes contain asymmetry between one side and the other, much remains unknown about the origin, physical effect, biological function, and consequences of its disruption. Here we have investigated a number of different types of membrane asymmetry using both atomistic and CG models. Our results illustrate how lipid membranes are complex and cooperative aggregates. Changes on one side of the membrane affect the opposite leaflet. In general, the two leaflets attempt to balance asymmetric deviations. For example, if one leaflet becomes more ordered, the opposite leaflet will become more ordered. How the bilayer responds to the asymmetry depends on the composition and chemical structure of its lipids. We chose simulation boxes that were small, to preclude the effect of bilayer bending in response to asymmetry. Future work on the system size dependence of asymmetry would be of interest, or large-scale CG simulations of small unilamellar vesicles with asymmetric lipid distributions.

Here we have used both atomistic and CG simulations to study a number of asymmetric bilayer systems. While there are many detailed insights, the overall trend we observe is the coupling between the bilayers two leaflets. By running extensive simulations for CG systems with 3–4 lipid components and systematically changing the asymmetry, we show that the number lipids in each leaflet can be adjusted to minimize leaflet ST (Figures 1C,F). The resulting number asymmetry is somewhat higher than using estimates from the APL of the two symmetric bilayers. This effect had been shown previously (Doktorova and Weinstein, 2018; Hossein and Deserno, 2020), but here we have run extensive CG simulations to characterize the behavior with and without cholesterol. Hossein and Deserno (2020) also showed that in some cases a more complex balance between spontaneous curvature and tension is desired. Because the two leaflets are interacting, they will alter the others elastic properties and therefore shift the packing of the lipids. While this can have important consequences, we also note that many bilayer properties are quite similar over a wide range of asymmetries (Figure 1D), suggesting that getting the exact number asymmetry may not seriously impact many simulation results. This effect is also due to the lipid’s flexibility, where each leaflet compensates, so the difference between the two leaflets is minimized.

We have also addressed the influence of cholesterol on membrane asymmetry. Cholesterol has been hypothesized to be able to reduce membrane asymmetry by being able to flip quickly from one leaflet to the other and equilibrate membrane stress (Miettinen and Lipowsky, 2019). Our results suggest that for the systems tested here this is true for some properties, such as the overall membrane density and APL, cholesterol can also increase the effect on some other properties, such as the LPP. Cholesterol makes the situation of simulating membrane asymmetry more complicated, because it is able to flip-flop on a time scale easily accessible to CG simulations, but not for AA simulations of long saturated lipids. This is also likely the case for other fast flipping lipids such as fatty acids, ceramides and diacylglycerol (Bennett and Tieleman, 2012; Janke et al., 2014). Cholesterol’s distribution is also more complicated because it will be driven by its chemical potential for either leaflet, and its ability to influence the membrane’s spontaneous curvature (Allender et al., 2019; Hossein and Deserno, 2020). Considerable work in this area is needed to properly assess lipid distributions, such as free energy methods for determining their relative chemical potential in each leaflet.

By systematically investigating AA and CG models, we are able to compare specific difference for the effect of membrane asymmetry. We find that many of the properties and differences are captured using both levels of detail. While there are substantial differences in the shape and magnitude of the LPPs, many of the changes for asymmetric membranes are similar for the two models. This is crucial for future simulation work using multiscale approaches to study asymmetric membranes. These methods are attractive due to the fast sampling of the CG models, which is necessary for many slow processes involving lipids, and the fine chemical details provided by the AA models. Ensuring that the CG asymmetric lipid distribution matches the AA models is necessary, or the AA simulation will be biased by the initial CG simulation result. The electrostatic potentials for Martini are shown to be opposite to that of AA models, illustrating a fundamental limitation for the model. It is also likely that the large preference of cholesterol for DPPC over DIPC in the asymmetric bilayer could be due to Martini’s strong repulsion between cholesterol and poly unsaturated phospholipids. Future work comparing the lipid chemical potentials for both models might help provide more qualitative differences. The overall good agreement between the AA and CG results is promising for large scale multiscale membrane simulations. Future work investigating the effect of proteins inclusion in asymmetric membrane systems would also be of interest. Another aspect of asymmetry for future work are pH gradients, which require more sophisticated computational methods.

We examined how asymmetric ion distributions can affect the membrane structure, as well as how asymmetric lipid distributions influence the electrostatic potential across the membrane. Our results match previous results showing that changes in lipid headgroup can shift the overall membrane electrostatic potential (Gurtovenko and Vattulainen, 2007). We also show that changing the lipid tails can have a substantial effect on the membrane’s dipole potential. The fact that lipid distribution influences membrane electrostatic potentials has important implications on numerous biological processes, such as nerve cell activation. Another source of membrane asymmetry is the bilayer separating aqueous compartments with chemical and/or ionic gradients. We show that relatively large electrostatic potentials have a modest effect on the membranes overall structure. This is somewhat surprising, as these membranes are likely to form pores if simulated for long times, so the ions can equilibrate (Gurtovenko and Vattulainen, 2005). These large electrostatic gradients might promote structural fluctuations that promote pore formation, such as water wires across the membrane, which have been shown to be a pre-pore state for pore formation in lipid bilayers (Bennett et al., 2014). Our results also illustrate how changes in lipid tails have a large effect on the electrostatic potential resulting from ionic imbalances. These differences are important for many biological processes, including the activity of voltage-gated ion channels (Bezanilla, 2007), GPCR signaling (Rinne et al., 2015) and the insertion and penetration of charged peptides, such as antimicrobial peptides (Teixeira et al., 2012) and cell-penetrating peptides (Moghal et al., 2020).

Conclusion

Asymmetry in lipid bilayer simulations is a crucial form of complexity in biological systems. Cells create and exploit asymmetry for energy production, signaling, and general transport of molecules. Studying membrane asymmetry in model systems is challenging for experiments and simulations, as this is an inherently non-equilibrium situation. For simulations, finite system sizes and time scales makes membrane asymmetry hard to study. Our results show how changes in number of lipids, lipid headgroups, lipid tails, ionic imbalances, and the effect of cholesterol impact bilayer’s structural properties. Overall, we show that CG Martini simulations reproduce the effects observed with AA CHARMM36 simulations. One leaflets structure effects the opposite leaflet, as the individual lipids adapt to their local environment. Fast flip-flop of cholesterol creates additional challenges, for simulating asymmetric membrane simulations. Cholesterol reduces the impact of membrane asymmetry for some structural properties, but can create more imbalance, especially when it has a large preference for lipids in one leaflet over the other. There is much work in the future for membrane simulations of asymmetry, including assessing chemical potential differences and enhanced sampling methods.

Data Availability Statement

The datasets generated for this study are available on request to the corresponding author.

Author Contributions

MB, SH, ML, LM, and MU wrote the analysis scripts. PS and HI managed the project. All authors contributed to running simulations, analyzing the result, and conceiving the project.

Funding

This work was partially supported by the Laboratory Directed Research and Development at the Lawrence Livermore National Laboratory (18-ERD-035) and performed under the auspices of the U.S. DOE by the Lawrence Livermore National Laboratory under contract DE-AC52-07NA27344. Release LLNL-JRNL-801959.

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Acknowledgments

We thank Fred Streitz for his guidance and support, and Livermore Computing for access to computing time. MB, SH, ML, LM, and MU thank the Harvey Mudd Clinic Program for support.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fcell.2020.00575/full#supplementary-material

Footnotes

1. ^ https://github.com/llnl-hmc-clinic/Bilayer-Simulation-Pipeline

References

Abraham, M. J., Murtola, T., Schulz, R., Páll, S., Smith, J. C., Hess, B., et al. (2015). GROMACS: high performance molecular simulations through multi-level parallelism from laptops to supercomputers. SoftwareX 1–2, 19–25. doi: 10.1016/j.softx.2015.06.001

CrossRef Full Text | Google Scholar

Allender, D., Sodt, A., and Schick, M. (2019). Cholesterol-dependent bending energy is important in cholesterol distribution of the plasma membrane. Biophys. J. 116, 2356–2366. doi: 10.1016/j.bpj.2019.03.028

PubMed Abstract | CrossRef Full Text | Google Scholar

Antila, H., Buslaev, P., Favela-Rosales, F., Ferreira, T. M., Gushchin, I., Javanainen, M., et al. (2019). Headgroup structure and cation binding in phosphatidylserine lipid bilayers. J. Phys. Chem. B 123, 9066–9079. doi: 10.1021/acs.jpcb.9b06091

PubMed Abstract | CrossRef Full Text | Google Scholar

Bennett, W. F. D., MacCallum, J. L., Hinner, M. J., Marrink, S. J., and Tieleman, D. P. (2009a). Molecular view of cholesterol flip-flop and chemical potential in different membrane environments. J. Am. Chem. Soc. 131, 12714–12720. doi: 10.1021/ja903529f

PubMed Abstract | CrossRef Full Text | Google Scholar

Bennett, W. F. D., MacCallum, J. L., and Tieleman, D. P. (2009b). Thermodynamic analysis of the effect of cholesterol on dipalmitoylphosphatidylcholine lipid membranes. J. Am. Chem. Soc. 131, 1972–1978. doi: 10.1021/ja808541r

PubMed Abstract | CrossRef Full Text | Google Scholar

Bennett, W. F. D., Sapay, N., and Tieleman, D. P. (2014). Atomistic simulations of pore formation and closure in lipid bilayers. Biophys. J. 106, 210–219. doi: 10.1016/j.bpj.2013.11.4486

PubMed Abstract | CrossRef Full Text | Google Scholar

Bennett, W. F. D., and Tieleman, D. P. (2012). Molecular simulation of rapid translocation of cholesterol, diacylglycerol, and ceramide in model raft and nonraft membranes. J. Lipid Res. 53, 421–429. doi: 10.1194/jlr.M022491

PubMed Abstract | CrossRef Full Text | Google Scholar

Berendsen, H. J., van Postma, J., van Gunsteren, W. F., DiNola, A., and Haak, J. R. (1984). Molecular dynamics with coupling to an external bath. J. Chem. Phys. 81, 3684–3690. doi: 10.1063/1.448118

CrossRef Full Text | Google Scholar

Bezanilla, F. (2007). “Voltage-gated ion channels,” in Biological Membrane Ion Channels, eds S. H. Chung, O. S. Andersen, and V. Krishnamurthy (Cham: Springer), 81–118. doi: 10.1007/0-387-68919-2_3

CrossRef Full Text | Google Scholar

Böckmann, R. A., De Groot, B. L., Kakorin, S., Neumann, E., and Grubmüller, H. (2008). Kinetics, statistics, and energetics of lipid membrane electroporation studied by molecular dynamics simulations. Biophys. J. 95, 1837–1850. doi: 10.1529/biophysj.108.129437

PubMed Abstract | CrossRef Full Text | Google Scholar

Buchoux, S. (2016). FATSLiM: a fast and robust software to analyze MD simulations of membranes. Bioinformatics 33, 133–134. doi: 10.1093/bioinformatics/btw563

PubMed Abstract | CrossRef Full Text | Google Scholar

Bussi, G., Donadio, D., and Parrinello, M. (2007). Canonical sampling through velocity rescaling. J. Chem. Phys. 126:014101. doi: 10.1063/1.2408420

CrossRef Full Text | Google Scholar

Callan-Jones, A., Sorre, B., and Bassereau, P. (2011). Curvature-driven lipid sorting in biomembranes. Cold Spring Harb. Perspect. Biol. 3:a004648. doi: 10.1101/cshperspect.a004648

PubMed Abstract | CrossRef Full Text | Google Scholar

Collins, M. D., and Keller, S. L. (2008). Tuning lipid mixtures to induce or suppress domain formation across leaflets of unsupported asymmetric bilayers. Proc. Natl. Acad. Sci. U.S.A. 105, 124–128. doi: 10.1073/pnas.0702970105

PubMed Abstract | CrossRef Full Text | Google Scholar

Darden, T., York, D., and Pedersen, L. (1993). Particle mesh Ewald: an N⋅ log (N) method for Ewald sums in large systems. J. Chem. Phys. 98, 10089–10092. doi: 10.1063/1.464397

CrossRef Full Text | Google Scholar

De Jong, D. H., Baoukina, S., Ingólfsson, H. I., and Marrink, S. J. (2016). Martini straight: boosting performance using a shorter cutoff and GPUs. Comput. Phys. Commun. 199, 1–7. doi: 10.1016/j.cpc.2015.09.014

CrossRef Full Text | Google Scholar

Doktorova, M., and Weinstein, H. (2018). Accurate in silico modeling of asymmetric bilayers based on biophysical principles. Biophys. J. 115, 1638–1643. doi: 10.1016/j.bpj.2018.09.008

PubMed Abstract | CrossRef Full Text | Google Scholar

Durell, S. R., Brooks, B. R., and Ben-Naim, A. (1994). Solvent-induced forces between two hydrophilic groups. J. Phys. Chem. 98, 2198–2202. doi: 10.1021/j100059a038

CrossRef Full Text | Google Scholar

Essmann, U., Perera, L., Berkowitz, M. L., Darden, T., Lee, H., and Pedersen, L. G. (1995). A smooth particle mesh Ewald method. J. Chem. Phys. 103, 8577–8593. doi: 10.1063/1.470117

CrossRef Full Text | Google Scholar

Esteban-Martín, S., Risselada, H. J., Salgado, J., and Marrink, S. J. (2009). Stability of asymmetric lipid bilayers assessed by molecular dynamics simulations. J. Am. Chem. Soc. 131, 15194–15202. doi: 10.1021/ja904450t

PubMed Abstract | CrossRef Full Text | Google Scholar

Gruenberg, J. (2020). Life in the lumen: the multivesicular endosome. Traffic 21, 76–93. doi: 10.1111/tra.12715

PubMed Abstract | CrossRef Full Text | Google Scholar

Gurtovenko, A. A., Anwar, J., and Vattulainen, I. (2010). Defect-mediated trafficking across cell membranes: insights from in silico modeling. Chem. Rev. 110, 6077–6103. doi: 10.1021/cr1000783

PubMed Abstract | CrossRef Full Text | Google Scholar

Gurtovenko, A. A., and Lyulina, A. S. (2014). Electroporation of asymmetric phospholipid membranes. J. Phys. Chem. B 118, 9909–9918. doi: 10.1021/jp5028355

PubMed Abstract | CrossRef Full Text | Google Scholar

Gurtovenko, A. A., and Vattulainen, I. (2005). Pore formation coupled to ion transport through lipid membranes as induced by transmembrane ionic charge imbalance: atomistic molecular dynamics study. J. Am. Chem. Soc. 127, 17570–17571. doi: 10.1021/ja053129n

PubMed Abstract | CrossRef Full Text | Google Scholar

Gurtovenko, A. A., and Vattulainen, I. (2007). Lipid transmembrane asymmetry and intrinsic membrane potential: two sides of the same coin. J. Am. Chem. Soc. 129, 5358–5359. doi: 10.1021/ja070949m

PubMed Abstract | CrossRef Full Text | Google Scholar

Hess, B. (2008). P-LINCS: a parallel linear constraint solver for molecular simulation. J. Chem. Theory Comput. 4, 116–122. doi: 10.1021/ct700200b

PubMed Abstract | CrossRef Full Text | Google Scholar

Hess, B., Bekker, H., Berendsen, H. J., and Fraaije, J. G. (1997). LINCS: a linear constraint solver for molecular simulations. J. Comput. Chem. 18, 1463–1472. doi: 10.1002/(sici)1096-987x(199709)18:12<1463::aid-jcc4>3.0.co;2-h

CrossRef Full Text | Google Scholar

Hossein, A., and Deserno, M. (2020). Spontaneous curvature, differential stress, and bending modulus of asymmetric lipid membranes. Biophys. J. 118, 624–642.

Google Scholar

Huber, R. G., Carpenter, T. S., Dube, N., Holdbrook, D. A., Ingólfsson, H. I., Irvine, W. A., et al. (2019). “Multiscale modeling and simulation approaches to lipid–protein interactions,” in Lipid-Protein Interactions, ed. J. Kleinschmidt (Cham: Springer), 1–30. doi: 10.1007/978-1-4939-9512-7_1

CrossRef Full Text | Google Scholar

Ingólfsson, H. I., Carpenter, T. S., Bhatia, H., Bremer, P.-T., Marrink, S. J., and Lightstone, F. C. (2017). Computational lipidomics of the neuronal plasma membrane. Biophys. J. 113, 2271–2280.

Google Scholar

Ingólfsson, H. I., Melo, M. N., van Eerden, F. J., Arnarez, C., Lopez, C. A., Wassenaar, T. A., et al. (2014). Lipid organization of the plasma membrane. J. Am. Chem. Soc. 136, 14554–14559. doi: 10.1021/ja507832e

PubMed Abstract | CrossRef Full Text | Google Scholar

Janke, J. J., Bennett, W. D., and Tieleman, D. P. (2014). Oleic acid phase behavior from molecular dynamics simulations. Langmuir 30, 10661–10667. doi: 10.1021/la501962n

PubMed Abstract | CrossRef Full Text | Google Scholar

Javanainen, M., and Martinez-Seara, H. (2019). Rapid diffusion of cholesterol along polyunsaturated membranes via deep dives. Phys. Chem. Chem. Phys. 21, 11660–11669. doi: 10.1039/c9cp02022e

PubMed Abstract | CrossRef Full Text | Google Scholar

Jo, S., Rui, H., Lim, J. B., Klauda, J. B., and Im, W. (2010). Cholesterol flip-flop: insights from free energy simulation studies. J. Phys. Chem. B 114, 13342–13348. doi: 10.1021/jp108166k

PubMed Abstract | CrossRef Full Text | Google Scholar

Jorgensen, W. L., Chandrasekhar, J., Madura, J. D., Impey, R. W., and Klein, M. L. (1983). Comparison of simple potential functions for simulating liquid water. J. Chem. Phys. 79, 926–935. doi: 10.1063/1.445869

CrossRef Full Text | Google Scholar

Khalid, S., Berglund, N. A., Holdbrook, D. A., Leung, Y. M., and Parkin, J. (2015). The membranes of Gram-negative bacteria: progress in molecular modelling and simulation. Biochem. Soc. Trans. 43, 162–167. doi: 10.1042/bst20140262

PubMed Abstract | CrossRef Full Text | Google Scholar

Klauda, J. B., Venable, R. M., Freites, J. A., O’Connor, J. W., Tobias, D. J., Mondragon-Ramirez, C., et al. (2010). Update of the CHARMM all-atom additive force field for lipids: validation on six lipid types. J. Phys. Chem. B 114, 7830–7843. doi: 10.1021/jp101759q

PubMed Abstract | CrossRef Full Text | Google Scholar

Lange, Y., and Steck, T. L. (2008). Cholesterol homeostasis and the escape tendency (activity) of plasma membrane cholesterol. Prog. Lipid Res. 47, 319–332. doi: 10.1016/j.plipres.2008.03.001

PubMed Abstract | CrossRef Full Text | Google Scholar

Lorent, J. H., Levental, K. R., Ganesan, L., Rivera-Longsworth, G., Sezgin, E., Doktorova, M., et al. (2020). Plasma membranes are asymmetric in lipid unsaturation, packing and protein shape. Nat. Chem. Biol. 16:710. doi: 10.1038/s41589-020-0564-3

PubMed Abstract | CrossRef Full Text | Google Scholar

Marquardt, D., Heberle, F. A., Miti, T., Eicher, B., London, E., Katsaras, J., et al. (2017). 1H NMR shows slow phospholipid flip-flop in gel and fluid bilayers. Langmuir 33, 3731–3741. doi: 10.1021/acs.langmuir.6b04485

PubMed Abstract | CrossRef Full Text | Google Scholar

Marrink, S. J., Corradi, V., Souza, P. C., Ingólfsson, H. I., Tieleman, D. P., and Sansom, M. S. (2019). Computational modeling of realistic cell membranes. Chem. Rev. 119, 6184–6226. doi: 10.1021/acs.chemrev.8b00460

PubMed Abstract | CrossRef Full Text | Google Scholar

Marrink, S. J., Risselada, H. J., Yefimov, S., Tieleman, D. P., and De Vries, A. H. (2007). The MARTINI force field: coarse grained model for biomolecular simulations. J. Phys. Chem. B 111, 7812–7824.

Google Scholar

Martyna, G. J., Klein, M. L., and Tuckerman, M. (1992). Nosé–Hoover chains: the canonical ensemble via continuous dynamics. J. Chem. Phys. 97, 2635–2643. doi: 10.1063/1.463940

CrossRef Full Text | Google Scholar

McNamee, M. G., and McConnell, H. M. (1973). Transmembrane potentials and phospholipid flip-flop in excitable membrane vesicles. Biochemistry 12, 2951–2958. doi: 10.1021/bi00740a001

PubMed Abstract | CrossRef Full Text | Google Scholar

Melcr, J., Martinez-Seara, H., Nencini, R., Kolafa, J., Jungwirth, P., and Ollila, O. S. (2018). Accurate binding of sodium and calcium to a POPC bilayer by effective inclusion of electronic polarization. J. Phys. Chem. B 122, 4546–4557. doi: 10.1021/acs.jpcb.7b12510

PubMed Abstract | CrossRef Full Text | Google Scholar

Melo, M., Ingólfsson, H., and Marrink, S. (2015). Parameters for Martini sterols and hopanoids based on a virtual-site description. J. Chem. Phys. 143:243152. doi: 10.1063/1.4937783

CrossRef Full Text | Google Scholar

Michaud-Agrawal, N., Denning, E. J., Woolf, T. B., and Beckstein, O. (2011). MDAnalysis: a toolkit for the analysis of molecular dynamics simulations. J. Comput. Chem. 32, 2319–2327. doi: 10.1002/jcc.21787

PubMed Abstract | CrossRef Full Text | Google Scholar

Miettinen, M. S., and Lipowsky, R. (2019). Bilayer membranes with frequent flip-flops have tensionless leaflets. Nano Lett. 19, 5011–5016. doi: 10.1021/acs.nanolett.9b01239

PubMed Abstract | CrossRef Full Text | Google Scholar

Miyamoto, S., and Kollman, P. A. (1992). Settle: an analytical version of the SHAKE and RATTLE algorithm for rigid water models. J. Comput. Chem. 13, 952–962. doi: 10.1002/jcc.540130805

CrossRef Full Text | Google Scholar

Moghal, M. M. R., Islam, M. Z., Hossain, F., Saha, S. K., and Yamazaki, M. (2020). Role of membrane potential on entry of cell-penetrating peptide transportan 10 into single vesicles. Biophys. J. 118, 57–69. doi: 10.1016/j.bpj.2019.11.012

PubMed Abstract | CrossRef Full Text | Google Scholar

Mondal, M., Mesmin, B., Mukherjee, S., and Maxfield, F. R. (2009). Sterols are mainly in the cytoplasmic leaflet of the plasma membrane and the endocytic recycling compartment in CHO cells. Mol. Biol. Cell 20, 581–588. doi: 10.1091/mbc.e08-07-0785

PubMed Abstract | CrossRef Full Text | Google Scholar

Neria, E., Fischer, S., and Karplus, M. (1996). Simulation of activation free energies in molecular systems. J. Chem. Phys. 105, 1902–1921. doi: 10.1063/1.472061

CrossRef Full Text | Google Scholar

Park, S., Beaven, A. H., Klauda, J. B., and Im, W. (2015). How tolerant are membrane simulations with mismatch in area per lipid between leaflets? J. Chem. Theory Comput. 11, 3466–3477. doi: 10.1021/acs.jctc.5b00232

PubMed Abstract | CrossRef Full Text | Google Scholar

Parrinello, M., and Rahman, A. (1981). Polymorphic transitions in single crystals: a new molecular dynamics method. J. Appl. Phys. 52, 7182–7190. doi: 10.1063/1.328693

CrossRef Full Text | Google Scholar

Perlmutter, J. D., and Sachs, J. N. (2011). Interleaflet interaction and asymmetry in phase separated lipid bilayers: molecular dynamics simulations. J. Am. Chem. Soc. 133, 6563–6577. doi: 10.1021/ja106626r

PubMed Abstract | CrossRef Full Text | Google Scholar

Rinne, A., Mobarec, J. C., Mahaut-Smith, M., Kolb, P., and Bünemann, M. (2015). The mode of agonist binding to a G protein–coupled receptor switches the effect that voltage changes have on signaling. Sci. Signal. 8:ra110. doi: 10.1126/scisignal.aac7419

PubMed Abstract | CrossRef Full Text | Google Scholar

Sachs, J. N., Crozier, P. S., and Woolf, T. B. (2004). Atomistic simulations of biologically realistic transmembrane potential gradients. J. Chem. Phys. 121, 10847–10851.

Google Scholar

Shahane, G., Ding, W., Palaiokostas, M., and Orsi, M. (2019). Physical properties of model biological lipid bilayers: insights from all-atom molecular dynamics simulations. J. Mol. Model. 25:76. doi: 10.1007/s00894-019-3964-0

PubMed Abstract | CrossRef Full Text | Google Scholar

Siggel, M., Bhaskara, R. M., and Hummer, G. (2019). Phospholipid scramblases remodel the shape of asymmetric membranes. J. Phys. Chem. Lett. 10, 6351–6354. doi: 10.1021/acs.jpclett.9b02531

PubMed Abstract | CrossRef Full Text | Google Scholar

Singharoy, A., Maffeo, C., Delgado-Magnero, K. H., Swainsbury, D. J., Sener, M., Kleinekathöfer, U., et al. (2019). Atoms to phenotypes: molecular design principles of cellular energy metabolism. Cell 179, 1098–1111.

Google Scholar

Skotland, T., Sandvig, K., and Llorente, A. (2017). Lipids in exosomes: current knowledge and the way forward. Prog. Lipid Res. 66, 30–41. doi: 10.1016/j.plipres.2017.03.001

PubMed Abstract | CrossRef Full Text | Google Scholar

Steck, T. L., and Lange, Y. (2018). Transverse distribution of plasma membrane bilayer cholesterol: picking sides. Traffic 19, 750–760. doi: 10.1111/tra.12586

PubMed Abstract | CrossRef Full Text | Google Scholar

Teixeira, V., Feio, M. J., and Bastos, M. (2012). Role of lipids in the interaction of antimicrobial peptides with membranes. Prog. Lipid Res. 51, 149–177. doi: 10.1016/j.plipres.2011.12.005

PubMed Abstract | CrossRef Full Text | Google Scholar

Tian, J., Nickels, J., Katsaras, J., and Cheng, X. (2016). Behavior of bilayer leaflets in asymmetric model membranes: atomistic simulation studies. J. Phys. Chem. B 120, 8438–8448. doi: 10.1021/acs.jpcb.6b02148

PubMed Abstract | CrossRef Full Text | Google Scholar

Tieleman, D. P., and Marrink, S.-J. (2006). Lipids out of equilibrium: energetics of desorption and pore mediated flip-flop. J. Am. Chem. Soc. 128, 12462–12467. doi: 10.1021/ja0624321

PubMed Abstract | CrossRef Full Text | Google Scholar

Torres-Sánchez, A., Vanegas, J. M., and Arroyo, M. (2015). Examining the mechanical equilibrium of microscopic stresses in molecular simulations. Phys. Rev. Lett. 114:258102.

Google Scholar

Van Der Spoel, D., Lindahl, E., Hess, B., Groenhof, G., Mark, A. E., and Berendsen, H. J. (2005). GROMACS: fast, flexible, and free. J. Comput. Chem. 26, 1701–1718. doi: 10.1002/jcc.20291

PubMed Abstract | CrossRef Full Text | Google Scholar

Van Meer, G., Voelker, D. R., and Feigenson, G. W. (2008). Membrane lipids: where they are and how they behave. Nat. Rev. Mol. Cell Biol. 9, 112–124. doi: 10.1038/nrm2330

PubMed Abstract | CrossRef Full Text | Google Scholar

Vanegas, J. M., Torres-Sánchez, A., and Arroyo, M. (2014). Importance of force decomposition for local stress calculations in biomembrane molecular simulations. J. Chem. Theory Comput. 10, 691–702. doi: 10.1021/ct4008926

PubMed Abstract | CrossRef Full Text | Google Scholar

Wassenaar, T. A., Ingólfsson, H. I., Böckmann, R. A., Tieleman, D. P., and Marrink, S. J. (2015). Computational lipidomics with insane: a versatile tool for generating custom membranes for molecular simulations. J. Chem. Theory Comput. 11, 2144–2155. doi: 10.1021/acs.jctc.5b00209

PubMed Abstract | CrossRef Full Text | Google Scholar

Wassenaar, T. A., Pluhackova, K., Böckmann, R. A., Marrink, S. J., and Tieleman, D. P. (2014). Going backward: a flexible geometric approach to reverse transformation from coarse grained to atomistic models. J. Chem. Theory Comput. 10, 676–690. doi: 10.1021/ct400617g

PubMed Abstract | CrossRef Full Text | Google Scholar