Weak carbohydrate–carbohydrate interactions in membrane adhesion are fuzzy and generic

Abstract

Carbohydrates such as the trisaccharide motif Le^X are key constituents of cell surfaces. Despite intense research, the interactions between carbohydrates of apposing cells or membranes are not well understood. In this article, we investigate carbohydrate–carbohydrate interactions in membrane adhesion as well as in solution with extensive atomistic molecular dynamics simulations that exceed the simulation times of previous studies by orders of magnitude. For Le^X, we obtain association constants of soluble carbohydrates, adhesion energies of lipid-anchored carbohydrates, and maximally sustained forces of carbohydrate complexes in membrane adhesion that are in good agreement with experimental results in the literature. Our simulations thus appear to provide a realistic, detailed picture of Le^X–Le^X interactions in solution and during membrane adhesion. In this picture, the Le^X–Le^X interactions are fuzzy, i.e. Le^X pairs interact in a large variety of short-lived, bound conformations. For the synthetic tetrasaccharide Lac 2, which is composed of two lactose units, we observe similarly fuzzy interactions and obtain association constants of both soluble and lipid-anchored variants that are comparable to the corresponding association constants of Le^X. The fuzzy, weak carbohydrate–carbohydrate interactions quantified in our simulations thus appear to be a generic feature of small, neutral carbohydrates such as Le^X and Lac 2.

---

Introduction

Carbohydrates are omnipresent at cell surfaces as constituents of glycolipids and glycoproteins.^1–3 During cell adhesion, these carbohydrates come into contact with proteins and carbohydrates on apposing cell surfaces. While specific interactions between carbohydrates and proteins are known to play important roles in cell adhesion events, the role of carbohydrate–carbohydrate interactions in these events is less clear.^4–8 About three decades ago, homophilic carbohydrate–carbohydrate interactions of the trisaccharide Lewis^X (Le^X) have been reported to be involved in embryonal cell compaction and aggregation,^9–11 and interactions between long carbohydrate chains have been linked to the species-specific aggregation of marine sponges.¹² In the following decades, carbohydrate–carbohydrate interactions in adhesion have been investigated in a variety of reconstituted or synthetic systems including nanoparticles and surfaces functionalized with carbohydrates,^13–15 atomic force microscopy setups,^16–19 and reconstituted vesicles^20,21 or membranes^22–24 containing glycolipids. While some carbohydrate–carbohydrate interactions have been reported to be strong,^17,25,26 interactions of small, neutral carbohydrates are typically considered to be weak.^27,28 However, the binding association constants, in particular at membrane interfaces, and the structural binding mechanisms are often not known.

In this article, we present detailed results from atomistic molecular dynamics simulations of carbohydrate–carbohydrate interactions in membrane adhesion and in solution for Le^X and the synthetic saccharide Lac 2, which is composed of two lactose units²⁹ (see Fig. 1). Our simulations employ a recent carbohydrate force field³⁰ that allows a more faithful representation of carbohydrate–carbohydrate interactions^30–32 and exceed the times and system sizes in previous simulation studies of carbohydrate–carbohydrate interactions in solution^33–35 by orders of magnitude. Le^X has been investigated extensively as a model system for carbohydrate–carbohydrate interactions,^{13,14,16,19–21,23} and experimental data available from these investigations are central to corroborate our simulation results. In our Le^X glycolipids, the Le^X trisaccharide is connected via a lactose disaccharide and a glycerol linker to lipids tails (see Fig. 1). In our Lac 2 glycolipids, the Le^X trisaccharide is replaced by another lactose disaccharide, which allows to compare the carbohydrate–carbohydrate interactions of Le^X to those of the common saccharide lactose. From simulations of soluble pairs of Le^X and Lac 2, we obtain association constants K_a of the order of 10 M⁻¹, which agrees with a K_a value of Le^X derived from weak affinity chromatography experiments.^16,36 From simulations of pairs of Le^X and Lac 2 glycolipids at apposing membrane surfaces, we obtain comparable association constants K_trans for the Le^X and Lac 2 glycolipids that strongly decrease with increasing membrane separation. For the membrane separation and thermal roughness of membrane multilayers with 10 mol% Le^X glycolipids measured in neutron scattering experiments,²³ we determine an adhesion energy per area of the order of 10 μJ m⁻² from our K_trans values, in agreement with the adhesion energy per area reported for vesicles that contain 10 mol% of Le^X glycolipids.²⁰ The average force on bound Le^X glycolipid complexes determined in our simulations increases with increasing membrane separation up to a maximum value of about 20 pN, which agrees with the Le^X–Le^X unbinding force obtained from atomic force microscopy experiments.¹⁶ The agreement with experimental results indicates that our simulations provide a realistic, detailed picture of weak carbohydrate–carbohydrate interactions in solution as well as in membrane adhesion. A striking feature is that the carbohydrate–carbohydrate interactions are fuzzy, i.e. both soluble and lipid-anchored variants of Le^X and Lac 2 interact in our simulations via a large variety of diverse, bound conformations.

Fig. 1 Structures of the soluble and lipid-anchored saccharides investigated in our simulations.

Results

Interactions of soluble carbohydrates

We first consider the interaction of two Le^X trisaccharides in solution and compare this Le^X–Le^X pair interaction to the interaction of two Lac 2 tetrasaccharides, which are composed of two lactose units²⁹ (see Fig. 1). Standard carbohydrate force fields lead to osmotic pressures for solutions of neutral carbohydrates that are systematically too low compared to experimental values. This underestimation of the osmotic pressure of the carbohydrate solutions results from an overestimation of attractive carbohydrate–carbohydrate interactions.^30,31 To avoid unrealistically attractive carbohydrate–carbohydrate interactions, we have used the GLYCAM06^TIP5P_OSMOr14 force field, in which the van der Waals parameters for saccharide–saccharide interactions of the standard force field GLYCAM06 have been reparametrized to correctly reproduce experimentally measured osmotic pressures.³⁰ The GLYCAM06^TIP5P_OSMOr14 force field employs the TIP5P water model because this water model leads to more reliable carbohydrate–carbohydrate interactions in GLYCAM06, compared to the standard TIP3P water model.^32,37 Using graphics processing units (GPUs) and the software AMBER GPU,^38,39 we have generated 50 simulation trajectories with a length of 2.0 μs for two Le^X molecules in a periodic simulation box of volume V = 131.5 nm³, and 40 trajectories with a length of 1 μs or close to 1 μs for two Lac 2 molecules in a simulation box of volume V = 260.5 nm³, at the simulation temperature 30 °C. Our total simulation times are 100 μs for the Le^X pair and 39.5 μs for the Lac 2 pair, which greatly exceed the total simulation times up to 40 ns in previous simulation studies of Le^X–Le^X pair interactions in solution^33,34 and the total simulation time of a few ns for pair interactions of trisaccharide epitopes from marine sponges.³⁵

In our simulations, we observe thousands of interaction events in which the two Le^X molecules or the two Lac 2 molecules are in contact. These interaction events are separated by longer or shorter trajectory parts in which the two molecules are not in contact. Fig. 2(a) and (b) display pair conformations of Le^X and Lac 2 in which the two molecules exhibit at least 20 or 50 contacts of non-hydrogen atoms, respectively. The shown pair conformations are randomly selected from the simulation frames of our trajectories. One of the carbohydrate molecules is aligned in the pair conformations and represented in blue colors, while the other molecule is represented in red/yellow colors. The clouds of red/yellow molecules around the aligned blue molecules in these conformations illustrate that the carbohydrate–carbohydrate interactions are ‘fuzzy’,^40,41i.e. the two molecules interact in broad ensembles of conformations, rather than via a single binding conformation. For both Le^X and Lac 2, the ensembles of conformations with at least 50 contacts are narrower than the ensembles of conformations with at least 20 contacts. In conformations with at least 50 contacts, the two Le^X molecules tend to stack above each other in different orientations, and the two Lac 2 molecules tend to align parallel or anti-parallel. However, the probability distributions of contact numbers in Fig. 2(c) illustrate that pair conformations with 50 or more contacts of non-hydrogen atoms are rather rare and not typical. The probability distributions decrease monotonously with increasing number of contacts.

Fig. 2 (a) and (b) Randomly selected pair conformations of two Le^X and two Lac 2 molecules with at least 20 or at least 50 contacts between non-hydrogen atoms within a distance less than 0.45 nm, respectively. One of the molecules is aligned in the 50 pair conformations and represented in blue colors, while the other molecule is represented in red/yellow colors. In the aligned Le^X molecules, fucose is represented in dark blue, galactose in light blue, and N-acytylglucosamine in cyan. In the other Le^X molecules, these monosaccharide units are represented in red, orange, and yellow, respectively. In the aligned Lac 2 molecules, the terminal galactose is represented in dark blue, the adjacent glucose in light blue, and the remaining galactose and glucose in cyan. In the other Lac 2 molecules, these monosaccharides are shown in red, orange, and yellow. (c) Probability distributions of the number of contacts between non-hydrogen atoms obtained from our simulations of two soluble Le^X or two soluble Lac 2 molecules. (d) Average lifetime of interaction events as a function of the maximum number of contacts of the interaction events. Interaction events are consecutive stretches of simulation frames at intervals of 0.1 ns with nonzero contacts of the two molecules. The error bars represent the standard deviations of the observed lifetimes. (e) Radial distribution functions g(r) of two soluble Le^X or Lac 2 molecules with center-of-mass distance r.

The interaction events of the two Le^X molecules or the two Lac 2 molecules can be characterized by their lifetime and by the maximum number of contacts of the events. These interaction events are obtained from our simulation trajectories as consecutive stretches of frames at intervals of 0.1 ns with nonzero contacts of the two molecules. Fig. 2(d) shows that the average lifetime of the interaction events increases with the maximum number of contacts observed during the event. With average lifetimes in the nanoseconds range, the interactions of the two Le^X or the two Lac 2 molecules are rather short-lived. Nonetheless, the radial distribution functions in Fig. 2(e) indicate that the interactions are attractive. The maxima of the radial distribution functions at center-of-mass distances of about 0.8 nm for Le^X and 0.6 nm for Lac 2 are significantly larger than the value 1 for a non-interacting ideal solution.

Quantifying the attractive interactions of the two Le^X or two Lac 2 molecules requires distinguishing bound and unbound states. This distinction is somewhat arbitrary because of the fuzzy interactions of the carbohydrates. The probability distributions of carboyhydrate–carbohydrate contact numbers in Fig. 2(c) are monotonously decreasing and, thus, not bimodal as required for a clear distinction of two states. Table 1 presents association constants of two Le^X or two Lac 2 molecules calculated for different cutoffs n_c of the maximum number of contacts of interaction events. In these calculations, only interaction events with a maximum number of contacts larger or equal to the cutoff n_c are taken to be binding events. The probability P_b that the two Le^X or two Lac 2 molecules are bound has been determined from the total duration of the binding events, and the association constants from K_a = VP_b/P_u where P_u = 1 − P_b is the probability that the molecules are unbound, and V is the volume of the simulation box. The K_a values in Table 1 slightly decrease with increasing contact cutoff n_c for binding events. For Le^X, a K_a value of 10 M⁻¹ has been obtained from weak affinity chromatography experiments,³⁶ which is of the same order of magnitude as the values derived from our simulations.

Table 1 Association constants K_a in units of M⁻¹ for different cutoffs n_c for the contact number of binding events

	n c = 5	n c = 10	n c = 20
Le^X	6.4 ± 0.3	5.7 ± 0.3	4.5 ± 0.3
Lac 2	13.2 ± 1.0	12.3 ± 1.0	10.7 ± 0.9

---

Interactions of lipid-anchored carbohydrates

To investigate the interactions of two lipid-anchored Le^X or two lipid-anchored Lac 2 molecules, we have performed simulations of Le^X and Lac 2 glycolipids embedded in POPC lipid membranes. Our Le^X and Lac 2 glycolipids have the same lipid tails as POPC, and carbohydrate tips that are connected to these lipid tails by a glycerol linker group (see Fig. 1). The carbohydrate tip of the Le^X glycolipid consists of the Le^X trisaccharide and an additional lactose disaccharide as spacer between Le^X and the glycerol linker. The Lac 2 glycolipid has the linear Lac 2 tetrasaccharide as carbohydrate tip. The force field of our simulations combines the GLYCAM06^TIP5P_OSMOr14 carbohydrate force field^30,42 for the TIP5P water model with the AMBER Lipid14 force field⁴³ for lipid membranes. Because simulations of AMBER Lipid14 POPC membranes in TIP5P water lead to an unreasonably small area per lipid, we have rescaled the Lennard-Jones interactions between the TIP5P water molecules and the lipid headgroup atoms to obtain the same area per lipid as in standard AMBER Lipid14 simulations with the TIP3P water model (see Methods).

We quantify the interactions of two Le^X or two Lac 2 glycolipids at apposing membrane surfaces in a system that consists of a single lipid bilayer with one glycolipid anchored in each monolayer (see Fig. 3). In this system, the two glycolipids in the different monolayers interact due to the periodic boundary conditions of the simulation box, and the separation of the membrane monolayers can be adjusted by varying the number of water molecules in the simulation box. The values for the membrane separation l given in Fig. 3 correspond to the separation from membrane midplane to membrane midplane and, thus, to the height of the simulation box. At each membrane separation, we have generated 10 trajectories with a length of 3 μs for the Le^X system and a length of 1 μs for the Lac 2 system at the temperature 30 °C. The total simulation times at each membrane separation thus are 30 μs and 10 μs for the Le^X and Lac 2 systems, respectively. The membranes contain in each monolayer 35 lipids besides the single glycolipid and have an area A of 23.3 nm². The height of the simulation box l increases with the number of water molecules n_w as l ≃ 3.8 nm + 0.013n_w nm. The thickness of the water layer in the simulations thus is about l − 3.8 nm.

Fig. 3 (Top) Membrane conformations with two unbound or bound Le^X glycolipids. The Le^X glycolipids are anchored in the different monolayers of the membrane and interact because of the periodic boundary conditions of the simulation box. The height of the simulation box corresponds to the membrane separation from bilayer midplane to midplane. Each membrane monolayer contains 35 POPC lipids, which have the same lipid tails as the Le^X glycolipids. The fucose and galactose at the branched tip of the Le^X glycolipids are represented in red and orange, and the remaining three monosaccharide units in yellow. (Bottom) 50 randomly selected complexes of the carboyhydrate tips of the Le^X glycolipids at different membrane separations. The selected complexes exhibit at least 10 contacts between non-hydrogen atoms of the two carbohydate tips. The carbohydrate tip of the lower Le^X glycolipid is aligned in the 50 complexes and represented in blue colors, while the carbohydrate tip of the upper glycolipid is represented in red/yellow colors. The Le^X motif of the carbohydrate tips are represented in the same colors as in Fig. 2(a). The lactose disaccharides of the carbohydrate tips, which are located between the Le^X trisaccharide and the linker group of the glycolipid, are represented in light blue and light yellow, respectively.

The interactions of the glycolipids strongly depend on the membrane separation. For the membrane separations l = 5.5, 6.0, 6.5, and 7.0 nm, 50 randomly selected complexes of the Le^X glycolipid tips with at least 10 contacts of non-hydrogen atoms are displayed at the bottom of Fig. 3. The carbohydrate tip of the lower Le^X glycolipid is aligned in the 50 complexes and represented in blue colors, while the carbohydrate tip of the upper glycolipid is represented in red/yellow colors. The clouds of red/yellow carbohydrates illustrate that the interactions of lipid-anchored Le^X are fuzzy, similar to soluble Le^X and Lac 2 (see Fig. 2). The overlap of the cloud of the upper, red/yellow carbohydrates with the lower, blue carbohydrate decreases with increasing membrane separation. At the membrane separation 5.5 nm, the Le^X glycolipids interact via their entire carbohydrate tips. At the separation 6.0 nm, the interactions are limited to the Le^X trisaccharide of the glycolipid tip, and at the membrane separations 6.5 nm and 7.0 nm, the interactions are further restricted to the galactose and fucose monosaccharides at the branched end of the Le^X glycolipid. The decrease of interactions with increasing separation is also reflected in the probability distributions of contact numbers shown in Fig. 4(a) and in the average lifetime of the interaction events for different maximum numbers of contacts in Fig. 4(b). At the smallest membrane separation 5.5 nm, complexes of Le^X glycolipids can exhibit up to 60 and more contacts of non-hydrogen atoms (see Fig. 4(a)), and average lifetimes up to 50 ns for interaction events with a maximum number of 60 contacts (see inset of Fig. 4(b)), which are about one order of magnitude larger than the average lifetimes for interaction events of soluble Le^X molecules with the same maximum number of contacts. At the membrane separations 6.0 and 6.5 nm, the overall contact numbers and lifetimes of interaction events are significantly smaller.

Fig. 4 (a) Probability distributions of the number of contacts between two Le^X glycolipids at different membrane separations. (b) Average lifetime of interaction events at different membrane separations as a function of the maximum contact number of the events.

Analogous to soluble carbohydrates, the binding association constants K_trans = AP_b/(1 − P_b) of the glycolipids in the different membrane monolayers can be determined from the probability P_b that the two Le^X or two Lac 2 glycolipids are bound. The binding constants shown in Fig. 5 are calculated for binding events with a maximum number of at least n_c = 5 contacts of non-hydrogen atoms. For the larger binding cutoff n_c = 10, the K_trans values of the two Le^X glycolipids are about 10% smaller than the values in Fig. 5 at the membrane separations 5.5 and 6.0 nm, and the values of the Lac 2 glycolipids are about 15% smaller at these separations. The K_trans values decrease with increasing membrane separation. For membrane separations larger than about 7.5 nm, the glycolipids cannot form contacts.

Fig. 5 Binding constant K_trans of two Le^X and two Lac 2 glycolipids versus membrane separation, calculated for binding events with a maximum number of at least n_c = 5 contacts of non-hydrogen atoms.

The binding constant K_trans can be related to membrane adhesion energies, which have been measured for membrane vesicles that contain 10 mol% of Le^X glycolipids.^20,27 For two apposing, large membrane surfaces of area A that contain a total number of N_t glycolipids, the free energy difference for forming the nth bond of the glycolipids is (see Methods)

ΔG_n = −k_BT ln[K_trans(N_t − n + 1)²/nA] (1)

The free energy differences ΔG_n are negative and, thus, favourable, from bond 1 until the equilibrium number n_eq of bonds. For bond numbers n > n_eq, the free energy difference ΔG_n is positive and, thus, unfavorable for binding. The adhesion free energy g_ad per area now can be calculated by summing up the free energy differences ΔG_n from bond 1 to bond n_eq:For an area per lipid of 0.65 nm² measured in our simulations, the area of a membrane surface that contains N_t glycolipids at a concentration of 10 mol% is A ≃ 6.5N_t nm². From eqn (1) and (2) and the values of K_trans for the Le^X glycolipids in Fig. 5, we obtain the adhesion free energies g_ad = 320 ± 60, 150 ± 20, 28 ± 5, and 5 ± 2 μJ m⁻² at the membrane separations l = 5.5, 6.0, 6.5, and 7.0 nm respectively. For lipid vesicles that contain 10 mol% of Le^X glycolipids, an adhesion free energy per area of 27 ± 2 μJ m⁻² has been reported,²⁰ which is comparable to the adhesion free energy obtained from our simulations with membrane separation 6.5 nm.

Forces on lipid-anchored carbohydrates in trans-direction

The binding of glycolipids in our simulations is associated with deviations of the glycolipids relative to the surrounding lipids. These deviations in the trans-direction perpendicular to the membrane surface result from forces on bound glycolipid complexes. Fig. 6(a) illustrates distributions of trans-deviations between the center of mass of the linker group of a Le^X glycolipid (see Fig. 1) and the center of mass of all lipid head groups in the same monolayer as the glycolipid. The trans-deviations d are calculated from the simulation frames of our trajectories at intervals of 0.1 ns. We obtain two values of d per simulation frame for the two glycolipids relative to the monolayer in which they are embedded. An increase in d indicates glycolipid motion away from the membrane midplane. With increasing membrane separation, the distributions for bound Le^X glycolipids deviate more and more from the distribution for unbound Le^X, which reflects increasing forces. The distribution of trans-deviations d of unbound Le^X glycolipids shown in Fig. 6(a) is calculated from our simulation trajectories at the membrane separation 8.0 nm, at which Le^X bonds do not occur, and can be approximated by a Gaussian distribution exp[−V(d)/k_BT] with . The trans-deviations d of unbound Le^X glycolipids thus can be described by a harmonic potential V(d) with force constant k and mean extension d_u, which can be determined from the standard deviation σ and mean d̄ of the Gaussian as k = k_BT/σ² = 94 ± 4 pN nm⁻¹ and d_u = d̄ = −0.31 ± 0.10 nm. The distributions of trans-deviations of bound Le^X glycolipids in Fig. 6(a) are calculated from our simulation trajectories at the membrane separations 5.5, 6.0, 6.5, and 7.0 nm, for binding events with a maximum number of at least n_c = 5 contacts of non-hydrogen atoms. The average force f = k(d_b − d_u) on bound Le^X glycolipids at the membrane separations l = 5.5, 6.0, 6.5, and 7.0 nm then can be calculated from the difference between the mean trans-deviations d_b = −0.26 ± 0.01, −0.22 ± 0.01, −0.16 ± 0.01, and −0.08 ± 0.02 nm of the bound glycolipids at these membrane separations and the mean trans-deviation d_u of the unbound glycolipids. The force f on bound Le^X glycolipids increases with increasing membrane separation up to a value of 21.7 ± 2.4 pN at the separation 7.0 nm (see Fig. 6(b)). This maximal force value agrees with the unbinding force 20 ± 4 pN of two Le^X molecules obtained from atomic force microscopy experiments.¹⁶ For bound Lac 2 glycolipids, we obtain a maximal force of 14.7 ± 3.5 pN at the separation 7.0 nm, which is about of the same magnitude as the maximal force sustained by the Le^X complexes.

Fig. 6 (a) Distributions of trans-deviations of Le^X linker groups relative to the surrounding lipids. The trans-deviations are calculated as the difference between the center of mass of the Le^X glycolipid linker group (see Fig. 1) and the center of mass of all lipid head groups in the same monolayer as the glycolipid. These trans-deviations of Le^X in the direction perpendicular to the membrane plane are determined from the simulation trajectories of the system illustrated in Fig. 3. (b) Forces on bound and unbound Le^X glycolipids at the different membrane separations. The trans-deviations and forces of bound glycolipids are obtained from the simulation frames of binding events with a maximum number of at least n_c = 5 contacts of non-hydrogen atoms. Deviations to force values obtained for the cutoff n_c = 10 are smaller than the error bars. Forces an unbound glycolipids are calculated from simulation frames with zero contacts between the glycolipids.

The forces on bound Le^X glycolipids lead to an adhesion pressure between the membranes. Fig. 7 illustrates the adhesion pressure p of membranes that contain 10 mol% of Le^X glycolipids as a function of the membrane separation. The adhesion pressure is estimated as p = P_b f/A where P_b is the probability that a Le^X glycolipid is bound at the concentration 10 mol%, f is the average force on the bound glycolipid, and A ≃ 6.5 nm² is the average membrane area of membrane patch with a single glycolipid at this concentration (see above). The negative pressure values for membrane separations l of 7.0 nm and smaller, at which the glycolipids can bind, indicate membrane attraction. From integration of the pressure profile along the dashed interpolation line shown in Fig. 7, we obtain adhesion energies for l = 6.0 nm and g_ad ≃ 30 μJ m⁻² for l = 6.5 nm. These adhesion energies per area agree with values g_ad = 150 ± 20 μJ m⁻² and 28 ± 5 μJ m⁻² obtained directly from the binding constants K_trans at the membrane separations l = 6.0 and 6.5 nm (see above), which indicates that average forces f on bound Le^X glycolipids of Fig. 6(b) are consistent with the binding constants K_trans shown in Fig. 5.

Fig. 7 Adhesion pressure p of membranes with 10 mol% of Le^X glycolipids obtained for the force values f on bound Le^X of Fig. 6(b). The dashed interpolation line is added as a guide for the eye and used to estimate adhesion energies via integration (see text). In this integration, the pressure p is taken to be zero at separations l ≥ 7.5 nm.

Discussion and conclusions

The membranes in our simulation systems are essentially planar because of the small size of the membranes, and because the glycolipid in one monolayer interacts with the glycolipid in the other monolayer across the periodic boundary of the simulation box. In larger, experimental systems, in contrast, the membranes exhibit thermally excited shape fluctuations, which lead to a steric repulsion between adjacent membranes.⁴⁴ During membrane adhesion, this steric repulsion needs to be overcome by attractive interactions.⁴⁵ The average separation and thermal roughness of the adhering membranes is determined by the interplay of the attractive interactions and the steric repulsion.⁴⁶ From neutron scattering experiments of DPPC membrane multilayers that contain 10 mol% of Le^X glycolipids,²³ an average membrane separation of l̄ = 7.7 ± 0.1 nm and a relative membrane roughness of ξ_⊥ = 0.73 ± 0.03 nm has been obtained.† Because of the periodicity of the membrane multilayers, the distribution of the local membrane separations l between adjacent membranes can be approximated by the symmetric Gaussian distribution with mean l̄ and standard deviation ξ_⊥. The average membrane separation l̄ obtained from neutron scattering is larger than the membrane separations at which the Le^X glycolipids interact in our simulations. Trans-binding of the glycolipids therefore requires local membrane separations of the fluctuating membranes that are smaller than the average separation of the membranes. The average adhesion energy per area of adjacent membranes can be estimated as , where g_ad(l) is the adhesion energy as a function of the local membrane separation l. From the four values of g_ad(l) at the membrane separations l = 5.5, 6.0, 6.5, and 7.0 nm determined in the section “Interactions of lipid-anchored carbohydrates”, we obtain the estimate ḡ_ad = 7 ± 3 μJ m⁻² for the average separation l̄ and relative membrane roughness ξ_⊥ of the neutron scattering experiments. This estimate of the average adhesion energy per area is comparable in magnitude to the adhesion free energy per area of 27 ± 2 μJ m⁻² reported for adhering membrane vesicles that contain 10 mol% of Le^X glycolipids.²⁰ The Le^X glycolipids embedded in the vesicles have the same carbohydrate tip as the Le^X glycolipids of the neutron scattering experiments and of our simulations. However, the carbohydrate tip of the vesicle system is connected to a ceramide, which contains a different linker between the carbohydrate tip and the lipid tails. Another difference is that the neutron scattering experiments have been performed at the temperature 60 °C to ensure that the DPPC membranes in these experiments are fluid.²³ The Le^X glycolipids of our simulations differ from those of the neutron scattering experiments only in the lipid tails. We have focused on POPC membranes and corresponding glycolipid tails to be able to run simulations of fluid membranes at the temperature 30 °C, which is close to the calibration temperature of the force fields. In principle, membrane tension suppresses shape fluctuations of the membranes and can lead to stronger adhesion. However, the suppression of fluctuations occurs only on lateral length scales larger than the characteristic length ,⁴⁷ which adopts values between 100 and 400 nm for typical membrane tensions σ of a few μN m⁻¹ (ref. 48–50) and typical membrane bending rigidities κ between 10 and 40k_BT.^51,52 These values are significantly larger than the lateral correlation length ξ_∥ of membranes adhering via Le^X glycolipids, which is only a few nanometers for a relative membrane roughness ξ_⊥ of about 0.7 nm.⁵³ On these small length scales, the membrane shape fluctuations are dominated by the bending energy of the membranes and the adhesion energies of the glycolipids, and are not affected by membrane tension.

Table 2 Membrane thickness d_m and lipid diffusion coefficient D from simulations with TIP5P water for different values of the scaling factor α, from simulations with TIP3P water, and from experiments on POPC lipid membranes

α	d m [nm]	D [μm^2 s^−1]
1.2	3.82 ± 0.03	3.27 ± 0.10
1.3	3.67 ± 0.01	4.15 ± 0.08
1.4	3.51 ± 0.01	5.57 ± 0.16
1.45	3.43 ± 0.01	5.85 ± 0.15
1.5	3.36 ± 0.01	6.77 ± 0.11
1.55	3.31 ± 0.01	7.31 ± 0.19
TIP3P	3.54 ± 0.01	5.45 ± 0.19
Exp.	3.68 (ref. 68)	10.7 (ref. 69)

---

The fuzzy interactions and comparable magnitude of the association constants of Le^X and Lac 2 obtained in our simulations indicate that the interactions of small, neutral carbohydrates such as Le^X and Lac 2 are rather generic and not dependent on specific, structural aspects of the carbohydrates. The good agreement to experimental results for the association constant of soluble Le^X,³⁶ adhesion energies of membranes with Le^X glycolipids,²⁰ and maximally sustained forces of Le^X complexes¹⁶ suggests that our simulations provide a realistic, detailed picture of weak carbohydrate–carbohydrate interactions in solution as well as in membrane adhesion. The fuzzy binding reduces the loss of rotational and translational entropy of the molecules during binding,⁴⁰ because binding can occur for a large variety of different relative orientations of the saccharides, in contrast to e.g. binding via specific hydrogen-bond patterns as suggested previously for Le^X based on simulations on short timescales up to 40 ns.^33,34 The fuzzy binding results from a subtle interplay between the rotational and translational entropy of the saccharides and the van der Waals, hydrogen bond, and hydrophobic interactions of the saccharides in the various binding conformations.

We have investigated the binding of Le^X in the absence of Ca²⁺. Several groups have reported that Le^X binding depends on Ca²⁺,^{13,14,19–21,54,55} whereas other groups have observed no dependence on Ca²⁺.^16,23 As pointed out by Kunze et al.,²¹ the Ca²⁺ concentrations used by most groups are of the order of 10 mM and, thus, greatly beyond physiological Ca²⁺ concentrations. In vesicle adhesion experiments, Kunze et al.²¹ observed a rather small increase of the number of bound vesicles for a physiological Ca²⁺ concentration of 0.9 mM, compared to experiments in the absence of Ca²⁺. However, a strong increase of the number of bound vesicles in the experiments occurred for a Ca²⁺ concentration of 10 mM. In atomic force microscopy experiments of Le^X unbinding,¹⁶ in contrast, the same unbinding force of about 20 ± 4 pN has been obtained both in the absence of Ca²⁺ and for a Ca²⁺ concentration of 10 mM. Overall, these experimental results suggest that the binding of Le^X is not strongly affected at least by physiological concentrations of Ca²⁺.

Methods

Simulations of soluble carbohydrates

System setup. We have used the GLYCAM06^TIP5P_OSMOr14 carbohydrate force field^30,42 in our simulations of soluble pairs of Le^X and Lac 2 in water. Initial structures of the Le^X trisaccharides and Lac 2 tetrasacchardies were created with the Glycam Carbohydrate Builder program⁵⁶ and solvated in truncated octahedral simulation boxes with 4287 TIP5P water molecules for the Le^X pair and with 8504 TIP5P water molecules for the Lac 2 pair. In the initial conformations, the two saccharides were placed in the simulation boxes such that they were not in contact. We have subsequently minimized the simulation systems in 5000 minimization steps of steepest decent and additional 5000 steps of the conjugent gradient algorithm. The systems were then heated from the temperature 0 K to 303 K at constant volume in 50 000 integration time steps of 2 fs with temperature control by a Langevin thermostat⁵⁷ with collision frequency γ = 1.0 ps⁻¹.

Production runs. After equilibration for 2 ns at 303 K, we have generated 50 independent trajectories for the Le^X pair and 40 trajectories for the Lac 2 pair with a 2 fs integration step in AMBER 14 and 16 GPU^38,39 using the Monte-Carlo barostat⁵⁸ and a Langevin thermostat with collision frequency γ = 1.0 ps⁻¹ to keep the temperature at 303 K and the pressure at 1 bar. On these trajectories, the lengths of bonds that contain hydrogens were restrained with the SHAKE algorithm,^59,60 non-bonded interactions were truncated at a cutoff value of 1 nm, and the Particle Mesh Ewald algorithm (PME)^61,62 was used to treat all electrostatic interactions. The 50 simulation trajectories for the Le^X pair have a length of 2.0 μs, and the 40 trajectories for the Lac 2 pair have a length of 1 μs or close to 1 μs. The total simulation times of these trajectories are 100 μs for the Le^X system and 39.5 μs for the Lac 2 system.

Analysis of trajectories. We have identified interactions events of the two Le^X or two Lac 2 molecules along the simulation trajectories as consecutive stretches of simulation frames at intervals of 0.1 ns with nonzero contacts of the molecules. These interaction events are separated by stretches of simulation frames with zero contacts and can be characterized by their lifetime and by the maximum number of contacts during the events. The contacts are defined as contacts between non-hydrogen atoms of the two molecules within a distance of less than 0.45 nm. We consider interaction events with a maximum number of contacts that is larger or equal to a cutoff number n_c as binding events. For the cutoff numbers n_c = 5, 10, and 20, we have obtained 7253, 4820, and 2331 binding events of the two Le^X molecules on all trajectories, and 2369, 1573, and 823 binding events of the two Lac 2 molecules. We have thus observed dozens of binding and unbinding events on each trajectory, with binding and unbinding times that are significantly smaller than the trajectory lengths (see also Fig. 2(d)). To ensure independence from the initial, unbound conformations of the trajectories, we have discarded the first 100 ns on all trajectories in our calculations of the binding probablity P_b of the two molecules, which is defined as the fraction of simulation frames belonging to binding events. We have calculated P_b for each trajectory and have determined the overall value and error of P_b as mean and error of the mean of the values for all trajectories. The association constants K_a reported in Table 1 were calculated from these binding probabilities via the relation K_a = VP_b/(1 − P_b) where V is the simulation box volume.⁶³ The errors of K_a are calculated by error propagation from the errors of P_b. The errors of the probability distributions and radial distribution functions in Fig. 2(c) and (e) are calculated as error of the mean of the corresponding quantities for the individual trajectories.

Simulations of lipid-anchored saccharides

System setup. We have generated the initial structures of the POPC lipid membranes with the CHARMM-GUI program.⁶⁵ For our simulations with glycolipids, one lipid in each monolayer has been replaced by a Le^X or Lac 2 glycolipid, which have the same lipids tails as POPC (see Fig. 1). Following ref. 43, we have performed the initial minimization and equilibration steps of all membrane systems as follows: we have first performed a minimization of the water molecules for fixed lipids and glycolipids in 2500 minimization steps of steepest descent and subsequent 2500 steps of the conjugent gradient algorithm. The lipids and glycolipids have been fixed by harmonic constraints with a force constant of 500 kcal mol⁻¹ Å⁻¹ in this minimization. We have next removed the harmonic constraints, and have repeated the same minimization steps for the complete systems. The subsequent heating of the systems has been performed in three steps: (1) heating from 0 K to 100 K at constant volume with harmonic constraints on lipids and glycolipids with a force constant of 20 kcal mol⁻¹ Å⁻¹; (2) heating from 100 K to 200 K with a reduced force constant of 10 kcal mol⁻¹ Å⁻¹ of the harmonic constraints on lipids and glycolipids; and (3) heating from 200 K to 303 K at constant pressure and a membrane tension of zero with the same harmonic constraints as in the second step using a semi-isotropic pressure coupling and the Berendsen barostat⁶⁶ with a pressure relaxation time of 3 ps. Each heating step consist of 10 000 MD integration steps of length 2 fs with temperature control by a Langevin thermostat with a collision frequency of 5.0 ps⁻¹.

Rescaling of Lennard-Jones interactions between water and lipid headgroups. We have used the GLYCAM06^TIP5P_OSMOr14 carbohydrate force field^30,42 for the carbohydrates and the AMBER Lipid14 force field⁴³ for the lipids of our MD simulations of POPC membranes with glycolipids. Simulations of AMBER Lipid14 POPC membranes in TIP5P water lead to an unreasonably small area per lipid of 0.514 ± 0.002 nm² (see Fig. 8) and to density profiles that deviate significantly from profiles obtained from simulations in the standard TIP3P water model (see Fig. 9), which has been used in the parametrization of the AMBER Lipid14 force field.⁴³ We have therefore rescaled the well depth of the Lennard-Jones interactions between the TIP5P water molecules and the Lipid14 lipid headgroup atoms by a scaling factor α in order to obtain the same area per lipid as in simulations of POPC membranes with TIP3P water. We chose to rescale the Lennard-Jones interactions between water and lipid headgroups because the density profiles of AMBER Lipid14 POPC membranes in TIP5P water show a smaller overlap between the water and lipid head group regions, compared to TIP3P water (see Fig. 9). This smaller overlap likely results from weaker Lennard-Jones interactions, and not from different atom sizes, because the Lennard-Jones radius 3.502 Å of the TIP5P oxygen atom is in fact smaller than the radius 3.53 Å of the TIP3P oxygen atom. Therefore, we have only rescaled the Lennard-Jones well-depth ε for the interaction between TIP5P water and the lipid head group atoms by a scaling factor α.

Fig. 8 Area per lipid for a Lipid14 POPC membrane in TIP5P water as a function of the scaling factor α for well-depth of the Lennard-Jones interactions between TIP5P water molecule and the lipid head group atoms. The dashed horizontal lines represent the area per lipid values of POPC membranes for Lipid14 in TIP3P water and from experiments.⁶⁴ The dashed line through the data points is a guide for the eye. Errors or the simulation data are smaller than the point sizes. The error of the experimental value is indicated by they shaded region. The temperature of the simulations and experiments is 30 °C.

Fig. 9 Electron density profiles for a Lipid14 POPC membrane in TIP3P and TIP5P water at the temperature 303 K. The membrane is composed of 128 lipids.

Fig. 8 illustrates simulation results for the area per lipid as a function of the scaling factor α. The membranes in these simulations consists of 128 POPC lipids, and the number of water molecules is 6400. For each value of α, we have generated 10 independent trajectories of length 150 ns with semi-isotropic pressure coupling at a membrane tension of zero and a temperature of 303 K using the same barostat and thermostat settings as in the last heating step of the system setup (see above). We have determined the area per lipid from the last 100 ns of these trajectories, with errors calculated as error of the mean of the values for the individual trajectories. The value α = 1.4 leads to an area per lipid in TIP5P simulations that is close to the area per lipid both in TIP3P simulations and in experiments (see Fig. 8). We have therefore used α = 1.4 in our simulations of lipid-anchored saccharides. For this value of α, the density profile of AMBER Lipid14 POPC membranes in TIP5P water (not shown) is practically identical to density profile in TIP3P water, and the membrane thickness d_m and lateral diffusion coefficient D of the lipids are identical within errors or close to the values obtained in TIP3P simulations (see Table 2). We have determined the bilayer thickness as the distance between the electron density peaks of the lipid head groups, and the lateral diffusion constant from the relation D = MSD(t)/(4t) where MSD(t) is the mean-squared-displacement of a lipid molecule at time t. To obtain MSD(t), we have first removed the center of mass motion of each leaflet to eliminate the ‘caterpillar effect’⁶⁷ and have divided our trajectories into 20 ns fragments. We have then calculated MSD(t) from the MSD profiles of single lipids by averaging over all lipids and all trajectory fragments. The diffusion coefficients in Table 2 are calculated from linear fits in the time intervals from t = 10 ns to 20 ns in which MSD(t) approaches a constant slope.

Production runs. The membranes of our simulations with two Le^X or two Lac 2 glycolipids are composed of 35 POPC lipids and one glycolipid in each monolayer. By varying the number of water molecules in the simulation box, we have created several membrane systems that differ in simulation box height. In our simulations with Le^X glycolipids, we have obtained the average box heights l = 5.49, 6.00, 6.51, 7.01, 7.50, and 8.00 nm for the numbers 1264, 1722, 2107, 2493, 2878, and 3263 of TIP5P water molecules, respectively. In our simulations with Lac 2 glycolipids, we have obtained the average box heights l = 5.58, 6.05, 6.54, 7.03, and 7.53 nm for the numbers 1373, 1753, 2132, 2512, and 2891 of water molecules. The height l of the rectangular simulation box corresponds to the separation from membrane midplane to membrane midplane across the periodic boundary of the box in the direction perpendicular to the membrane. After equilibration for 100 ns, we have produced 10 independent trajectories for each system with the software AMBER 16 GPU.^38,39 The trajectories have a length of 3 μs for the Le^X systems and a length of 1 μs for the Lac 2 systems. We have regulated the simulation temperature of 303 K using a Langevin thermostat⁵⁷ with a collision frequency of 5.0 ps⁻¹, and have employed a semi-isotropic pressure coupling with a pressure of 1 bar in all directions, which corresponds to a membrane tension of zero. We have used the Berendsen barostat⁶⁶ with relaxation time τ = 3 ps for the pressure regulation because of the stability of the semi-isotropic pressure coupling in AMBER 16 GPU in combination with this barostat. For large systems as considered here, the weak-coupling scheme of the Berendsen barostat can be expected to lead to results that are essentially equivalent to other barostats.⁷⁰ We have constrained the bond lengths for hydrogen atoms with the SHAKE algorithm^59,60 and have used an integration timestep of 2 fs in all simulations. A cutoff length of 1.0 nm was used in calculating the non-bonded interactions with the Particle Mesh Ewald (PME) algorithm.^61,62

Analysis of trajectories. We have identified interactions events between the carbohydrate tips of the two Le^X or two Lac 2 glycolipids in the same way as described above for the soluble saccharides. For two Le^X glycolipids, we have obtained 1490, 1609, 588, 141 binding events with a maximum contact number of at least n_c = 5 on the trajectories at the membrane separations 5.5, 6.0, 6.5, and 7.0 nm, respectively. For two Lac 2 glycolipids, we have obtained 609, 413, 183, and 34 such binding events on the trajectories at the corresponding membrane separations.

To ensure independence from the initial conformation of the trajectories, we have discarded the first 10% of each trajectory in our calculations of the binding probablity P_b of the two molecules. In analogy to soluble carbohydrates, we have determined P_b and its error as mean and error of the mean of the values for the 10 trajectories at a given membrane separation. The binding constant then follows as K_trans = AP_b/(1 − P_b) where A is the membrane area.⁶³ We have calculated the errors of the probability distributions in Fig. 4(a) and 6(a) and of the forces in Fig. 6(b) as error of the mean of the corresponding quantities for the individual trajectories.

Calculation of adhesion free energies from trans-binding constants of membrane-anchored molecules

The binding constant K_trans of molecules anchored to two apposing membrane surfaces 1 and 2 of area A is related to the on- and off-rate constants of these molecules via

K_trans = k_on/k_off (3)

If the total numbers of the molecules at the two surfaces are N₁ and N₂, up to n ≤ min(N₁,N₂) trans-bonds can be formed. The effective rate for going from a state with n − 1 trans-bonds to a state with n bonds is⁶³

k₊ = k_on(N₁ − n + 1)(N₂ − n + 1)/A (4)

and the effective rate for going back from n bonds to n − 1 is

k₋ = nk_off (5)

The condition of detailed balance implies

P_n−1k₊ = k₋P_n (6)

where P_n is the equilibrium probability of the state with n trans-bonds. The free-energy difference ΔG_n between the states with n and n − 1 bonds is related to the equilibrium probabilities via

exp[−ΔG_nk_BT] = P_n/P_n−1 (7)

From these equations, we obtainThe adhesion free energy g_ad per area then can be calculated by summing up the free energy differences ΔG_n from bond 1 to bond n_eq where n_eq is the equilibrium number of bonds at which ΔG_n changes sign (see eqn (2)).

Conflicts of interest

There are no conflicts to declare.

Acknowledgements

Financial support of the International Max Planck Research School (IMPRS) on Multiscale Bio-Systems and by the German Research Foundation (DFG) via Emmy Noether grant SCHN 1396/1 is gratefully acknowledged. We would like to thank Mark Santer for helpful discussions. Open Access funding provided by the Max Planck Society

References

1. B. Alberts, A. Johnson, J. Lewis, D. Morgan, M. Raff, K. Roberts, P. Walter, J. Wilson and T. Hunt, Molecular Biology of the Cell, Garland Science, New York, 6th edn, 2014.
2. A. Varki, Nature, 2007, 446, 1023–1029.
3. J. W. Dennis, I. R. Nabi and M. Demetriou, Cell, 2009, 139, 1229–1241.
4. R. L. Schnaar, Arch. Biochem. Biophys., 2004, 426, 163–172.
5. F. T. Liu and G. A. Rabinovich, Nat. Rev. Cancer, 2005, 5, 29–41.
6. J. Arnaud, A. Audfray and A. Imberty, Chem. Soc. Rev., 2013, 42, 4798–4813.
7. A. Varki, Glycobiology, 2017, 27, 3–49.
8. J. Poole, C. J. Day, M. von Itzstein, J. C. Paton and M. P. Jennings, Nat. Rev. Microbiol., 2018, 16, 440–452.
9. I. Eggens, B. Fenderson, T. Toyokuni, B. Dean, M. Stroud and S. Hakomori, J. Biol. Chem., 1989, 264, 9476–9484.
10. N. Kojima, B. A. Fenderson, M. R. Stroud, R. I. Goldberg, R. Habermann, T. Toyokuni and S.-I. Hakomori, Glycoconjugate J., 1994, 11, 238–248.
11. K. Handa, T. Takatani-Nakase, L. Larue, M. P. Stemmler, R. Kemler and S.-I. Hakomori, Biochem. Biophys. Res. Commun., 2007, 358, 247–252.
12. G. N. Misevic, J. Finne and M. M. Burger, J. Biol. Chem., 1987, 262, 5870–5877.
13. J. M. de la Fuente, A. G. Barrientos, T. C. Rojas, J. Rojo, J. Cañada, A. Fernández and S. Penadés, Angew. Chem., Int. Ed., 2001, 40, 2257–2261.
14. M. J. Hernaiz, J. M. de la Fuente, A. G. Barrientos and S. Penades, Angew. Chem., Int. Ed., 2002, 41, 1554–1557.
15. J. M. de la Fuente, P. Eaton, A. G. Barrientos, M. Menendez and S. Penades, J. Am. Chem. Soc., 2005, 127, 6192–6197.
16. C. Tromas, J. Rojo, J. M. de la Fuente, A. G. Barrientos, R. Garcia and S. Penades, Angew. Chem., Int. Ed., 2001, 40, 3052–3055.
17. I. Bucior, S. Scheuring, A. Engel and M. M. Burger, J. Cell Biol., 2004, 165, 529–537.
18. B. Lorenz, L. Álvarez de Cienfuegos, M. Oelkers, E. Kriemen, C. Brand, M. Stephan, E. Sunnick, D. Yüksel, V. Kalsani, K. Kumar, D. B. Werz and A. Janshoff, J. Am. Chem. Soc., 2012, 134, 3326–3329.
19. H. Witt, F. Savic, M. Oelkers, S. I. Awan, D. B. Werz, B. Geil and A. Janshoff, Biophys. J., 2016, 110, 1582–1592.
20. C. Gourier, F. Pincet, E. Perez, Y. M. Zhang, Z. Y. Zhu, J. M. Mallet and P. Sinay, Angew. Chem., Int. Ed., 2005, 44, 1683–1687.
21. A. Kunze, M. Bally, F. Höök and G. Larson, Sci. Rep., 2013, 3, 1452.
22. Z. W. Yu, T. L. Calvert and D. Leckband, Biochemistry, 1998, 37, 1540–1550.
23. E. Schneck, B. Deme, C. Gege and M. Tanaka, Biophys. J., 2011, 100, 2151–2159.
24. V. M. Latza, B. Deme and E. Schneck, Biophys. J., 2020, 118, 1602–1611.
25. C. J. Day, E. N. Tran, E. A. Semchenko, G. Tram, L. E. Hartley-Tassell, P. S. K. Ng, R. M. King, R. Ulanovsky, S. McAtamney, M. A. Apicella, J. Tiralongo, R. Morona, V. Korolik and M. P. Jennings, Proc. Natl. Acad. Sci. U. S. A., 2015, 112, E7266–E7275.
26. Y. Yu, T. Tyrikos-Ergas, Y. Zhu, G. Fittolani, V. Bordoni, A. Singhal, R. J. Fair, A. Grafmueller, P. H. Seeberger and M. Delbianco, Angew. Chem., Int. Ed., 2019, 58, 13127–13132.
27. F. Pincet, T. Le Bouar, Y. M. Zhang, J. Esnault, J. M. Mallet, E. Perez and P. Sinay, Biophys. J., 2001, 80, 1354–1358.
28. T. R. Patel, S. E. Harding, A. Ebringerova, M. Deszczynski, Z. Hromadkova, A. Togola, B. S. Paulsen, G. A. Morris and A. J. Rowe, Biophys. J., 2007, 93, 741–749.
29. M. F. Schneider, G. Mathe, M. Tanaka, C. Gege and R. R. Schmidt, J. Phys. Chem. B, 2001, 105, 5178–5185.
30. J. Sauter and A. Grafmüller, J. Chem. Theory Comput., 2016, 12, 4375–4384.
31. W. K. Lay, M. S. Miller and A. H. Elcock, J. Chem. Theory Comput., 2016, 12, 1401–1407.
32. R. J. Woods, Chem. Rev., 2018, 118, 8005–8024.
33. C. Gourmala, Y. Luo, F. Barbault, Y. Zhang, S. Ghalem, F. Maurel and B. Fan, J. Mol. Struct.: THEOCHEM, 2007, 821, 22–29.
34. Y. Luo, F. Barbault, C. Gourmala, Y. Zhang, F. Maurel, Y. Hu and B. T. Fan, J. Mol. Model., 2008, 14, 901–910.
35. J. I. Santos, A. Carvalho de Souza, F. J. Cañada, S. Martín-Santamaría, J. P. Kamerling and J. Jiménez-Barbero, ChemBioChem, 2009, 10, 511–519.
36. J. M. de la Fuente and S. Penades, Glycoconjugate J., 2004, 21, 149–163.
37. J. Sauter and A. Grafmüller, J. Chem. Theory Comput., 2015, 11, 1765–1774.
38. R. Salomon-Ferrer, A. W. Goetz, D. Poole, S. Le Grand and R. C. Walker, J. Chem. Theory Comput., 2013, 9, 3878–3888.
39. S. Le Grand, A. W. Goetz and R. C. Walker, Comput. Phys. Commun., 2013, 184, 374–380.
40. P. Tompa and M. Fuxreiter, Trends Biochem. Sci., 2008, 33, 2–8.
41. V. N. Uversky and A. K. Dunker, Anal. Chem., 2012, 84, 2096–2104.
42. K. N. Kirschner, A. B. Yongye, S. M. Tschampel, J. Gonzalez-Outeirino, C. R. Daniels, B. L. Foley and R. J. Woods, J. Comput. Chem., 2008, 29, 622–655.
43. C. J. Dickson, B. D. Madej, A. A. Skjevik, R. M. Betz, K. Teigen, I. R. Gould and R. C. Walker, J. Chem. Theory Comput., 2014, 10, 865–879.
44. W. Helfrich, Z. Naturforsch., A: Phys., Phys. Chem., Kosmophys., 1978, 33, 305–315.
45. R. Lipowsky and S. Leibler, Phys. Rev. Lett., 1986, 56, 2541–2544.
46. J. Steinkühler, B. Rozycki, C. Alvey, R. Lipowsky, T. R. Weikl, R. Dimova and D. E. Discher, J. Cell Sci., 2019, 132, jcs216770.
47. T. R. Weikl, Annu. Rev. Phys. Chem., 2018, 69, 521–539.
48. R. Simson, E. Wallraff, J. Faix, J. Niewohner, G. Gerisch and E. Sackmann, Biophys. J., 1998, 74, 514–522.
49. G. Popescu, T. Ikeda, K. Goda, C. A. Best-Popescu, M. Laposata, S. Manley, R. R. Dasari, K. Badizadegan and M. S. Feld, Phys. Rev. Lett., 2006, 97, 218101.
50. T. Betz, M. Lenz, J.-F. Joanny and C. Sykes, Proc. Natl. Acad. Sci. U. S. A., 2009, 106, 15320–15325.
51. J. F. Nagle, Faraday Discuss., 2013, 161, 11–29.
52. R. Dimova, Adv. Colloid Interface Sci., 2014, 208, 225–234.
53. G.-K. Xu, J. Hu, R. Lipowsky and T. R. Weikl, J. Chem. Phys., 2015, 143, 243136.
54. A. Geyer, C. Gege and R. R. Schmidt, Angew. Chem., Int. Ed., 2000, 39, 3246.
55. G. Nodet, L. Poggi, D. Abergel, C. Gourmala, D. Dong, Y. Zhang, J.-M. Mallet and G. Bodenhausen, J. Am. Chem. Soc., 2007, 129, 9080–9085.
56. Woods Group, GLYCAM Web, Complex Carbohydrate Research Center, University of Georgia, Athens, GA, 2005–2020, http://glycam.org.
57. R. Salomon-Ferrer, D. A. Case and R. C. Walker, Wiley Interdiscip. Rev.: Comput. Mol. Sci., 2013, 3, 198–210.
58. M. P. Allen and D. J. Tildesley, Computer simulation of liquids, Oxford University Press, 2017.
59. S. Miyamoto and P. A. Kollman, J. Comput. Chem., 1992, 13, 952–962.
60. J.-P. Ryckaert, G. Ciccotti and H. J. Berendsen, J. Comput. Phys., 1977, 23, 327–341.
61. U. Essmann, L. Perera, M. L. Berkowitz, T. Darden, H. Lee and L. G. Pedersen, J. Chem. Phys., 1995, 103, 8577–8593.
62. T. Darden, D. York and L. Pedersen, J. Chem. Phys., 1993, 98, 10089–10092.
63. J. Hu, R. Lipowsky and T. R. Weikl, Proc. Natl. Acad. Sci. U. S. A., 2013, 110, 15283–15288.
64. N. Kucerka, M.-P. Nieh and J. Katsaras, Biochim. Biophys. Acta, Biomembr., 2011, 1808, 2761–2771.
65. S. Jo, T. Kim, V. G. Iyer and W. Im, J. Comput. Chem., 2008, 29, 1859–1865.
66. H. J. Berendsen, J. v. Postma, W. F. van Gunsteren, A. DiNola and J. Haak, J. Chem. Phys., 1984, 81, 3684–3690.
67. C. Anezo, A. H. de Vries, H. D. Holtje, D. P. Tieleman and S. J. Marrink, J. Phys. Chem. B, 2003, 107, 9424–9433.
68. N. Kucerka, S. Tristram-Nagle and J. F. Nagle, J. Membr. Biol., 2005, 208, 193–202.
69. A. Filippov, G. Oradd and G. Lindblom, Langmuir, 2003, 19, 6397–6400.
70. T. Morishita, J. Chem. Phys., 2000, 113, 2976–2982.

---

Footnote

† The relative membrane roughness follows from eqn (2) of ref. 23 as with parameter values given in Table 2. Here, g₁(r) is the membrane displacement correlation function of adjacent membranes in the multilayer.

---