ORIGINAL RESEARCH article
Sec. Lipids, Membranes and Membranous Organelles
Volume 9 - 2022 | https://doi.org/10.3389/fmolb.2022.958537
Lipid/water interface of galactolipid bilayers in different lyotropic liquid-crystalline phases
- 1Faculty of Biochemistry, Biophysics and Biotechnology, Jagiellonian University, Krakow, Poland
- 2Faculty of Mathematics and Computer Science, Jagiellonian University, Krakow, Poland
In this study, carried out using computational methods, the organisation of the lipid/water interface of bilayers composed of galactolipids with both α-linolenoyl acyl chains is analysed and compared in three different lyotropic liquid-crystalline phases. These systems include the monogalactosyldiglyceride (MGDG) and digalactosyldiglyceride (DGDG) bilayers in the lamellar phase, the MGDG double bilayer during stalk phase formation and the inverse hexagonal MGDG phase. For each system, lipid-water and direct and water-mediated lipid-lipid interactions between the lipids of one bilayer leaflet and those of two apposing leaflets at the onset of new phase (stalk) formation, are identified. A network of interactions between DGDG molecules and its topological properties are derived and compared to those for the MGDG bilayer.
Biological membranes (biomembranes) surround each cell and cell organelle. Their fundamental structural element is a lipid matrix that also plays the role of a selective permeability barrier. The biological functions that biomembranes can fulfil depend on the lipid composition of the matrix. The composition can vary within a wide range and determines the types and strength of intermolecular interactions and the molecular dynamics of lipids. Subsequently, it determines the physicochemical, biophysical, mechanical and other properties of the matrix and thus of the biomembrane. As the lipid matrix is an intricate system, experimental and computational studies are carried out on much simpler model membranes. Model membranes are hydrated lipid bilayers of a controlled lipid composition typical of the specific biomembrane. Lipid bilayers have three distinct regions, namely the bulk water region, the polar interface consisting of the lipid heads and water molecules, and the nonpolar bilayer core consisting of the lipid hydrocarbon chains. The interfacial region separates the other two regions and constitutes the first barrier preventing free movement of molecules across the bilayer. Moreover, many important processes occur there (Bondar and Lemieux, 2019). The interface thus plays an essential role in the functioning of the biomembrane.
Even in simple model membranes the lipid/water interface is structurally and dynamically complex. Structurally, because it consists of different types of polar, nonpolar and charged chemical groups and water molecules; dynamically, because the groups are in constant motion and the interfacial water molecules, even though predominantly bound to the lipid head groups (Markiewicz et al., 2015; Calero and Franzese, 2019), undergo rotational and translational motion and exchange with bulk water fast (Rog et al., 2009). Intermolecular interactions, dynamics and spatial organisation of the lipid head groups and water molecules at the interface are strongly interrelated, and this mutual dependence regulates the properties of the interface and thus of the membrane (Disalvo and Disalvo, 2015; Nickels et al., 2015; Frias and Disalvo, 2021).
The lipid composition of the matrix depends on the type of biomembrane within the cell and the function of the cell within the organism. The matrix of the mammalian plasma membrane consists primarily of glycerophospholipids (PL), i.e. phosphatidylcholine (PC), phosphatidylethanolamine (PE), phosphatidylserine (PS), and sphingomyelin (SM), with one saturated and the other mono-cis-unsaturated acyl chains (van Meer et al., 2008). PE together with phosphatidylglycerol are the main lipid representatives of the inner bacterial membranes (Dowhan, 1997). The main constituents of thylakoid membranes of chloroplasts are glycolipids with the galactose moieties and the glycerol backbone as the head group, i.e. monogalactosyldiglyceride (MGDG) and digalactosyldiglyceride (DGDG), and both α-linolenoyl (di-18:3, cis) acyl chains (Dormann and Benning, 2002). Poly-unsaturation of galactolipid acyl chains is indispensable for proper functioning of thylakoid membranes as summarised in Bratek et al., 2019 (Bratek et al., 2019) and citations therein. Lipopolysaccharides and lipid A are the main constituents of the outer membrane of Gram-negative (G–) bacteria (Brandenburg et al., 2016).
Due to the importance of the mammalian plasma membrane and the fact that they are relatively straightforward to handle, single or binary mixed PL bilayers as well as those also containing other natural membrane components have been extensively studied and much is known about their interfaces. The lipid/water interface of PL bilayers has been studied using experimental methods e.g. (Gawrisch et al., 1978; Volkov et al., 2007a; Volkov et al., 2007b; Disalvo et al., 2008; Zhao et al., 2008; Beranova et al., 2012; Cheng et al., 2014; Pokorna et al., 2014), although more detailed information about its properties has been provided by computer modelling, e.g. (Rog et al., 2009; Berkowitz and Vacha, 2012; Disalvo et al., 2014; Nickels et al., 2015; Pasenkiewicz-Gierula et al., 2016; Laage et al., 2017; Elola and Rodriguez, 2018; Martelli et al., 2018; Srivastava and Debnath, 2018; Tian and Chiu, 2018; Calero and Franzese, 2019; Kucerka et al., 2019; Srivastava et al., 2019; Deplazes et al., 2020; Luo et al., 2020; Szczelina et al., 2020; Hande and Chakrabarty, 2022).
In contrast to phospholipids, glycolipids are relatively less frequently studied in spite of their widespread occurrence. Publications on the lipid/water interface of galactolipid bilayers either in the lamellar or non-lamellar phases are rather scarce. Experimental studies were carried out on bilayers consisting of galactolipids with 18:3, 18:2 18:1, 18:0 and 16:0 acyl chains, e.g., (Shipley et al., 1973; Marra, 1986; Mcdaniel, 1988; Webb and Green, 1991; Bottier et al., 2007) whereas computational studies were carried out on di-18:3 MGDG, e.g., (Markiewicz et al., 2015; Baczynski et al., 2018; Szczelina et al., 2020), 80% di-18:3 DGDG and 20% 18:3–16:0 DGDG (Kanduc et al., 2017) and di-16:0 glucolipid and di-16:0 galactolipid (Rog et al., 2007) bilayers. Also publications on the interface of lipopolysaccharide and lipid A bilayers are not numerous, e.g. (Snyder et al., 1999; Wu et al., 2013; Murzyn and Pasenkiewicz-Gierula, 2015; Kim et al., 2016; Luna et al., 2021; Paracini et al., 2022).
Whereas di-18:3-cis DGDG is a bilayer-forming lipid (Deme et al., 2014), di-18:3-cis MGDG is not (Deme et al., 2014); due to the cone shape under ambient conditions it forms an inverse hexagonal (HII) phase in water spontaneously (Sanderson and Williams, 1992). In this study, the organisation of the lipid/water interface of di-18:3-cis MGDG bilayers in different lyotropic phases and of the di-18:3-cis DGDG lamellar bilayer are analysed and compared. In particular, lipid-water as well as direct and water-mediated lipid-lipid interactions are identified. These interactions take place within the same bilayer interface but also between lipids belonging to the interfaces of apposing leaflets when a new phase begins to form. The strength and branching of inter-lipid interactions at the DGDG bilayer interface are analysed using a formal network analysis approach. The analysis demonstrates that the interactions together with the lipid head groups form a dynamic but stable and extended network. The topological properties of the network are determined and compared with those of the MGDG bilayer (Szczelina et al., 2020).
2 Systems and methods
2.1 Simulation systems
In this molecular modelling study, the lipid/water interface of galactolipid systems in three different lyotropic liquid-crystalline phases is investigated. The galactolipids used to build the systems are monogalactosyldiglyceride (MGDG) and digalactosyldiglyceride (DGDG), each with both α-linolenoyl (di-18:3, cis) acyl chains (Figure 1). The investigated phases are MGDG and DGDG lamellar bilayers (Figure 2); the MGDG double bilayer, which forms the stalk phase (Figure 2); and the MGDG inverse hexagonal phase (HII) (Bratek et al., 2019).
FIGURE 1. Molecular structures of di-18:3-cis (A) MGDG and (B) DGDG. The numbering of the acyl chains and glycerol backbone atoms is according to Sundaralingam’s nomenclature (Sundaralingam, 1972), with an exception for the C1″ and C3″ carbon atoms that are swapped here. The numbering of the galactose ring atoms is according to IUPAC convention (McNaught, 1996). The numbers of the carbon and oxygen atoms of the β ring are marked with ′ and the C1, C2 and C3 atoms of the glycerol backbone are marked with ″ to distinguish the atoms of the α and β galactose rings and of the glycerol. The chemical symbol for carbon atoms, C, is omitted and the hydrogen atoms are not shown except for the polar ones shown as empty circles. Oxygen (O) atoms are dark and the carbon atoms are light grey circles, respectively.
FIGURE 2. Final structures of the (A) MGDG and (B) DGDG bilayers after 320 and 1,050 ns of MD simulations, respectively. (C) Initial and (D) final (after ∼1,800 ns of MD simulation) structures of the double MGDG bilayer (W15). (E) The image of W15 in (D) was replicated along the x- and z-axis over periodic boundaries. The atoms are represented in standard colours, except for acyl chain carbon atoms, which are dark blue. The water is shown as a transparent blue surface. The hydrogen atoms are not shown.
The construction and conformational analysis of the computer model of the di-18:3-cis MGDG molecule (Figure 1) are described in Refs (Baczynski et al., 2015; Baczynski et al., 2018) and those of di-18:3-cis DGDG molecule (Figure 1) are described in Supporting Information (SI). MGDG and DGDG with poly-unsaturated acyl chains were chosen because of their widespread occurrence in Nature and the role they play in photosynthetic membranes (Dormann and Benning, 2002). MGDG and DGDG lamellar bilayers were built from scratch using the Packmol package (Martinez et al., 2009). The MGDG bilayer consisted of 450 MGDG molecules (15 × 15 in each leaflet); the DGDG bilayer consisted of 200 DGDG molecules (10 × 10 in each leaflet). The bilayers were hydrated with 30 H2O/lipid, i.e. with 13,500 and 6,000 H2O molecules, respectively, and MD simulated; the MGDG bilayer for 320 ns (Figure 2A), and the DGDG bilayer for 1,050 ns (Figure 2B). 30 H2O/lipid in the MGDG and DGDG bilayers is more than required for their full hydration. The hydration of DGDG in the lamellar phase and MGDG in hexagonal phase is similar (Brentel et al., 1985) and is 16 H2O/DGDG (Crowe et al., 1990) and ∼19 H2O/DGDG (Rand and Parsegian, 1989) and 13–14 H2O/MGDG (Brentel et al., 1985; Crowe et al., 1990; Selstam et al., 1990). This hydration corresponds to about 22% water by weight in both bilayers (Shipley et al., 1973).
The MGDG bilayer was validated in Ref. (Baczynski et al., 2015). Because a pure di-18:3-cis MGDG bilayer does not form in water spontaneously, e.g. (Dormann and Benning, 2002) there are no experimental data for this bilayer. Therefore, the MGDG bilayer was validated indirectly by comparing its structural properties with those of a pure di-16:0 MGDG bilayer (Rog et al., 2007; Lopez et al., 2013) and binary di-18:3 MGDG-DMPC bilayers (Kapla et al., 2012); also, with the help of a well-studied dioleoylPC (DOPC) bilayer, e.g., (Liu and Nagle, 2004; Kucerka et al., 2005; Kucerka et al., 2008). The DGDG bilayer is validated in the Results section.
A stalk is a crucial intermediate in the membrane fusion mechanism (Kozlovsky et al., 2004; Kasson and Pande, 2007). This is a local connection of lipids that belong to the inner leaflets of two bilayers which come into close contact (Figures 2D,E) as a result of their partial dehydration, and involves lipid mixing between these leaflets (Kozlov et al., 1989; Ohta-Iino et al., 2001; Salditt and Aeffner, 2016). The MGDG stalk was generated in MD simulation of the double bilayer system. The double bilayer was constructed by duplicating the MGDG bilayer after 300 ns of MD simulation and placing one bilayer on top of the other (Figure 2C and film Supplementary Video S1). The intra-bilayer water layer contained 6,750 H2O molecules (15 H2O/MGDG) and the outer water layer contained 13,500 H2O molecules (30 H2O/MGDG). The double bilayer (W15 system) was MD simulated for nearly 2,000 ns (Figure 2D). The first vertical connection between the head groups of two galactolipid molecules across the “inner” water layer formed spontaneously (Kozlovsky et al., 2004) within 1 ns of MD simulation (Supplementary Figure S3). A detailed description of the simulation and stalk formation will be presented elsewhere.
The construction of the MGDG HII phase, its MD simulation and validation were described in detail in Ref. (Bratek et al., 2019). Several structures of the phase were tested before its stable computer model was achieved. The stable HII phase consisted of sixteen cylinders, each containing 5,400 water molecules, which corresponded to 30 H2O/MGDG. The MGDG HII phase generated in 3-µs MD simulation was analysed to obtain such basic structural parameters as hexagonal lattice constant, circular and effective radii of the water channel, surface area/MGDG, and order parameter profiles for the MGDG acyl chains (Bratek et al., 2019). To validate the computer model, experimental data for mainly di-18:3 MGDG (Shipley et al., 1973), mainly di-18:2 MGDG (Bottier et al., 2007), dioleoylPE (DOPE) (Rand and Fuller, 1994) and palmitoyloleoylPE (POPE) (Rappolt et al., 2003) HII phases, were used. The comparison with experimental data was performed by extrapolating the linear dependence of the experimental parameters on the HII phase hydration level, to a hydration level of 30 H2O/MGDG in the computer-generated HII phase.
2.2 Simulation parameters and conditions
Force field parameters for the α-linolenic chain, and the glycerol moiety of the galactolipids, were taken directly from the all-atom optimised potentials for liquid simulations (OPLS-AA) force field associated with the software package GROMACS 5.05 (Abraham et al., 2015) and supplemented with the partial charges on the glycerol backbone of galactolipids from Ref. (Maciejewski et al., 2014).
The head group of MGDG comprises a single β-D-galactose and the glycerol backbone to which the galactose is attached by an O-glycosidic bond called here β-1′-1″ linkage (cf. Figure 1). The head group of DGDG comprises two galactose moieties, α-D-galactose and β-D-galactose linked by an O-glycosidic bond called here α-1-6 linkage (cf. Figure 1), also attached to the glycerol backbone by the β-1′-1″ linkage. For the galactose moieties of MGDG and DGDG, OPLS-AA parameters for carbohydrates (Damm et al., 1997) were used. The parameters for MGDG and other mono-glycoglycerolipids have been tested successfully in previous atomistic MD-simulation studies (Rog et al., 2007; Baczynski et al., 2015; Markiewicz et al., 2015; Baczynski et al., 2018; Bratek et al., 2019; Szczelina et al., 2020). The DGDG bilayer is validated in the Results section. For water, the transferable intermolecular potential three-point model (TIP3P) was used (Jorgensen et al., 1983).
MD simulations of the lamellar galactolipid bilayers were carried out in the NPT ensemble, under a pressure of 1 atm and at a temperature of 295 K (22°C) using the software package GROMACS (Abraham et al., 2015). To control the temperature and pressure, for the first 20 ns of MD simulation, the Berendsen thermostat and barostat (Berendsen et al., 1984) were used, and then the Nosé-Hoover (Nose, 1984; Hoover, 1985) and the Parrinello-Rahman (Parrinello and Rahman, 1981) methods were used, respectively. The relaxation time for the temperature was 0.6 ps and for the pressure 1.0 ps. The temperatures of the solute and solvent were controlled independently, and the pressure was controlled anisotropically.
The linear constraint solver (LINCS) algorithm (Hess et al., 1997) was used to preserve the length of any covalent bond with a hydrogen atom, and the time step was set to 2 fs. The van der Waals interactions were cut-off at 1.0 nm. The long-range electrostatic interactions were evaluated using the particle-mesh Ewald summation method with a β spline interpolation order of 5, and a direct sum tolerance of 10–5 (Essmann et al., 1995). For the real space, a cut-off of 1.0 nm, three-dimensional periodic boundary conditions (PBC), and the usual minimum image convention, were used (Essmann et al., 1995). The list of non-bonded pairs was updated every 5 time steps.
The W15 system was MD simulated at 295 K for 320 ns (film SF1, SI). Then, to speed up the process of the MGDG stalk structure formation, the following ∼1.4-µs simulation was carried out at 333 K (60°C) with the time step of 1.5 fs. After that, the temperature was gradually lowered to 295 K and after reaching this temperature, MD simulation was continued for 100 ns with a 2-fs time step. The temperature profile of this simulation is shown in Supplementary Figure S4. All other simulation parameters and conditions as well as the simulation programme were the same as in the case of the lamellar systems.
All trajectories analysed in this paper were recorded every 1 ps.
The MGDG HII phase was generated in a 3-µs MD simulation, also at 295 K (22°C) and the trajectory was recorded every 2 ps (Bratek et al., 2019).
2.3 Network analysis
The methodology used to analyse the interaction network at the bilayer interface is described in detail in Ref. (Szczelina et al., 2020). The basis of network analysis and the main weighted network parameters are summarised here only briefly. Mathematically, a network can be described and modelled by means of graph theory. In the following, the terms “network” and “graph” are interchangeable. Here, the objects of the graph are lipid molecules (centres-of-mass) in one bilayer leaflet (nodes) and intermolecular interactions between them are the graph edges. Consecutive pairs of nodes connected by edges form a path. A cluster is a set of interconnected nodes where each node has a path to all other nodes. A cluster size is the number of nodes which make a particular cluster. A graph is connected when it is made of only one cluster. The node degree is the number of edges connecting this node to other nodes. The node strength is determined by the number of individual interactions that account for each edge of the node, the average energy of each type of interaction and lifetime of the edge. A network bridge is an edge removing which disconnects the graph.
Network analysis was carried out using NetworkX (Hagberg et al., 2008), a Python language software package for creating, manipulating, and studying the structure, dynamics, and functions of complex networks. Network bridges were identified using a bridge-finding algorithm that employs the chain decompositions described in Ref. (Schmidt, 2013). The network at the DGDG bilayer was visualised using Cytoscape (Shannon et al., 2003).
3.1 Systems equilibration and validation
Time profiles of the potential energy (Ep), the average surface area per lipid (AL) and the bilayer width (DRR) of the MGDG and DGDG bilayers are shown in Figure 3. The average AL was obtained by dividing the simulation box surface area by the number of lipids in one bilayer leaflet. The average DRR was defined and calculated as the distance between the average positions of the centres-of-mass of the single galactose rings (MGDG) or of the double galactose rings (DGDG), in the opposite bilayer leaflets, in a similar fashion to Ref. (Baczynski et al., 2015). Each of the three bilayer parameters converged to some constant value (Figure 3). The time profiles for the MGDG bilayer (Figures 3A,C,E) indicated that the bilayer equilibrated within ∼50 ns of MD simulation. However, the equilibration time of the DGDG bilayer was difficult to assess on the basis of the time profiles in Figure 3; therefore, the average values of Ep, AL and DRR were calculated at three time segments (400–500, 500–600, 600–1,050, Supplementary Table S2), which indicated that the DGDG bilayer equilibrated within ∼500 ns of MD simulation. The average values of Ep, AL and DRR for the MGDG and the DGDG bilayer are given in Table 1 and are marked in Figure 3 as straight red lines. Additionally, to compare with some experimental values, the average DCC bilayer width was calculated as the distance between the average positions of the C2″ atoms (cf. Figure 1) in the opposite leaflets of the MGDG and DGDG bilayers (Table 1). For the MGDG bilayer the averages presented below were obtained for the time range 200–300 ns, while for the DGDG bilayer they were obtained for the range 900–1,000 ns of the respective MD simulations. Errors in the average values derived are standard deviation estimates.
FIGURE 3. Equilibration of the lamellar MGDG and DGDG bilayers and the W15 system. Time profiles of the (A) MGDG and (B) DGDG potential energy (Ep); (C) MGDG and (D) DGDG average surface area per lipid (AL); (E) MGDG and (F) DGDG average bilayer width (DRR), during 320 and 1,050 ns, respectively, of MD simulations at 295 K. (A–F) The time (20 ns) when the T and p control methods were switched (cf. sec. 2.2) is marked with a green arrow; the red line shows the average value of a given parameter. Time profiles of Ep (G) for the whole ∼1,800-ns and (H) for the last 100-ns (when the system’s temperature was 295 K) of MD simulation of the W15 system.
The equilibrated lamellar MGDG and DGDG bilayers are shown in Figures 2A,B, respectively.
The values for AL of 61.77 ± 0.50 Å2 and DRR of 41.60 ± 0.28 Å obtained in this study for the MGDG bilayer are very close to those obtained in Ref. (Baczynski et al., 2015), where the MGDG bilayer was validated (cf. sections 2.1).
The values for AL of 65.84 ± 0.63 Å2 and DRR of 42.04 ± 0.33 Å obtained in this study for the DGDG bilayer MD simulated for 1 µs can be compared with those published in the literature—for AL they range between 63 and 78 Å2 and for DRR between 41 and 44 Å. Most of the published values of AL and DRR were obtained from MD simulations, either coarse grained of 64 ± 1 Å2 and of 41 Å, respectively, for the di-16:0 DGDG bilayer (Lopez et al., 2013) and 63 ± 1 Å2 and 44 ± 3 Å, respectively, for the di-18:3 DGDG bilayer (Navarro-Retamal et al., 2018), atomistic united-atom of ∼78 Å2 and 41.7 Å, respectively, for the mixed DGDG bilayer containing 80% di-18:3 DGDG and 20% 18:3–16:0 (Kanduc et al., 2017) or atomistic all-atom of 67 Å2 and 42 Å, respectively, for the di-16:0 DGDG bilayer (Lopez et al., 2013). The experimental values of AL and DRR obtained using X-ray diffraction are ∼75 Å2 and 41.6 Å, respectively (Shipley et al., 1973), and of DRR obtained using neutron diffraction is 41 Å (Mcdaniel, 1988), for the bilayers consisting of DGDG with 18:3, 16:0, 18:4, 18:2, and smaller amounts of 16:1, 18:0, 18:1 acyl chains, but mainly of di-18:3 DGDG (Shipley et al., 1973).
The results of MD simulations show some dependence of the acyl chain unsaturation on the bilayer structural properties, although it should be remembered that the computer models of the DGDG molecule used in those studies had different resolutions, thus the structural parameters derived may somewhat differ from one another. Other differences in the results may stem from the differing acyl chain compositions of the bilayers, e.g. (Shipley et al., 1973; Kanduc et al., 2017). Besides, the surface area/DGDG in the bilayer in Ref. (Kanduc et al., 2017), after initial equilibration, was kept constant, so its value may be not accurate.
The above comparisons demonstrate that the DGDG bilayer generated in this MD simulation study is effective in reproducing the basic bilayer properties determined in previous studies. Furthermore, the entries in Supplementary Table S1 imply that the conformational states of the DGDG head group concur well with experimental, e.g. (Hirotsu and Higuchi, 1976; Wormald et al., 2002; Ziolkowska et al., 2007; Nakae et al., 2018) and computer simulation, e.g. (Peric-Hassler et al., 2010) data for other disaccharides (cf. SI). Thus, the conclusion that the DGDG bilayer is positively validated is justified.
The initial and final structures of the W15 system are shown in Figures 2C,D. W15 is in the process of stalk phase formation (film Supplementary Video S1) and thus is not at equilibrium. Nevertheless, its energy profile (Figure 3G) was calculated for the whole MD simulation time of ∼1,800 nsas well as for the last 100 ns (Figure 3H), when the temperature, after lowering from 333 to 295 K, was stable at 295 K (cf. Methods). The average value of the whole system’s Ep, as well as the values of AL, DRR and DCC for its “flat” part (marked with a black frame in Figure 4A) calculated over the last 60 ns of MD simulation of the W15 systems equilibrated at 295 K, are given in Table 1.
FIGURE 4. Three regions of W15 and top view of the HII phase. (A) The concave region is not marked; the larger connect-15 region formed through the “inner” water layer (15 H2O/MGDG) is marked with a black frame and the smaller connect-30 region newly forming over the PBC through the “outer” water layer (30 H2O/MGDG) is marked with a red frame. The atoms are represented in standard colours, except for the acyl chain carbon atoms which are dark blue. The water is shown as a transparent blue surface. The hydrogen atoms are not shown. (B) Top view of the HII phase generated Ref. (Bratek et al., 2019) after 3 µs of MD simulation. The water-filled cylinders are in different colours to show the HII phase structure better. The water molecules are not shown.
3.2 Lipid-water H-bonds
MGDG has four OH groups that are both donors and acceptors of hydrogen bonds (H-bond) and six O atoms that are only acceptors of H-bonds (Figure 1) can thus make numerous H-bonds with water molecules. The average numbers of particular types of the MGDG-water interactions are given in Table 2.
DGDG has seven OH groups that are both donors and acceptors of H-bonds and eight O atoms that are only acceptors of H-bonds (Figure 1). Accordingly, the average numbers of DGDG-water interactions (Table 2) are greater than those of MGDG, although somewhat smaller than expected.
Water molecules bind preferentially to the MGDG and DGDG rings and are 50% more often H-bond donors than acceptors (Table 2). In the water-glycerol H-bonding, water is the only H-bond donor. The number of water-glycerol H-bonds is only slightly smaller in the DGDG than the MGDG bilayer (Table 2).
The smaller than expected number of H-bonds with water and H-bonded water molecules in the DGDG than the MGDG bilayer is to some extent compensated by the larger number of water bridges (WB) (Pasenkiewicz-Gierula et al., 1997), which link pairs of galactolipid head groups. WBs form predominantly between galactolipid rings, although the number of glycerol-ring WBs in both bilayers is also quite significant (Table 2).
From Figures 2D, 4A it is apparent that the W15 system has two distinct regions—“region of full hydration” and “region of reduced hydration” (Figure 4A). The “full hydration region”, which is called concave is not marked in Figure 4A. The larger “reduced hydration” region forms as a result of the local cross-water connection of two inner leaflets of the double bilayer that were originally separated by the thinner “inner” water layer (15 H2O/MGDG). It is marked with a black frame in Figures 4A and is called connect-15. The smaller “reduced hydration” region is still forming over the PBC as a result of the local cross-water connection of two outer leaflets of the double bilayer that were originally separated by the thicker “outer” water layer (30 H2O/MGDG). It is marked with a red frame in Figures 4A and is called connect-30. The connections are more visible in Figure 2E.
The average number of each type of MGDG-water interaction in a specified region of the W15 system and in the MGDG HII channels is given in Table 3. The numbers of MGDG-water H-bonds and water bridges (horizontal) as well as H2O molecules H-bonded by MGDG in the concave region of W15 are very similar to those in the MGDG lamellar bilayer, but those in the connect regions, particularly in the connect-15, are smaller (Tables 2, 3). The numbers of H-bonds and H-bonded H2O molecules in the MGDG HII phase are smaller than those in the concave region, but are similar to those in the connect regions of W15 (Table 3).
In the connect-15 and connect-30 regions, the horizontal and vertical WBs can be distinguished. The horizontal WBs are between lipids of the same bilayer leaflet (Figure 5A) and the vertical are between lipids of the apposing leaflets (Figure 5B). In both regions, the number of horizontal WBs is larger than that of the vertical ones but somewhat smaller than the number of them in the MGDG lamellar bilayer and the concave region. This could indicate that there is some competition between the horizontal and vertical WBs. Nevertheless, the total number of WBs (horizontal and vertical) in each connect region is larger than the number of those in the concave region and the lamellar bilayer. The number of water bridges (horizontal) in the MGDG HII phase is higher than in the MGDG bilayer and any W15 region (Table 3).
FIGURE 5. Examples of MGDG-MGDG H-bonds and water bridges at the interface of the W15 connect-15 region. (A) Horizontal (between lipids belonging to the same bilayer leaflet) interactions; (B) vertical (between lipids belonging to the apposing leaflets) interactions. The molecules are shown as sticks in standard colours (acyl chains are cut off). The dotted blue lines represent intermolecular interactions. In (A) the water molecule is an acceptor of two H-bonds, in (B) the water molecule is a donor of two H-bonds.
3.3 Lipid-lipid interactions
The MGDG and DGDG heads have both H-bond acceptor and donor groups. Therefore, they can be linked by direct inter-lipid H-bonds at the lipid/water interface, in addition to water bridges, which are water-mediated lipid-lipid interactions (Table 2, 3). The numbers of direct H-bonds in the MGDG and DGDG lamellar bilayers, the W15 system and its specified regions, as well as in the MGDG HII phase, are given in Table 4. The number of H-bond acceptor and donor groups of DGDG is 50% larger than that of MGDG; however, the number of DGDG-water H-bonds is only ∼37% greater, whereas the numbers of WBs and direct H-bonds are 70 and ∼130%, respectively, greater than those of MGDG (Tables 2, 4). This indicates that at the bilayer interface the head groups of DGDG interact preferentially with one another, rather than with water, whereas interactions between the MGDG head groups and water are relatively numerous.
To obtain a better insight, the numbers of ring-water and ring-ring interactions in the DGDG bilayer were calculated for the α and β galactose rings separately. The results given in Tables 2, 4 show that the number of intermolecular interactions of each DGDG ring is different.
In the connect regions of W15, the MGDG head groups form both horizontal and vertical inter-lipid H-bonds (Figure 5), as in the case of WBs. The average numbers of horizontal H-bonds in the three regions of W15 are similar to each other and also similar to the number of them in the MGDG lamellar bilayer. This implies that, in contrast to WBs, different H-bond donor and acceptor groups of MGDG are involved in formation of the horizontal and the vertical direct H-bonds.
The number of lipid-lipid H-bonds in the MGDG HII phase is higher than in the MGDG bilayer or any W15 region (Table 4), as is the number of WBs.
3.4 Orientation of the galactolipid head group
The orientation of the MGDG head group in the bilayer is determined here, as in Ref. (Baczynski et al., 2018), by angle ω between the MGDG head group vector, which connects the C2″ atom in the glycerol and the O4’ atom in the galactose ring (Figure 1A), and the bilayer normal. Ref. (Baczynski et al., 2018) correlations between angle ω and the conformation of the torsion angles of the glycosidic linkages, and also between angle ω and the numbers of head-water and head-head interactions at the MGDG bilayer interface, were calculated. Significant conclusions of those calculations were that there was virtually no correlation between the orientation of the MGDG head group and the conformation of its glycosidic linkage, and that there was only a weak correlation between the MGDG head group orientation and the number of intermolecular interactions of the head.
The most probable (preferred) orientation, called here tilt, of the head group is angle ω, for which the ω distribution has the main maximum. The distributions of ω and tilts for MGDG in the bilayer and in the connect-15 region of W15 obtained in this study are shown in Supplementary Figure S5A, B and are given in Table 1, respectively. In the MGDG bilayer the ω distribution is smooth and the vector tilt is 32°. Both are similar to those in Ref. (Baczynski et al., 2018). In the connect-15 region the ω distribution has a long tail and the vector tilt is 42°. These results can possibly be linked to the somewhat uneven surface of the connect-15 region (Figure 4A) and to the onset of the rotation of some of the MGDG molecules in the process of formation of the stalk structure.
The tilt of the DGDG α ring in the bilayer was obtained from the distribution of angle ω between the α ring vector (C2″-O4 vector, Figure 1B) and the bilayer normal (Supplementary Figure S5C), and that of the DGDG β ring was obtained from the distribution of the ω angle between the β ring vector (C2″-O4’ vector, Figure 1B) and the bilayer normal (Supplementary Figure S5D); both tilts are given in Table 1. The tilt of the DGDG α ring of 30° is almost the same as that of the MGDG β ring of 32°. However, their ω distributions differ. In addition to the main maximum in the ω distribution of the α ring vector at 30°, there are smaller but clear maxima at ∼60°, ∼80° and the last one at ∼140°. These maxima indicate that the α ring may have three additional less populated but stable orientations.
The tilt of the DGDG β ring is 36°. Even though the tilts of the DGDG α and β ring vectors are similar, the rings belong to different planes (Supplementary Figure S6). The distribution of the angle between the planes of the α and β rings shown in Supplementary Figure S7 has two maxima. The higher, relatively narrow maximum is at 82° and the significantly lower one is at 162°. The angles are most likely determined by the preferred populations of the torsion angles of the α-1-6 and β-1′-1″ glycosidic linkages. On the basis of the results of Ref. (Baczynski et al., 2018) it is justified to assume the angles between the α and β ring planes do not depend on the α and β ring tilts.
3.5 Density profile of the terminal CH3 groups of galactolipid acyl chains
The density profiles of the terminal CH3 groups of the poly-cis-unsaturated α-linolenoyl acyl chains of MGDG and DGDG across the bilayer were calculated to estimate the probability of finding the groups in the interfacial region of each bilayer, and to compare this probability with the results of previous experimental, e.g. (Feix et al., 1984; Mihailescu et al., 2011) and computer simulation, e.g. (Mihailescu et al., 2011) studies. The profiles across the MGDG and DGDG bilayers and connect-15 region are shown in Supplementary Figure S8. The probability was calculated for each leaflet of the bilayers from the area under the fragment of the CH3 profile where the electron density of water is non-zero (Supplementary Figure S8). The estimated probability (averaged over both leaflets) of finding a CH3 group at the bilayer interface is ∼15% in the MGDG bilayer and ∼16% in the DGDG bilayer.
3.6 Network analysis of the DGDG bilayer
At the bilayer interface, the galactolipid head groups and interactions (H-bonds and water bridges) between them create a network of interactions. In Ref. (Szczelina et al., 2020) several topological properties of the network in the MGDG bilayers were determined. Here, the same methodology (cf. section 2.3) is used to analyse the interaction network in the DGDG bilayer. The values of the network parameters (cf. section 2.3), averaged over the last 50 ns of MD simulation of the DGDG bilayer, are given in Table 5, together with those obtained for the MGDG bilayers in Ref. (Szczelina et al., 2020). As the interacting groups of MGDG and DGDG head groups are the same, in this calculation the average energies of H-bonding and water bridging of DGDG are assumed to be the same as those of MGDG calculated in Ref. (Szczelina et al., 2020). The interaction network in the DGDG bilayer is presented in Figure 6 and its time changes are shown in film Supplementary Video S2.
FIGURE 6. The networks of lipid interconnections at the interfaces of the DGDG bilayer at two time frames 1 ps apart. (A) The network is connected and (B) the connected network is broken into two clusters. The black rectangle depicts the basic simulation box. Nodes (centres-of-mass of the lipids) are presented as dots in the x,y-plane, and edges as lines connecting respective nodes. The largest cluster is in blue and a single-node cluster is in orange. To avoid problems with edges crossing PBC, 9 copies of each node are presented (strong colour tone for the cluster, soft colour tones for its 8 PBC copies); the edges are drawn in the basic simulation box and all its copies.
The values given in Table 5 are graphically presented in Supplementary Figures S9,S10. These figures also provide additional information. The distribution of the node degrees shown in Supplementary Figure S9A indicates that the most probable degree in the DGDG bilayer is 5. This means that two DGDG head groups are most often connected by five individual inter-lipid interactions, and the smallest number of such interactions is two. The results shown in Supplementary Figure S9C demonstrate that only clusters of sizes 1, 2, 98, 99 and 100 have non-zero probability of forming and the probability of forming a cluster of size 100 is at least two orders of magnitude larger than that of the remaining ones. Time profiles of the average number of clusters and the average sizes of the smallest and the largest clusters are shown in Supplementary Figure S10; the averages are over two networks, each in one bilayer leaflet. Supplementary Figure S10A reveals that in each bilayer leaflet the network is connected for most of the time. However, from time to time one of the connected networks breaks for a short while into two clusters and the average number of clusters is then 1.5. Only in one case does the network break into three clusters (of sizes 1, 1 and 98) and the average number of clusters is then 2 (Supplementary Figure S10A). The average size of the smallest cluster is either 100 or nearly 50 (Supplementary Figure S10B) and that of the largest is either 100 or nearly 100 (Supplementary Figure S10C). For most of the time, the sizes of the smallest and the largest clusters are 100. The time profile of the average node strength plotted in Supplementary Figure S10D only fluctuates around the average value, and this indicates that the average node strength is stable during the simulation time.
4.1 MGDG and DGDG bilayers
An interesting result of this study is that the DGDG bilayer equilibrated after a much longer time than the MGDG bilayer. This effect was previously noticed by Kanduč et al. (Kanduc et al., 2017) and explained as arising from the “pronounced hydrogen-bonding capabilities” of DGDG (Kanduc et al., 2017). Our results indicate that direct DGDG-DGDG H-bonds as well as WBs at the DGDG bilayer interface are indeed numerous (Table 4). The detailed spatial organisation of these inter-lipid links is revealed by network analysis (Hagberg et al., 2008) (cf. Section 3.6). The network of the connections (H-bonds and WBs) is shown in Figure 6 and its dynamics in film Supplementary Video S2. The figure, the film and the network parameters in Table 5 as well as Supplementary Figures S9,S10 demonstrate that the connections are not only numerous, but also extended, and branched and for most of the time encompass all DGDG molecules in each bilayer leaflet.
The interaction network in the DGDG bilayer is qualitatively similar to those in the MGDG bilayers (Szczelina et al., 2020), but is considerably more stable. It also differs in some other aspects. The most probable node degree in the DGDG bilayer is 5 (Supplementary Figure S9A), whereas in the MGDG bilayer it is 3 (Szczelina et al., 2020). The average number of clusters is ∼1 (Table 5; Supplementary Figures S9B,S10A), whereas in the MGDG bilayer, depending on its size, it is ∼1.5 or ∼3 (Table 5). The number of network bridges (0.21) is much smaller than in the MGDG bilayers (Table 5). The small number of network bridges in the DGDG bilayer indicates that the connected network disconnects rarely (Supplementary Figure S10A) and much less often than in the 4 MGDG bilayer, cf. Supporting Information of Ref. (Szczelina et al., 2020).
Probably due to the greater AL (65.8 vs. 61.8 Å2), the lifetimes of network edges consisting of only H-bond interactions and of only water bridge interactions in the DGDG bilayer are about half the length of those in the MGDG bilayers (Table 5). However, because the numbers of individual inter-lipid interactions are about twice as great and the energies of the interactions are the same, the node strength in the DGDG and MGDG bilayers are similar.
In the DGDG bilayer, as in the MGDG bilayers (Szczelina et al., 2020), the average node strength is stable during the analysis time (Supplementary Figure S10D). This indicates that the pattern of lipid interconnections at the bilayer interface is also stable in spite of the relatively short lifetimes of the network edges (Table 5), their fast rearrangements (film Supplementary Video S2) and the overall dynamics of the bilayer interfacial region.
Despite the fact that the network of interconnections at the DGDG bilayer interface is more stable, extended and branched than that at the MGDG bilayer interface, the large head group and cylindrical shape of the DGDG molecule prevent formation of non-bilayer phases, as is the case with the MGDG aggregates.
An apparent disproportion in the number of inter-lipid links at the bilayer interface between the DGDG and the MGDG bilayer (cf. section 3.3) can possibly be explained as follows. The tilt of the MGDG β galactose rings relative to the bilayer normal is 32° (Table 1), thus their polar groups are quite exposed to the water phase. In contrast, the DGDG β rings are screened from water by the α rings so they make fewer H-bonds with water than those of MGDG (Table 2). Due to smaller hydration, the polar groups of the DGDG β ring make fewer ring-ring WBs but more direct ring-ring H-bonds than those of the MGDG β ring and the DGDG α ring (Table 4). As a consequence of the hydration disparity of the DGDG α and β rings, the β-α WBs are more numerous and the β-α H-bonds are less numerous than those of the α-α rings (Tables 2, 4).
The tilt of the DGDG α galactose rings (O4-C2″ vector) is 30° and is practically the same as that of the MGDG β galactose rings, which is 32° (Table 1). However, the average numbers of ring-water H-bonds and ring-ring H-bonds and WBs made by the DGDG α ring are smaller than those made by the MGDG β ring (Tables 2, 4). This is because in addition to the α ring-water and α-α rings interactions, the α ring also interacts with the DGDG β ring. These interactions can be seen in the ω angle distribution in Supplementary Figure S5 as secondary maxima. They indicate that the DGDG α ring, on top of the preferred one, has three other less populated but stable orientations. Two of these orientations (∼60° and ∼80°) are possibly stabilised by its interactions with the β ring (particularly WBs) and one (∼140°) with the glycerol backbone (Tables 2, 4). The values obtained for the DGDG α ring orientation are only in partial agreement with conclusions derived on the basis of experimental data that “the polar head group of this lipid (DGDG) is oriented parallel to the plane of the bilayer” (Marra, 1986; Mcdaniel, 1988). The distribution of angle ω in Supplementary Figure S5 indicates that only a small fraction of the DGDG α rings is oriented parallel to the bilayer plane.
Previous experimental, e.g. (Feix et al., 1984; Mihailescu et al., 2011) and computer simulation, e.g. (Mihailescu et al., 2011) studies have revealed that the terminal CH3 groups of PL acyl chains can locate in the bilayer interfacial region. Using the electron-electron double-resonance methodology Felix et al. (Feix et al., 1984) showed that the probability of finding the CH3 group of a saturated acyl chain at the interface is 14%. While using the neutron diffraction methodology Mihailescu et al. (Mihailescu et al., 2011) showed that the probability of finding the CH3 group of a mono-cis-unsaturated acyl chain at the bilayer interface is 20%. This location of the CH3 groups was also found in an MD simulation study (Mihailescu et al., 2011). In this study the probability of finding the CH3 group of a poly-cis-unsaturated acyl chain was calculated from the electron densities of the CH3 groups and the water of each bilayer leaflet (Supplementary Figure S8), and was found to be ∼15% in the MGDG bilayer and ∼16% in the DGDG bilayer. Thus, these results predict that the terminal CH3 group not only of a saturated and a mono-unsaturated acyl chain of PC but also of a poly-unsaturated acyl chain of MGDG and DGDG can be found in the interfacial bilayer region.
4.2 W15 system
In the W15 system, two MGDG bilayers were initially placed parallel to each other and separated by two water layers, namely, the thinner “inner” water layer containing 6,750 H2O molecules (15 H2O/MGDG) and the thicker “outer” water layer containing 13,500 H2O molecules (30 H2O/MGDG) (Figure 2). In the course of MD simulation, the distance between the bilayers decreased in some places and increased in some others, indicating the onset of stalk structure formation (film Supplementary Video S1). The first vertical lipid-lipid contact across the “inner” water layer formed within the initial 1 ns of MD simulation at 295 K (Supplementary Figure S3). Local partial dehydration leading to formation of the connect-15 region can be seen in film SF1 (SI). In this dehydration process, each MGDG molecule loses on average approximately one H-bond with water and ∼1.5 H-bonded water molecules, but gains ∼0.6 WBs and ∼0.4 inter-lipid H-bonds (Table 3). Thus, interactions with water are replaced by lipid-lipid interaction. Water molecules move from the connect to the concave regions; this process is relatively fast as can be deduced from film Supplementary Video S1.
Hydration of MGDG molecules as well as the number of inter-lipid interactions in the W15 concave region are practically the same as in the MGDG bilayer. This might be because the W15 is in the process of stalk structure formation. The stalk structure involves lipid mixing between apposing leaflets (Kozlov et al., 1989; Ohta-Iino et al., 2001; Salditt and Aeffner, 2016), which requires rotation of lipid molecules so as to transform the concave region into a water-filled tube, whose inner surface consists of lipid heads. This rotation only started in the connect-15 region of W15, so the properties of the concave region are more like those of the MGDG bilayer than of the MGDG HII phase. The onset of rotation of some of the MGDG molecules in W15 can be deduced from the long tail in the ω distribution for the MGDG ring vector in the connect-15 region (Supplementary Figure S5B), and the relatively high value of the MGDG tilt, being 42° (Table 1).
4.3 HII phase
Details of the construction and MD simulation of the MGDG HII phase are described in Ref. (Bratek et al., 2019). The MGDG HII phase consisted of sixteen cylinders (Figure 4B), each filled with 5,400 water molecules (30 H2O/MGDG). In that paper, the basic structural parameters of the HII phase were identified (cf. section 2.1). In this paper, the average hydration and the average number of inter-lipid interactions of the MGDG head groups are calculated. Compared to the MGDG bilayer and the concave region of W15, the MGDG head groups in the HII phase are less hydrated than in the bilayer, but their hydration is similar to that in the connect regions of W15, whereas their mutual interactions are more numerous than in the bilayer and W15.
4.4 Effect of acyl chains on the HII phase structure
Combined X-ray, neutron scattering and MD simulation studies indicate that the length and mono- and poly-unsaturation of PL acyl chains have an impact, among others, on the lipid surface area in PC, e.g. (Pabst et al., 2010; Marquardt et al., 2020) and PE, e.g. (Kucerka et al., 2015) bilayers. Using a different experimental approach, the effect of the length and unsaturation of the acyl chains on lipid hydration in PC and PE monolayers is revealed (Maltseva et al., 2022). Yet, our previous MD simulation study on cis- and trans-mono-unsaturated PC bilayers indicates that the conformation of the double bond does not have much impact on the lipid surface area (Murzyn et al., 2001).
The results for the HII phase are in contrast with those for lamellar PL bilayers. The experimentally derived structural parameters such as hexagonal lattice constant (dhex) and radius of the water channel (r), as a function of the hydration level for mainly di-18:3 MGDG, mainly di-18:2 MGDG, di-18:1 DOPE and 16:0–18:1 POPE HII phases, either in the case of dhex or r, lie on one straight line, irrespectively of the degree of the acyl chains (Ref. (Bratek et al., 2019) and explanations therein). Moreover, the values of the average surface area/MGDG in the cylinders of the HII phases as a function of the hydration level also lie on one curve, irrespectively of the type of acyl chains (Ref. Bratek et al., 2019), SI). Thus, these structural parameters depend predominantly on the hydration level of the phase and not on the degree of the unsaturation of the acyl chains.
The analyses presented in this paper revealed:
(1) In the interfacial region of the MGDG and DGDG bilayers, the galactolipid and water molecules interact via direct H-bonds and water bridges.
(2) At the bilayer/water interface MGDG interacts with water more readily than DGDG.
(3) At the bilayer/water interface the lipid-lipid interactions are more readily formed in the DGDG than the MGDG bilayer.
(4) The disproportionally higher number of DGDG-DGDG interactions relative to the number of the DGDG H-bond donor and acceptor groups can be explained by screening the DGDG β rings from the water by the α rings. This screening results in the hydration disparity of the DGDG α and β rings and the different preferences of the rings to interact via H-bonds and water bridges.
(5) The network of inter-lipid interactions at the DGDG bilayer interface is more stable and extended than that in the MGDG bilayer. Nevertheless, a DGDG aggregate under ambient conditions does not form HII phase in water spontaneously; this is most likely due to the cylindrical shape of the DGDG molecule and its large head group.
(6) In the system consisting of two MGDG bilayers separated by a water layer containing 6,750 H2O molecules (15 H2O/MGDG) a MGDG stalk structure begins to form; the structure is visible as local vertical contacts of MGDG head groups from the apposing bilayer leaflets separated by water-filled tunnels (W15 system).
(7) The number of lipid-lipid and lipid-water interactions at the interface of the water-filled tunnel of the MGDG stalk structure is similar to that of the MGDG bilayer.
(8) The number of lipid-water interactions in the locally connected regions of the MGDG stalk structure is smaller than that at the interface of the MGDG bilayer.
(9) In the locally connected regions of the MGDG stalk structure both horizontal (between lipids from the same bilayer leaflet) and vertical (between lipids from apposing bilayer leaflets) inter-lipid H-bonds and water bridges are formed.
(10) The total number of lipid-lipid horizontal and vertical interactions in the locally connected regions of the MGDG stalk structure is greater than the number of horizontal lipid-lipid interactions in the MGDG bilayer.
(11) The number of lipid-water interactions in the MGDG HII phase is similar to that in the locally connected regions of the stalk structure and smaller than that in the MGDG bilayer.
(12) The number of inter-lipid H-bonds (horizontal) in the MGDG HII phase is greater than in the MGDG bilayer (horizontal) and similar to the total (horizontal and vertical) number of inter-lipid H-bonds in the locally connected regions of the MGDG stalk structure. The number of water bridges in the MGDG HII phase is greater than in the MGDG bilayer and moderately similar to the total (horizontal and vertical) number of water bridges in the locally connected regions of the MGDG stalk structure.
(13) From 11 to 12 one can conclude that when the nonlamellar phase is formed, the lipid-water interactions are, to some extent, replaced by lipid-lipid interactions.
Data availability statement
The original contributions presented in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding author.
MP-G conceived the study. JH, MM and MPG designed the study. JH carried out all MD simulations and analyses. RS performed all network analyses. MM supervised the computations. JH and RS created all figures and films. MP-G wrote the final version of the manuscript with input from all authors. All authors discussed the results, read, and approved the submitted version.
This research was supported in part by PL-Grid Infrastructure. The work has been partly performed under the Project HPC-EUROPA3 (INFRAIA-2016-1–730,897), with the support of the EC Research Innovation Action under the H2020 Programme. The open-access publication of this article was funded by the Department of Computational Biophysics and Bioinformatics, Jagiellonian University, Krakow, Poland.
J.H. gratefully acknowledges the support of Waldemar Kulig, University of Helsinki, and the computer resources and technical support provided by CSC, also the support by the project for PhD students and Young Scientists FBBB N19/MNW/000014.
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.
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fmolb.2022.958537/full#supplementary-material
Abraham, M. J., Murtola, T., Schulz, R., Pall, 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
Baczynski, K., Markiewicz, M., and Pasenkiewicz-Gierula, M. (2018). Is the tilt of the lipid head group correlated with the number of intermolecular interactions at the bilayer interface? FEBS Lett. 592 (9), 1507–1515. doi:10.1002/1873-3468.13048
Beranova, L., Humpolickova, J., Sykora, J., Benda, A., Cwiklik, L., Jurkiewicz, P., et al. (2012). Effect of heavy water on phospholipid membranes: Experimental confirmation of molecular dynamics simulations. Phys. Chem. Chem. Phys. 14 (42), 14516–14522. doi:10.1039/c2cp41275f
Berendsen, H. J. C., Postma, J. P. M., Vangunsteren, W. F., Dinola, A., and Haak, J. R. (1984). Molecular-dynamics with coupling to an external bath. J. Chem. Phys. 81 (8), 3684–3690. doi:10.1063/1.448118
Bottier, C., Gean, J., Artzner, F., Desbat, B., Pezolet, M., Renault, A., et al. (2007). Galactosyl headgroup interactions control the molecular packing of wheat lipids in Langmuir films and in hydrated liquid-crystalline mesophases. Biochim. Biophys. Acta 1768 (6), 1526–1540. doi:10.1016/j.bbamem.2007.02.021
Brandenburg, K., Heinbockel, L., Correa, W., Fukuoka, S., Gutsmann, T., Zahringer, U., et al. (2016). Supramolecular structure of enterobacterial wild-type lipopolysaccharides (LPS), fractions thereof, and their neutralization by Pep19-2.5. J. Struct. Biol. 194 (1), 68–77. doi:10.1016/j.jsb.2016.01.014
Bratek, L., Markiewicz, M., Baczynski, K., Jalocha-Bratek, J., and Pasenkiewicz-Gierula, M. (2019). Inverse hexagonal phase of poly-unsaturated monogalactolipid: A computer model and analysis. J. Mol. Liq. 290, 111189. doi:10.1016/j.molliq.2019.111189
Brentel, I., Selstam, E., and Lindblom, G. (1985). Phase-equilibria of mixtures of plant galactolipids - the formation of a bicontinuous cubic phase. Biochimica Biophysica Acta - Biomembr. 812 (3), 816–826. doi:10.1016/0005-2736(85)90277-9
Damm, W., Frontera, A., TiradoRives, J., and Jorgensen, W. L. (1997). OPLS all-atom force field for carbohydrates. J. Comput. Chem. 18 (16), 1955–1970. doi:10.1002/(sici)1096-987x(199712)18:16<1955::aid-jcc1>3.0.co;2-l
Deme, B., Cataye, C., Block, M. A., Marechal, E., and Jouhet, J. (2014). Contribution of galactoglycerolipids to the 3-dimensional architecture of thylakoids. Faseb J. 28 (8), 3373–3383. doi:10.1096/fj.13-247395
Deplazes, E., Sarrami, F., and Poger, D. (2020). Effect of H3O+ on the structure and dynamics of water at the interface with phospholipid bilayers. J. Phys. Chem. B 124 (8), 1361–1373. doi:10.1021/acs.jpcb.9b10169
Disalvo, E. A., Lairion, F., Martini, F., Tymczyszyn, E., Frias, M., Almaleck, H., et al. (2008). Structural and functional properties of hydration and confined water in membrane interfaces. Biochim. Biophys. Acta 1778 (12), 2655–2670. doi:10.1016/j.bbamem.2008.08.025
Disalvo, E. A., Martini, M. F., Bouchet, A. M., Hollmann, A., and Frias, M. A. (2014). Structural and thermodynamic properties of water-membrane interphases: Significance for peptide/membrane interactions. Adv. Colloid Interface Sci. 211, 17–33. doi:10.1016/j.cis.2014.05.002
Disalvo, E. A. (2015). “Membrane hydration: A hint to a new model for biomembranes,” in Membrane hydration the role of water in the structure and function of biological membranes. Subcellular biochemistry. 71. Editor E. A. Disalvo (Cham: Springer), 1–16.
Feix, J. B., Popp, C. A., Venkataramu, S. D., Beth, A. H., Park, J. H., Hyde, J. S., et al. (1984). An electron-electron double-resonance study of interactions between [14N]- and [15N]stearic acid spin-label pairs: Lateral diffusion and vertical fluctuations in dimyristoylphosphatidylcholine. Biochemistry 23 (10), 2293–2299. doi:10.1021/bi00305a032
Hagberg, A. A., Schult, D. A., and Swart, P. J. (Editors) (2008). “Exploring network structure, dynamics, and function using NetworkX,” Proceedings of the 7th Python in science conference (SciPy2008) (Pasadena, CA USA.
Hande, V. R., and Chakrabarty, S. (2022). How far is “bulk water” from interfaces? Depends on the nature of the surface and what we measure. J. Phys. Chem. B 126, 1125–1135. doi:10.1021/acs.jpcb.1c08603
Hess, B., Bekker, H., Berendsen, H. J. C., and Fraaije, J. G. E. M. (1997). Lincs: A linear constraint solver for molecular simulations. J. Comput. Chem. 18 (12), 1463–1472. doi:10.1002/(sici)1096-987x(199709)18:12<1463::aid-jcc4>3.0.co;2-h
Hirotsu, K., and Higuchi, T. (1976). The conformations of oligosaccharides. III. The crystal and molecular structure of melibiose monohydrate. Bull. Chem. Soc. Jpn. 49 (5), 1240–1245. doi:10.1246/bcsj.49.1240
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 (2), 926–935. doi:10.1063/1.445869
Kanduc, M., Schlaich, A., de Vries, A. H., Jouhet, J., Marechal, E., Deme, B., et al. (2017). Tight cohesion between glycolipid membranes results from balanced water-headgroup interactions. Nat. Commun. 8, 14899. doi:10.1038/ncomms14899
Kapla, J., Stevensson, B., Dahlberg, M., and Maliniak, A. (2012). Molecular dynamics simulations of membranes composed of glycolipids and phospholipids. J. Phys. Chem. B 116 (1), 244–252. doi:10.1021/jp209268p
Kasson, P. M., and Pande, V. S. (2007). Control of membrane fusion mechanism by lipid composition: Predictions from ensemble molecular dynamics. PLoS Comput. Biol. 3 (11), e220. doi:10.1371/journal.pcbi.0030220
Kim, S., Patel, D. S., Park, S., Slusky, J., Klauda, J. B., Widmalm, G., et al. (2016). Bilayer properties of lipid A from various gram-negative bacteria. Biophys. J. 111 (8), 1750–1760. doi:10.1016/j.bpj.2016.09.001
Kozlov, M. M., Leikin, S. L., Chernomordik, L. V., Markin, V. S., and Chizmadzhev, Y. A. (1989). Stalk mechanism of vesicle fusion - intermixing of aqueous contents. Eur. Biophys. J. 17 (3), 121–129. doi:10.1007/BF00254765
Kozlovsky, Y., Efrat, A., Siegel, D. A., and Kozlov, M. M. (2004). Stalk phase formation: Effects of dehydration and saddle splay modulus. Biophys. J. 87 (4), 2508–2521. doi:10.1529/biophysj.103.038075
Kucerka, N., Nagle, J. F., Sachs, J. N., Feller, S. E., Pencer, J., Jackson, A., et al. (2008). Lipid bilayer structure determined by the simultaneous analysis of neutron and x-ray scattering data. Biophys. J. 95 (5), 2356–2367. doi:10.1529/biophysj.108.132662
Kucerka, N., Tristram-Nagle, S., and Nagle, J. F. (2005). Structure of fully hydrated fluid phase lipid bilayers with monounsaturated chains. J. Membr. Biol. 208 (3), 193–202. doi:10.1007/s00232-005-7006-8
Kucerka, N., van Oosten, B., Pan, J. J., Heberle, F. A., Harroun, T. A., Katsaras, J., et al. (2015). Molecular structures of fluid phosphatidylethanolamine bilayers obtained from simulation-to-experiment comparisons and experimental scattering density profiles. J. Phys. Chem. B 119 (5), 1947–1956. doi:10.1021/jp511159q
Liu, Y. F., and Nagle, J. F. (2004). Diffuse scattering provides material parameters and electron density profiles of biomembranes. Phys. Rev. E Stat. Nonlin. Soft Matter Phys. 69 (4), 040901. doi:10.1103/PhysRevE.69.040901
Lopez, C. A., Sovova, Z., van Eerden, F. J., de Vries, A. H., and Marrink, S. J. (2013). Martini force field parameters for glycolipids. J. Chem. Theory Comput. 9 (3), 1694–1708. doi:10.1021/ct3009655
Luo, Y., Zhou, G. B., Li, L., Xiong, S. Y., Yang, Z., Chen, X. S., et al. (2020). Hydrogen bond-induced responses in mid- and far-infrared spectra of interfacial water at phospholipid bilayers. Fluid Phase Equilibr 518, 112626. doi:10.1016/j.fluid.2020.112626
Maciejewski, A., Pasenkiewicz-Gierula, M., Cramariuc, O., Vattulainen, I., and Rog, T. (2014). Refined OPLS all-atom force field for saturated phosphatidylcholine bilayers at full hydration. J. Phys. Chem. B 118 (17), 4571–4581. doi:10.1021/jp5016627
Maltseva, D., Gonella, G., Ruysschaert, J-M., and Bonn, M. (2022). Phospholipid acyl tail affects lipid headgroup orientation and membrane hydration. J. Chem. Phys. 156, 234706–234710. doi:10.1063/5.0092237
Markiewicz, M., Baczynski, K., and Pasenkiewicz-Gierula, M. (2015). Properties of water hydrating the galactolipid and phospholipid bilayers: A molecular dynamics simulation study. Acta Biochim. Pol. 62 (3), 475–481. doi:10.18388/abp.2015_1077
Marquardt, D., Heberle, F. A., Pan, J. J., Cheng, X. L., Pabst, G., Harroun, T. A., et al. (2020). The structures of polyunsaturated lipid bilayers by joint refinement of neutron and X-ray scattering data. Chem. Phys. Lipids 229, 104892. doi:10.1016/j.chemphyslip.2020.104892
Marra, J. (1986). Direct measurements of attractive vanderwaals and adhesion forces between uncharged lipid bilayers in aqueous-solutions. J. Colloid Interface Sci. 109 (1), 11–20. doi:10.1016/0021-9797(86)90276-6
Martelli, F., Ko, H. Y., Borallo, C. C., and Franzese, G. (2018). Structural properties of water confined by phospholipid membranes. Front. Phys. (Beijing). 13 (1), 136801. doi:10.1007/s11467-017-0704-8
Martinez, L., Andrade, R., Birgin, E. G., and Martinez, J. M. (2009). Packmol: A package for building initial configurations for molecular dynamics simulations. J. Comput. Chem. 30 (13), 2157–2164. doi:10.1002/jcc.21224
Mihailescu, M., Vaswani, R. G., Jardon-Valadez, E., Castro-Roman, F., Freites, J. A., Worcester, D. L., et al. (2011). Acyl-chain methyl distributions of liquid-ordered and -disordered membranes. Biophys. J. 100 (6), 1455–1462. doi:10.1016/j.bpj.2011.01.035
Murzyn, K., and Pasenkiewicz-Gierula, M. (2015). Structural properties of the water/membrane interface of a bilayer built of the E. coli lipid A. J. Phys. Chem. B 119 (18), 5846–5856. doi:10.1021/jp5119629
Murzyn, K., Rog, T., Jezierski, G., Takaoka, Y., and Pasenkiewicz-Gierula, M. (2001). Effects of phospholipid unsaturation on the membrane/water interface: A molecular simulation study. Biophys. J. 81 (1), 170–183. doi:10.1016/S0006-3495(01)75689-5
Navarro-Retamal, C., Bremer, A., Ingolfsson, H. I., Alzate-Morales, J., Caballero, J., Thalhammer, A., et al. (2018). Folding and lipid composition determine membrane interaction of the disordered protein COR15A. Biophys. J. 115 (6), 968–980. doi:10.1016/j.bpj.2018.08.014
Ohta-Iino, S., Pasenkiewicz-Gierula, M., Takaoka, Y., Miyagawa, H., Kitamura, K., and Kusumi, A. (2001). Fast lipid disorientation at the onset of membrane fusion revealed by molecular dynamics simulations. Biophys. J. 81 (1), 217–224. doi:10.1016/S0006-3495(01)75693-7
Pabst, G., Kucerka, N., Nieh, M. P., Rheinstadter, M. C., and Katsaras, J. (2010). Applications of neutron and X-ray scattering to the study of biologically relevant model membranes. Chem. Phys. Lipids 163 (6), 460–479. doi:10.1016/j.chemphyslip.2010.03.010
Paracini, N., Schneck, E., Imberty, A., and Micciulla, S. (2022). Lipopolysaccharides at solid and liquid interfaces: Models for biophysical studies of the gram-negative bacterial outer membrane. Adv. Colloid Interfac. 301, 102603. doi:10.1016/j.cis.2022.102603
Pasenkiewicz-Gierula, M., Baczynski, K., Markiewicz, M., and Murzyn, K. (2016). Computer modelling studies of the bilayer/water interface. Biochim. Biophys. Acta 1858 (10), 2305–2321. doi:10.1016/j.bbamem.2016.01.024
Pasenkiewicz-Gierula, M., Takaoka, Y., Miyagawa, H., Kitamura, K., and Kusumi, A. (1997). Hydrogen bonding of water to phosphatidylcholine in the membrane as studied by a molecular dynamics simulation: Location, geometry, and lipid-lipid bridging via hydrogen-bonded water. J. Phys. Chem. A 101 (20), 3677–3691. doi:10.1021/jp962099v
Peric-Hassler, L., Hansen, H. S., Baron, R., and Hunenberger, P. H. (2010). Conformational properties of glucose-based disaccharides investigated using molecular dynamics simulations with local elevation umbrella sampling. Carbohydr. Res. 345 (12), 1781–1801. doi:10.1016/j.carres.2010.05.026
Pokorna, S., Jurkiewicz, P., Vazdar, M., Cwiklik, L., Jungwirth, P., Hof, M., et al. (2014). Does fluoride disrupt hydrogen bond network in cationic lipid bilayer? Time-dependent fluorescence shift of laurdan and molecular dynamics simulations. J. Chem. Phys. 141 (22), 22D516. doi:10.1063/1.4898798
Rand, R. P., and Fuller, N. L. (1994). Structural dimensions and their changes in a reentrant hexagonal-lamellar transition of phospholipids. Biophys. J. 66 (6), 2127–2138. doi:10.1016/S0006-3495(94)81008-2
Rappolt, M., Hickel, A., Bringezu, F., and Lohner, K. (2003). Mechanism of the lamellar/inverse hexagonal phase transition examined by high resolution X-ray diffraction. Biophys. J. 84 (5), 3111–3122. doi:10.1016/S0006-3495(03)70036-8
Rog, T., Murzyn, K., Milhaud, J., Karttunen, M., and Pasenkiewicz-Gierula, M. (2009). Water isotope effect on the phosphatidylcholine bilayer properties: A molecular dynamics simulation study. J. Phys. Chem. B 113 (8), 2378–2387. doi:10.1021/jp8048235
Rog, T., Vattulainen, I., Bunker, A., and Karttunen, M. (2007). Glycolipid membranes through atomistic simulations: Effect of glucose and galactose head groups on lipid bilayer properties. J. Phys. Chem. B 111 (34), 10146–10154. doi:10.1021/jp0730895
Sanderson, P. W., and Williams, W. P. (1992). Low-temperature phase behaviour of the major plant leaf lipid monogalactosyldiacylglycerol. Biochim. Biophys. Acta 1107 (1), 77–85. doi:10.1016/0005-2736(92)90331-f
Selstam, E., Brentel, I., and Lindblom, G. (1990). “The phase-structure of galactolipids and their role in the formation of the prolamellar body,” in Current research in photosynthesis. Editor M. Baltscheffsky (Dordrecht: Springer), 2749–2754.
Shannon, P., Markiel, A., Ozier, O., Baliga, N. S., Wang, J. T., Ramage, D., et al. (2003). Cytoscape: A software environment for integrated models of biomolecular interaction networks. Genome Res. 13 (11), 2498–2504. doi:10.1101/gr.1239303
Shipley, G. G., Green, J. P., and Nichols, B. W. (1973). The phase behavior of monogalactosyl, digalactosyl, and sulphoquinovosyl diglycerides. Biochim. Biophys. Acta 311 (4), 531–544. doi:10.1016/0005-2736(73)90128-4
Snyder, S., Kim, D., and McIntosh, T. J. (1999). Lipopolysaccharide bilayer structure: Effect of chemotype, core mutations, divalent cations, and temperature. Biochemistry 38 (33), 10758–10767. doi:10.1021/bi990867d
Srivastava, A., Malik, S., and Debnath, A. (2019). Heterogeneity in structure and dynamics of water near bilayers using TIP3P and TIP4P/2005 water models. Chem. Phys. 525, 110396. doi:10.1016/j.chemphys.2019.110396
Sundaralingam, M. (1972). Discussion paper: Molecular structures and conformations of the phospholipids and sphingomyelins. Ann. N. Y. Acad. Sci. 195 (20), 324–355. doi:10.1111/j.1749-6632.1972.tb54814.x
Szczelina, R., Baczynski, K., Markiewicz, M., and Pasenkiewicz-Gierula, M. (2020). Network of lipid interconnections at the interfaces of galactolipid and phospholipid bilayers. J. Mol. Liq. 298, 112002. doi:10.1016/j.molliq.2019.112002
Tian, C. A., and Chiu, C. C. (2018). Importance of hydrophilic groups on modulating the structural, mechanical, and interfacial properties of bilayers: A comparative molecular dynamics study of phosphatidylcholine and ion pair amphiphile membranes. Int. J. Mol. Sci. 19 (6), E1552. doi:10.3390/ijms19061552
Wormald, M. R., Petrescu, A. J., Pao, Y. L., Glithero, A., Elliott, T., Dwek, R. A., et al. (2002). Conformational studies of oligosaccharides and glycopeptides: Complementarity of NMR, X-ray crystallography, and molecular modelling. Chem. Rev. 102 (2), 371–386. doi:10.1021/cr990368i
Wu, E. L., Engstrom, O., Jo, S., Stuhlsatz, D., Yeom, M. S., Klauda, J. B., et al. (2013). Molecular dynamics and NMR spectroscopy studies of E. coli lipopolysaccharide structure and dynamics. Biophys. J. 105 (6), 1444–1455. doi:10.1016/j.bpj.2013.08.002
Zhao, W., Moilanen, D. E., Fenn, E. E., and Fayer, M. D. (2008). Water at the surfaces of aligned phospholipid multibilayer model membranes probed with ultrafast vibrational spectroscopy. J. Am. Chem. Soc. 130 (42), 13927–13937. doi:10.1021/ja803252y
Ziolkowska, N. E., Shenoy, S. R., O'Keefe, B. R., McMahon, J. B., Palmer, K. E., Dwek, R. A., et al. (2007). Crystallographic, thermodynamic, and molecular modeling studies of the mode of binding of oligosaccharides to the potent antiviral protein griffithsin. Proteins 67 (3), 661–670. doi:10.1002/prot.21336
Keywords: monogalactolipid, digalactolipid, bilayer, stalk structure, inverse hexagonal phase, inter-lipid interactions, interaction network
Citation: Hryc J, Szczelina R, Markiewicz M and Pasenkiewicz-Gierula M (2022) Lipid/water interface of galactolipid bilayers in different lyotropic liquid-crystalline phases. Front. Mol. Biosci. 9:958537. doi: 10.3389/fmolb.2022.958537
Received: 31 May 2022; Accepted: 11 July 2022;
Published: 15 August 2022.
Edited by:Isabel María López-Lara, National Autonomous University of Mexico, Mexico
Reviewed by:Christopher Garvey, Technical University of Munich, Germany
Ravi K. Shukla, DIT University, India
Copyright © 2022 Hryc, Szczelina, Markiewicz and Pasenkiewicz-Gierula. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.