The pH-dependent change in protonation of ionizable lipids is crucial for the success of lipid-based nanoparticles as mRNA delivery systems. Despite their widespread application in vaccines, the structural changes upon acidification are not well understood. Molecular dynamics simulations support structure prediction but require an a priori knowledge of the lipid packing and protonation degree. The presetting of the protonation degree is a challenging task in the case of ionizable lipids since it depends on pH and on the local lipid environment and often lacks experimental validation. Here, we introduce a methodology of combining all-atom molecular dynamics simulations with experimental total-reflection x-ray fluorescence and scattering measurements for the ionizable lipid Dlin-MC3-DMA (MC3) in POPC monolayers. This joint approach allows us to simultaneously determine the lipid packing and the protonation degree of MC3. The consistent parameterization is expected to be useful for further predictive modeling of the action of MC3-based lipid nanoparticles.

The pH-dependent change in protonation of ionizable lipids is fundamental for endosomal release and hence the success of lipid nano-particles (LNPs) as mRNA delivery systems.1 Ionizable lipids are already successfully used in FDA-approved drugs for the treatment of amyloidosis and in Moderna’s vaccines against SARS-CoV-2.2,3 Hence, the design of ionizable lipids with high transfection efficiencies is expected to play a key role in modern RNA therapeutics. Despite their importance, the microscopic structural changes upon acidification, preceding the endosomal escape of RNA, have not been resolved so far. Molecular dynamics (MD) simulations are particularly suited to provide atomistic insights, however, they require an a priori knowledge of the lipid packing and the protonation degree. The presetting of the protonation degree is not always a simple task since it depends sensitively on the environment. Three competing interactions have to be considered: (i) conformational rearrangements of the lipids, (ii) partial dehydration of the ionizable group, and (iii) interactions with polar and charged lipids or ions. Due to these contributions, it is not surprising that ionizable lipids can have significant pKa shifts. For example, the pKa value of Dlin-MC3-DMA (MC3) was found to be shifted from about 9.4 at infinite dilution4 to 6.4 inside a lipid nano-particle.5 

Several theoretical models are available for pKa predictions. Popular approaches are based on solving the linearized Poisson–Boltzmann equation6 or generalized Born models.7 Such simple models work well for simple systems but fail in more complicated cases where conformational rearrangements become important. MD simulations include these conformational rearrangements and, when combined with free energy methods such as free energy perturbation8 or thermodynamic integration,9 yield the free energy of deprotonation and hence the pKa. Using implicit or explicit solvent, free energy simulations have been successfully applied to calculate pKa shifts in proteins.10,11 In addition, constant-pH MD simulations have been developed to overcome limitations of fixed protonation states in conventional MD simulations.12,13 However, these methods are computationally expensive and their accuracy depends on the underlying force fields and the convergence of the sampling of all degrees of freedom. In particular in lipid systems, which show pH dependent structural phase transitions,14 such convergence is difficult to achieve. In addition, the correct packing of the lipids is crucial. To correctly reproduce the lateral structure, the force fields of all components including ionizable lipids, uncharged lipids, and the water model have to be accurate. Shortcomings of the force fields, for instance in reproducing the experimental surface tension of water, will lead to deviations in the packing such that simulations can fail to reproduce experimental surface pressure-area isotherms in simple monolayer systems.15–18 Consequently, consistent methods are required to correctly assign the protonation degree and the packing of lipids to yield robust and meaningful results. In this respect the combination of all-atom MD simulations and x-ray scattering techniques offers promising perspectives, as the latter are uniquely suited to determine the packing of lipid layers19–21 and the protonation degree22 by quantifying the interfacial counterion excess.23,24

The aim of this current work is to introduce a methodology of combining total reflection x-ray scattering and x-ray fluorescence experiments to yield a reliable prediction as benchmark for consistent modeling. We focus on a lipid monolayer composed of the phospholipid POPC and the transfection lipid MC3, one of the most widely used ionizable lipids due to its high transfection efficiency.25 These simple monolayers, bound to the air/water interface at controlled packing density, are an ideal starting point since they allow us to quantitatively compare the results from all-atom MD simulations and scattering experiments. The joint approach allows us to optimize the protonation degree and the lipid packing. In turn, the simulations yield a detailed insight into the layer structure at the molecular level.

The transfection lipid DLin-MC3-DMA (MC3, see Fig. 1) was purchased from MedChemExpress (Monmouth Junction, NJ, USA). The phospholipid 1-palmitoyl-2-oleoyl-glycero-3-phosphocholine (POPC, see Fig. 1) was purchased from Sigma-Aldrich (Merck KGaA, Germany). Potassium bromide (KBr), Hydrobromic acid (HBr) 48 wt. % in H2O and Tris(hydroxymethyl)-aminomethan (Tris-Base, Cl-free) were of analytical grade, purchased also from Sigma-Aldrich (Merck KGaA, Germany) and used without further purification. Chloroform (purity ≥ 99%, anhydrous) and methanol (purity ≥ 99.9%) were purchased from Sigma-Aldrich (Merck KGaA, Germany) and used as received.

FIG. 1.

Molecular structures of the phospholipid POPC and the ionizable lipid MC3.

FIG. 1.

Molecular structures of the phospholipid POPC and the ionizable lipid MC3.

Close modal

The mixture of MC3-POPC (20:80 mol. %) was prepared by mixing the two dry components and subsequently dissolving them at a concentration of 1 mg/ml in chloroform-methanol (7:3 v/v). This solution was then spread onto an air-aqueous interface. The resulting monolayer film was laterally compressed until reaching a final lateral pressure of π = 30 mN/m while in some cases recording isotherms (see Fig. S1 in the supplementary material). Two solutions were used as subphases for different experiments: (KBr 2 mM, Tris-Base 5 mM, pH 7.5) and (KBr 2 mM, Tris-Base 5 mM, pH 5.0). The corresponding pHs were reached by adding few drops of a 1M HBr solution. The resulting Br concentration after titration was 5.8 mM for pH 7.5 and 6.8 mM for pH 5.0 (see supplementary material for further details). We note that Tris-Base does not buffer very well anymore at pH 5.0, but sufficiently well to stabilize the pH over the time of the measurement (≈1 h). We worked at comparatively low Br bulk concentrations in order to minimize undesired x-ray fluorescence contributions from the illuminated aqueous bulk region close to the surface.

Grazing-incidence x-ray scattering experiments (GIXOS and TRXF, see below) were carried out at the beamline P08 at storage ring PETRA III of Deutsches Elektronen-Synchrotron (DESY, Hamburg, Germany), with the settings described in Refs. 26 and 27, from which the following text is partially reproduced. The Langmuir trough (Riegler and Kirstein, Potsdam, Germany) was located in a hermetically sealed container with Kapton windows, and the temperature was kept at 20 °C by a thermostat. The container was constantly flushed with a stream of humidified helium (He) to prevent air scattering and the generation of reactive oxygen species. The synchrotron x-ray beam was monochromatized to a photon energy of 15 keV, corresponding to a wavelength of λ = 0.827 Å. The incident angle was adjusted to θi = 0.07°, slightly below the critical angle of total reflection, θc = 0.082°. A ground glass plate was placed less than 1 mm beneath the illuminated area of the monolayer in order to reduce mechanically excited surface waves. The beam footprint on the water surface was 1 × 60 mm2 as imposed by the incident beam optics.

Under total-reflection conditions, an x-ray standing wave (SW) is formed at the air/water interface. The penetration depth Λ of its evanescent tail into the aqueous hemispace is a function of the angle of incidence θi28 and for the incident angle used is ≃8 nm. The exact shape ϕ(z) of the SW intensity along the vertical position z (for a given incident angle) follows from the interfacial electron density profile ρe(z), and can be calculated as described previously.21,29

1. Grazing-incidence x-ray off-specular scattering (GIXOS)

Analogous to conventional x-ray reflectometry, GIXOS allows reconstructing the interfacial electron density profile (i.e., the laterally-averaged structure of the lipid layer in the direction perpendicular to the surface) from the qz-dependent scattering intensity, however at fixed incident angle. The details of this technique are described in Refs. 20, 30, and 31. As explained more briefly in Kanduč et al.,32 the qz-dependence of the diffuse scattering intensity I(qxy ≠ 0, qz) recorded at low-enough yet finite qxy (“out of the specular plane”) with the help of a narrow slit contains information equivalent to that of the conventional reflectivity R(qz) and can be transformed as:
I(qxy0,qz)=I(qz)=R(qz)V(qz)RF(qz)
(1)
to good approximation, where RF(qz) is the reflectivity of an ideal interface between the two bulk media and V(qz) is the Vineyard function.33 The approximation is based on the assumption of conformal topographic roughness of all surfaces, which is justified for a lipid monolayer subject to capillary wave roughness. In the present work, the GIXOS signal was measured at qxy = 0.04 Å−1.

Like conventional reflectivity curves, GIXOS curves can be computed on the basis of an assumed interfacial electron density profile ρe(z), via the phase-correct summation of all reflected and transmitted partial waves occurring at the density gradients.21,29 Here, we either modeled ρe(z) as a set of rough slabs representing the tail and headgroup regions of the monolayer with adjustable thickness, electron density, and roughness parameters27,31,32 or extracted ρe(z) directly from the MD simulations. In all cases, ρe(z) was then discretized and the corresponding qz-dependent reflectivities R(qz) were calculated from the Fresnel reflection laws at each slab-slab interface using the iterative recipe of Parratt34 or the equivalent Abebe’s matrix method35 and subsequently multiplied with V(qz)/RF(qz) to obtain the theoretical GIXOS signal. For fits to the experimental data, a constant scale factor and a constant intensity background was applied to the modeled intensities, and the fitting range was limited to the consensus range of validity (qz>3qzC, where qzC is the qz-value at the critical angle of total reflection).31 

2. Total-reflection x-ray fluorescence (TRXF)

The fluorescence signal induced via photoelectric ionization by the x-ray beam under total reflection conditions was recorded with an Amptek X-123SDD detector (Amptek, Bedford, USA). The detector was placed almost parallel to the water surface and perpendicular to the x-ray beam axis, in order to keep elastic and Compton scattering into the detector as low as possible. The center of the detector view angle was set to coincide with the beam footprint position on the water surface.

The fluorescence intensity IF emitted by the Br counterions present at the interface is determined by their interfacial depth profile c(z).36 On a quantitative level, IF is proportional to the spatial integral over the product of c(z) and the known SW intensity profile ϕ(z) introduced above,
IF=Ac(z)ϕ(z)dz
(2)
where the prefactor A can be calibrated with a suitable reference. Here, we used the surface of a bare KBr solution as reference, for which c(z) is known and can be considered constant.21,36 Experimentally, IF was obtained by fitting the intensity peak associated with the respective Kα or Lβ emission lines of Br (at ≈12.0 keV and ≈1.5 keV) in the recorded fluorescence spectra with a Gaussian function. Constraints on the peak positions were imposed based on the tabulated line energies.37 The standing wave ϕ(z) was calculated with a slab-model representation of the interfacial electron density profile,21 whose parameters were previously obtained in fits to the GIXOS curves. Note that roughness can be neglected for the computation of ϕ(z) when qz is low, as is the case under total reflection.36 

1. Infinite dilution limit

The ionizable MC3 lipids are protonated according to the reaction
MC30+H3O+MC3++H2O
(3)
At infinite dilution of the lipids, the protonation degree is given by the Henderson–Hasselbalch equation
pH=log1αα+pKa.
(4)
for the protonation degree α, where we assume pKa = 9.47 at infinite dilution.4 

2. Self-consistent equation for the protonation degree

The pKa is related to the equilibrium constant Ka by
pKa=logKa
(5)
At the monolayer/water interface, the protonation of MC3 lipids are influenced by the distribution of H3O+ ions and the equilibrium constant Ka is given by
Ka=1αα[H3O+]surface
(6)
where [H3O+]surface is the concentration at the position of the ionizable group in the monolayer. The surface concentration is related to the bulk concentration via a Boltzmann distribution38,
[H3O+]surface=[H3O+]bulkeeΦ0/(kBT)
(7)
where e is the elementary charge, kBT the product of Boltzmann constant and temperature, and Φ0 the electrostatic potential at the position of the ionizable group. Combining Eqs. (6) and (7) and using the definitions for pH and pKa yields a self-consistent equation for the pH as a function of the protonation degree39 
pH=log1αα+pKaeΦ0kBTln10.
(8)
The last term is the electrostatic contribution. It takes the dependence of the surface potential on the H3O+ distribution into account. Without this term the ordinary equation for the acid deprotonation equilibrium in bulk solution is obtained [Eq. (4)].
The electrostatic potential Φ(z) is calculated from the Poisson Boltzmann equation at a planar interface
ε0εrd2Φ(z)dz2=i=±qic0eqiΦ(z)/kBT
(9)
Here, z is the distance perpendicular to the surface, q± = ±e is the ion charge in a 1:1 electrolyte, c0 is the bulk salt concentration, ɛ0 is the dielectric constant of vacuum, ɛr is the relative dielectric constant of water. For the constant charge boundary condition, we obtain the analytical solution
Φ(z)=2kBTeln1+Γ0eκz1Γ0eκz
(10)
where Γ0 is related to the surface potential via
Γ0=tanheΦ04kBT
(11)
and κ−1 is the Debye screening length
κ1=ε0εrkBTe2c0
(12)
The surface potential Φ0 = Φ(0) is obtained from
Φ0=2kBTearcsinhσ8kBTc0ε0εr
(13)
The surface charge density σ is related to the area per lipid via σ = αefMC3/Alip, where fMC3 is the fraction of MC3 among all lipids. The experimentally obtained ion concentrations of c0 = 6.8 mM for α = 67.5% and c0 = 5.8 mM for α = 15% were used for the theoretical predictions. In the self-consistent approach, the protonation degree depends on the surface potential Φ0. Clearly, Φ0 depends on σ and hence on the protonation degree α and the area per lipid Alip. In our approach, Φ0 and σ were obtained from the results of the final simulations with optimized α and Alip based on the experiments. σ was obtained by dividing the total number of protonated MC3 lipids in one leaflet of the monolayer by the area of the leaflet. The values of σ calculated at α = 67.5% and α = 15% were 0.029 and 0.006 Cm−2, respectively. With these values, we obtain surface potentials of Φ0 = 92.5 mV and Φ0 = 33.5 mV at α = 67.5% and α = 15%, respectively.
We simulated the MC3-POPC monolayers (20:80 mol. %) using the Gromacs simulation package (v-2021.5).40 The monolayers were created by translating the leaflets of an equilibrated bilayer which in turn was created using the Mem-Gen webserver.41 The AMBER Lipid 17 force field42 was used to describe the POPC lipids. For the cationic and neutral MC3 molecules we used our recently developed force fields.43 These parameters have the advantage that they closely reproduce the structure of MC3 in lipid layers as judged by neutron reflectometry experiments. In addition, the parameters are compatible with the AMBER force field family. Ions were described using the Mamatkulov–Schwierz force field parameters44 and the TIP3P water model was used.45 The systems contained 160 POPC and 40 MC3 lipids per monolayer, and 60 water molecules per lipid. Chloride ions were added to neutralize the positive charge of the protonated MC3 lipid molecules. The systems were energy minimized using a gradient descent algorithm. The temperature of the systems was maintained at 293.15 K using the velocity rescaling thermostat with stochastic term46 and a time constant of 1.0 ps. For the NPT simulations, semi-isotropic Parrinello–Rahman47 was used. The Lennard-Jones interaction potential was cut off and shifted to zero at 1.2 nm. Electrostatic interactions were cut-off at 1.2 nm and the particle mesh Ewald (PME) method48 with a Fourier grid spacing of 0.12 nm was used to evaluate the long-range electrostatics part. All bonds involving hydrogens were constrained using the LINCS algorithm.49 A timestep of 2.0 fs was used to integrate the equations of motion. Analysis was performed using Gromacs inbuilt routines and with MDAnalysis50 and trajectories were visualized with visual molecular dynamics (VMD).51 To obtain the desired area per lipid, the in-plane components Pxx and Pyy of the pressure tensor were varied in 11 steps from −50 to −34 bar. These simulations had durations of 100 ns. A semi-isotropic barostat was employed, such that Pxx = Pyy. The box size in z-direction, Lz = 20 nm, was kept constant. The surface tension γ can be deduced from the pressure tensor components52 as
γ=LznPzzPxx+Pyy2
(14)
where n = 2 is the number of interfaces in the simulation box.53 The value of γ was directly obtained using the Gromacs gmx energy routine. The covered range of in-plane pressure tensor components corresponds to a surface tension range of 34 mN/m ≤ γ ≤ 50 mN/m. The surface pressure is calculated from π = γ0γ, where γ0 = 58 mN/m is the surface tension of TIP3P water at 293.15 K as obtained through extrapolation of the data by Vega and de Miguel52 The last 50 ns of each of these simulations at given value of surface pressure were used to obtain the electron density, from which the GIXOS curves were evaluated using Eq. (1). The Refnx python package54 was used to obtain the reflectivity profiles [R(qz) and RF(qz) in Eq. (1)]. The simulations which best reproduced the experimental GIXOS curves at pH 5.0 and pH 7.5 for α = 1.0 had the parameters summarized in Table I. The simulations with experimentally validated values of Alip and α were performed for 500 ns in the NVT ensemble and had the parameters summarized in Table II.
TABLE I.

Simulation parameters that best reproduce the experimental GIXOS curves for α = 1.0.

pHFixed parametersCalculated parameters
αLz (nm)Pxx/Pyy (bar)Pzz (bar)Lx/Ly (nm)Alip (nm2)γ (mN/m)π (mN/m)
5.0 1.0 20.0 46.0 −0.36 ± 0.72 12.10 ± 0.02 0.73 ± 0.002 45.9 ± 0.9 12.1 ± 0.9 
7.5 1.0 20.0 48.0 1.26 ± 1.0 12.25 ± 0.02 0.75 ± 0.003 49.8 ± 1.2 8.2 ± 1.2 
pHFixed parametersCalculated parameters
αLz (nm)Pxx/Pyy (bar)Pzz (bar)Lx/Ly (nm)Alip (nm2)γ (mN/m)π (mN/m)
5.0 1.0 20.0 46.0 −0.36 ± 0.72 12.10 ± 0.02 0.73 ± 0.002 45.9 ± 0.9 12.1 ± 0.9 
7.5 1.0 20.0 48.0 1.26 ± 1.0 12.25 ± 0.02 0.75 ± 0.003 49.8 ± 1.2 8.2 ± 1.2 
TABLE II.

Parameters used for the simulations with experimentally validated area per lipid and protonation degree.

pHFixed parametersCalculated parameters
αLx/Ly (nm)Lz (nm)Alip (nm2)Pxx (bar)Pyy (bar)Pzz (bar)γ (mN/m)π (mN/m)
5.0 0.675 12.10 20.0 0.73 −47.8 ± 0.35 −48.2 ± 0.4 −0.36 ± 0.46 47.6 ± 0.6 10.4 ± 0.6 
7.5 0.150 12.25 20.0 0.75 −48.8 ± 0.6 −49.2 ± 0.8 −0.42 ± 0.56 49.1 ± 0.8 8.9 ± 0.8 
pHFixed parametersCalculated parameters
αLx/Ly (nm)Lz (nm)Alip (nm2)Pxx (bar)Pyy (bar)Pzz (bar)γ (mN/m)π (mN/m)
5.0 0.675 12.10 20.0 0.73 −47.8 ± 0.35 −48.2 ± 0.4 −0.36 ± 0.46 47.6 ± 0.6 10.4 ± 0.6 
7.5 0.150 12.25 20.0 0.75 −48.8 ± 0.6 −49.2 ± 0.8 −0.42 ± 0.56 49.1 ± 0.8 8.9 ± 0.8 

The lipid monolayers contain 80 mol. % of the charge-neutral zwitterionic phospholipid POPC and 20 mol. % of the ionizable lipid MC3 (Fig. 1). Dependent on the pH, the ionizable MC3 can be neutral (MC30) or positively charged (MC3+) such that the monolayer contains up to three lipid species.

Our aim is to elucidate the influence of pH on the properties of the monolayers at two different pH values, pH 7.5 and pH 5.0, representing the extracellular and the endosomal milieu, respectively. However, to perform simulations two important parameters have to be chosen. The first one is the lateral monolayer packing in terms of the average area per lipid, Alip. It should match the experimental one at a lateral pressure of π = 30 mN/m, a value representative of the packing in lipid bilayers.55 The second one is the protonation degree α of MC3, which should be consistent with the experimental value at given pH conditions. To correctly assign Alip and α, we obtain complementary information from two synchrotron-based grazing-incidence x-ray scattering techniques, namely grazing-incidence x-ray off-specular scattering (GIXOS) and total-reflection x-ray fluorescence (TRXF).36,56,57 While GIXOS, like conventional reflectometry, precisely reveals the monolayer thickness27 (albeit with limited insights into the conformations and lateral distributions of the different molecular species), TRXF is highly sensitive to the interfacial accumulation of counterions at charged monolayers23,24,36 and can therefore serve to determine protonation degrees.22 As will be shown later, the accuracy in the experimental determination of α is further enhanced by including the counterion distribution from the simulations.

Figure 2 outlines the methodology to obtain consistent values for Alip and α. Initially, α is estimated based on the Henderson–Hasselbalch equation which predicts α = 1 for both pH values (Fig. 3). Note that the self-consistent equation yields an improved estimate but requires the surface charge density, which is not known a priori. In a second step, we performed simulations with the assumed α value while systematically varying Alip to identify the values that match the layer thickness obtained in the GIXOS measurements at the two pH conditions. In the next step, the pH-dependent protonation degrees α are deduced from the TRXF measurements while imposing the known values of Alip as well as the shape of the counterion distribution found in the simulations. Finally, new simulations with the experimentally validated values of Alip and α are performed for both pH conditions. In the following, we discuss the results from each step in detail.

FIG. 2.

Methodology of combining all-atom molecular dynamics simulations with experimental total reflection x-ray scattering and x-ray fluorescence for consistent modeling of the area per lipid and the protonation degree.

FIG. 2.

Methodology of combining all-atom molecular dynamics simulations with experimental total reflection x-ray scattering and x-ray fluorescence for consistent modeling of the area per lipid and the protonation degree.

Close modal
FIG. 3.

Dependence of the protonation degree α on pH obtained from the Henderson–Hasselbalch Eq. (4) and the self-consistent Eq. (8) using pKa = 9.4.4 Dashed straight lines indicate the predictions for the two experimental pH values, 5.0 and 7.5.

FIG. 3.

Dependence of the protonation degree α on pH obtained from the Henderson–Hasselbalch Eq. (4) and the self-consistent Eq. (8) using pKa = 9.4.4 Dashed straight lines indicate the predictions for the two experimental pH values, 5.0 and 7.5.

Close modal

Figure 4(a) shows the GIXOS curves obtained experimentally with monolayers at π = 30 mN/m on aqueous subphases with pH 7.5 (top) and with pH 5.0 (bottom). The subphases additionally contained defined concentrations of Br and K+ ions, as required for the simultaneous TRXF measurements discussed further below. The curves contain information on the electron density profiles of the monolayers, where the overall layer thicknesses are mainly encoded in the qz-positions of the intensity minima. The lines superimposed to the experimental data points are the simulated GIXOS curves computed from the electron density profiles obtained in simulations. Based on the Henderson–Hasselbalch equation [Eq. (4)] the protonation degree is α ≈ 1 (Fig. 3). However, significant deviations from the experimental curve [dashed line in Fig. 4(a)] are obtained with simulations at the experimental surface pressure corresponding to a surface tension of γ = γ0π = 28 mN/m, with γ0 = 58 mN/m for TIP3P water (see Sec. II). Clearly, the packing of the lipids and hence Alip is not reproduced. This deviation is likely caused by the TIP3P water model which, like most water models, does not yield the correct surface tension of water.18,52 In order to provide improvement, we systematically vary Alip in the simulations to obtain the best-matching value for each pH, characterized by the minimum χ2 deviation (Figs. S4 and S5) from the experimental GIXOS data. The simulated GIXOS curves obtained after optimization (solid lines) reproduce the experimental data well, demonstrating that the simulations properly capture the overall structure of the monolayers.

FIG. 4.

(a) Comparison of experimental GIXOS curves with the simulated ones for the optimized values of Alip. The system consists of 20% MC3+ and 80% POPC. The curves at different pH values are vertically shifted for clarity. The dashed black line is obtained by using the same surface pressure as in the experiments (π = 30 mN/m) and corresponds to Alip = 0.63 nm2. The solid blue lines correspond to the simulations with optimized area per lipid (Alip = 0.75 nm2 for pH = 7.5 and Alip = 0.73 nm2 for pH = 5.0). (b) Snapshot from a simulation with Alip = 0.63 nm2 and α = 1 after 100 ns of equilibration. All other run parameters are summarized in Table I. MC3+ is represented in yellow with the head group nitrogen highlighted by yellow spheres. POPC is shown in white and chloride ions in cyan color. (c) Corresponding electron density profile ρe(z). The plot also shows the counterion density profile c(z) obtained in the simulations (data points). The solid line superimposed is an analytical fit with Eq. (15).

FIG. 4.

(a) Comparison of experimental GIXOS curves with the simulated ones for the optimized values of Alip. The system consists of 20% MC3+ and 80% POPC. The curves at different pH values are vertically shifted for clarity. The dashed black line is obtained by using the same surface pressure as in the experiments (π = 30 mN/m) and corresponds to Alip = 0.63 nm2. The solid blue lines correspond to the simulations with optimized area per lipid (Alip = 0.75 nm2 for pH = 7.5 and Alip = 0.73 nm2 for pH = 5.0). (b) Snapshot from a simulation with Alip = 0.63 nm2 and α = 1 after 100 ns of equilibration. All other run parameters are summarized in Table I. MC3+ is represented in yellow with the head group nitrogen highlighted by yellow spheres. POPC is shown in white and chloride ions in cyan color. (c) Corresponding electron density profile ρe(z). The plot also shows the counterion density profile c(z) obtained in the simulations (data points). The solid line superimposed is an analytical fit with Eq. (15).

Close modal

For pH 7.5 and pH 5.0, the best-matching values are Alip = 0.75 nm2 and Alip = 0.73 nm2, respectively, which is in the typical range for lipid bilayers the fluid phase.58 The associated surface pressures in the simulations are as low as π ≈ 12 mN/m (pH 5.0, see Table I) and π ≈ 8 mN/m (pH 7.5, see Table I) and deviate by a factor of ≈3 from the experimental value. Interestingly, the agreement is much better on the level of the surface tension, where the simulations yield γ ≈ 46 mN/m (pH 5.0, see Table I) and γ ≈ 50 mN/m (pH 7.5, see Table I), compared to the experimental value of γ = γ0π = 42 mN/m, where γ0 = 72 mN/m. As seen in Tables I and II, the protonation degree of MC3 only has a minor influence on γ and π for a fixed area per lipid. In conclusion, the lateral density in terms of Alip appears to be the most robust parameter for the experimental validation of MD simulations of lipid monolayers.

The electron density profile ρe(z) and the counterion density profile c(z) obtained in the simulations for Alip = 0.73 nm2 are shown in Fig. 4(c). The profiles are consistent with earlier reports on fluid lipid monolayers at comparable lateral pressures,19 however the electron density maximum representing the lipid headgroups is slightly less pronounced than for commonly studied phospholipid layers, because MC3 has a more compact headgroup without the electron-rich P atom.

The pH-dependent protonation degree of MC3 was determined experimentally through measurements of the interfacial excess of Br counterions, Γ, which coincides with the density of positive net charges in the monolayer.36 Br rather than the more common Cl was used as counterion in these experiments because the used TRXF setup was slightly more sensitive to heavier elements. Note, however, that Br is equally suited as indicator of the protonation degree. Figure 5(a) shows the x-ray fluorescence spectra in a narrow energy range around the peak associated with the Kα emission line of Br. The spectra were obtained under total reflection conditions with the bare aqueous subphases and with lipid monolayers at π = 30 mN/m both at pH 5.0 and pH 7.5. The complete spectra are shown in the supplementary material (Fig. S2). The intensities of the emission lines are much higher in the presence of the monolayer at both pH values, which reflects the accumulation of Br ions at the interface required to compensate the positive charges of MC3+. The Lβ emission lines behave consistently (see Fig. S2 in the supplementary material). The intensity is much higher at pH 5.0 than at pH 7.5, evidencing that the positive charge density is much higher at the lower pH and, thus, the protonation degree α. This observation immediately implies that α at pH 7.5 is much lower than 100%, in contradiction to the estimation from the Henderson–Hasselbalch equation. The slight difference in the intensities observed with the bare subphases at the two pH values can be attributed to the slightly different concentration of Br as a result of pH titration with HBr (see Sec. II).

FIG. 5.

(a) X-ray fluorescence spectra around the Kα emission peak of the Br counterions obtained with the bare aqueous subphases and with lipid monolayers at two pH values. (b) Illustration of the model used for data analysis (here pH 5.0): the two-slab representation of the experimentally obtained electron density ρe(z), the standing wave intensity ϕ(z), and the counterion profile c(z). For the sake of clarity, all curves are normalized to have maximal values of 1.

FIG. 5.

(a) X-ray fluorescence spectra around the Kα emission peak of the Br counterions obtained with the bare aqueous subphases and with lipid monolayers at two pH values. (b) Illustration of the model used for data analysis (here pH 5.0): the two-slab representation of the experimentally obtained electron density ρe(z), the standing wave intensity ϕ(z), and the counterion profile c(z). For the sake of clarity, all curves are normalized to have maximal values of 1.

Close modal
The following quantitative analysis takes the counterion profiles c(z) from the MD simulations into account. This approach can be viewed as a refinement with respect to commonly used treatments of the counterion distribution in terms of thin adsorption layers27,36 or related approximations.59 Here, the simulation-based distributions were modeled as a (truncated) Gaussian function on top of the constant experimental bulk concentration c0,
c(z)=c0+(cmaxc0)e(zzmax)2/(2σ2)zz0,
(15)
where z0 defines the interface between the tail and headgroup slabs. This functional form was fitted to the numerical counterion profiles extracted from the simulations [see solid line in Fig. 4(c)], yielding zmax = z0 + 0.56 nm; σ = 0.48 nm for pH 5.0 and zmax = z0 + 0.53 nm; σ = 0.47 nm for pH 7.5. With that, cmax remains as the only adjustable parameter, which is varied so as to match the measured fluorescence intensities IF in the presence of the lipid monolayer and at the given pH, according to Eq. (2). For illustration, Fig. 5(b) displays all ingredients required for this calculation: the two-slab model of the electron density ρe(z), obtained from the experimental GIXOS curves as described in the supplementary material (Fig. S3), the corresponding standing wave intensity ϕ(z), and the parameterized counterion profile c(z). The counterion excess then simply follows as
Γ=z0c(z)c0dz.
(16)
We obtain Γ = 0.185 nm−2 for pH 5.0 and Γ = 0.039 nm−2 for pH 7.5. For a known counterion excess, the protonation degree α of MC3 follows as:
α=ΓAlipfMC3
(17)
where fMC3 = 0.2 is the MC3 mole fraction. With that, we obtain α = 67.5% for pH 5.0 and α = 14.5% for pH 7.5.

The obtained protonation degrees are significantly lower compared to the prediction from the Henderson–Hasselbalch equation. This is not surprising since the equation is valid only in the limit of infinite dilution. In the monolayer system, interactions between polar lipids, charged lipids, ions and protons contribute to the protonation degree. The effect of the proton distribution is taken into account in the self-consistent Eq. (8). The protonation degrees predicted in this way are α = 99% for pH 5.0 and α = 68% for pH 7.5. Although these values provide a more meaningful estimate, they are still significantly higher than the ones experimentally observed and require an a priori estimate of the surface charge density. To provide improvement, the interactions between the lipids in the monolayer as well as conformational rearrangement and partial dehydration have to be taken into account. In particular the effect of partial dehydration is expected to have a significant effect since the uncharged headgroups are buried while the charged headgroups are solvent exposed [Figs. 6(a) and 6(b)]. In any case, reliable experimental measurements of the protonation degree are imperative for validation and improvement of theoretical predictions.

FIG. 6.

Simulations with experimentally validated Alip and α. (a) Side view snapshot of one of the monolayers with α = 67.5%. Colors: green (MC30), yellow (MC3+) and gray (POPC). (b) Number density profiles ρn(z) of water and of the head and tail groups of MC30, MC3+, and POPC. Solid lines represent head groups and dashed lines represent tail groups. Dotted line: water. The top scale applies to MC3, the bottom scale to POPC and water.

FIG. 6.

Simulations with experimentally validated Alip and α. (a) Side view snapshot of one of the monolayers with α = 67.5%. Colors: green (MC30), yellow (MC3+) and gray (POPC). (b) Number density profiles ρn(z) of water and of the head and tail groups of MC30, MC3+, and POPC. Solid lines represent head groups and dashed lines represent tail groups. Dotted line: water. The top scale applies to MC3, the bottom scale to POPC and water.

Close modal

Finally, it is important to confirm that the simulations with the experimentally validated protonation degree still reproduce the experimental GIXOS curves and counterion profiles that were used in the validation procedure (Fig. 2). Otherwise additional iterations would be required. The results show that after only one iteration the results converged and the simulated GIXOS curves do not significantly respond to the updated α values (see Fig. S6 in the supplementary material), and the peak in the counterion distribution still has the same position and width (see Fig. S7 in the supplementary material).

As the last step, MD simulations with the experimentally validated values of Alip and α for the two pH values were performed, with the parameters summarized in Table II. The correct protonation degrees were imposed by representing the MC3 molecules with suitable numbers of MC30 and MC3+ variants. In practice, 6 and 27 out of the 40 MC3 molecules in each monolayer were represented with MC3+ for pH 7.5 and pH 5.0, respectively, yielding α = 15.0% and α = 67.5%. These protonation degrees come as close as possible to the experimental values when considering that the numbers of lipids can only be varied in integer steps.

Figure 6(a) shows a snapshot from an equilibrated simulation with α = 67.5%, corresponding to pH 5.0. POPC matrix lipids are represented in white, MC30 in yellow, and MC3+ in green. Panel (b) shows the normalized number density profiles ρn of the tails and headgroups of all lipid species and of water in the z-direction, i.e., perpendicular to the membrane plane. The distribution of the positively charged headgroups of MC3+ largely overlaps with that of POPC’s zwitterionic headgroups, while the distribution of the uncharged and only weakly polar headgroups of MC30 is strongly shifted towards the hydrophobic region accommodating the chains. As shown in the supplementary material (Fig. S8), consistent observations can be made at α = 15.0%, corresponding to pH 7.5. This behavior has been noticed earlier43 and may provide the explanation for the small yet significant increase in Alip when increasing the pH from 5.0 to 7.5 (see above). Namely, in contrast to the naive expectation that a reduction in the density of charged lipids would decrease lateral repulsion and thus the lipids’ effective area requirement, the deep insertion of uncharged MC30 into the monolayer effectively increases the area requirement at high pH, where MC30 is more abundant.

We now turn to the in-plane (xy) distribution of the three lipid species. Figure 7(a) shows a top-view simulation snapshot of a monolayer after equilibration. Qualitatively, it is seen that both MC30 and MC3+, depicted in green and yellow, respectively, distribute rather evenly over the membrane plane, which maximizes their average distance. This behavior is quantified with the help of the normalized in-plane radial distribution functions (RDFs) of the headgroup N atoms, which are shown in panel (b). Indeed, while the RDF between two POPC headgroups on approach oscillates around unity and even exhibits a pronounced next-neighbor peak at r ≈ 0.9 nm due to steric and dipole interactions,60 the RDF between two MC3+ headgroups systematically drops below unity on approach, albeit with a few weak oscillations. This drop in the RDF can be attributed to the electrostatic repulsion between the positively charged headgroups. The RDF between two MC30 drops below unity on approach, too, but like-charge repulsion obviously cannot be the reason because MC30 carry no net charge. Calculation of the mean force reveals an oscillatory behavior (Fig. S10). Repulsion occurs at small distances since MC30 leads to a local arrangement of the POPC molecules which oppose to some extent the penetration of other MC30 molecules nearby.

FIG. 7.

Simulations with experimentally validated Alip and α. (a) Top view snapshot of one of the monolayers with α = 67.5%. Colors: green (MC30), yellow (MC3+) and gray (POPC). (b) In-plane radial distribution functions (RDFs) between the headgroup nitrogens of the same lipid species. (c) RDFs between headgroup nitrogens of different lipid species.

FIG. 7.

Simulations with experimentally validated Alip and α. (a) Top view snapshot of one of the monolayers with α = 67.5%. Colors: green (MC30), yellow (MC3+) and gray (POPC). (b) In-plane radial distribution functions (RDFs) between the headgroup nitrogens of the same lipid species. (c) RDFs between headgroup nitrogens of different lipid species.

Close modal

Another qualitative difference between the RDF of two MC30 and the other two RDFs shown in Fig. 7(b) is its finite value at r = 0. Apparently, two MC30 headgroups can share the same in-plane (xy) coordinates because they are somewhat mobile in the third dimension (z), even though their depth profile in Fig. 6(b) is not significantly broader than those of the other headgroups. To resolve this apparent discrepancy we have calculated the 3D-RDFs of the MC30 headgroups and compared them to the 3D-RDFs of the other headgroups (Fig. S9 in the supplementary material). It is seen that the radius of mutual steric exclusion is much smaller for MC30 than for the other headgroups, which explains why the probability for two MC30 headgroups to be at the same xy position is higher than that for the other headgroups.

Next we take a look at the mixed RDFs shown in Fig. 7(c). Strikingly, the mixed RDFs of MC30 with POPC and of MC30 with MC3+ both merely exhibit a weak dip near r = 0, demonstrating that there is little steric exclusion because the headgroups of MC30 roam at a different “altitude” than those of POPC and MC3+. In contrast, as seen from the mixed RDF of POPC with MC3+, their headgroups mutually exclude each other sterically for small r and exhibit a pronounced next-neighbor peak at r ≈ 0.9 nm due to attractive interactions between the positively charged MC3+ headgroups and the negatively charged phosphate groups in the PC headgroups.61 

Molecular dynamics simulations of ionizable lipids provide detailed insight into their structure but require an a priori knowledge of the lipid packing and protonation degree. We have introduced a methodology of combining all-atom molecular dynamics simulations with experimental total reflection x-ray scattering and x-ray fluorescence to make robust and reliable predictions of MC3 lipid area and degree of protonation at pH conditions that mimic the extracellular and endosomal milieus. The combined approach yields a protonation degree α = 14.5% and an average area per molecule Alip = 0.73 nm2 at pH 7.5. At pH 5.0, it yields α = 67.5% and Alip = 0.75 nm2. In both cases the value for the protonation degree differs significantly from the value expected from the Henderson–Hasselbalch equation. Hence, the presetting of the protonation degree in the simulations based on simple theories can fail due to the interactions between polar lipids, charged lipids, ions and protons in complex lipid interfaces.

The simulations with optimized Alip and α, obtained from the matching methodology, reveal that the headgroups of the uncharged MC3 enter deeply into the hydrophobic monolayer region while the charged MC3 headgroups are solvent-exposed, as seen from the density profiles in Fig. 6. Lipid dehydration upon deprotonation is therefore expected to have a significant effect on the protonation degree and is likely the reason why theoretical predictions based on the Henderson–Hasselbalch equation overestimate the protonation degree for these simple monolayer systems. In addition, the buried uncharged MC3 lipids are laterally almost uncorrelated from the more hydrophilic POPC headgroups and the charged MC3 lipids. This behavior is associated with a non-typical structural response of Alip to pH variations, which may play a significant role for the endosomal escape and subsequent biological action of mRNA.

The joint experimental and simulation approach presented here provides a benchmark for consistent modeling. The results can be used to test theoretical models or advanced simulation methods such as free energy simulations or constant-pH simulations. The methodology is also suited to investigate changes in protonation and lipid packing in response to interactions with nucleic acids. Another challenge is to extend this methodology for predictions in three-dimensional lipid mesophases and lipid nanoparticles.

Further details on pH titration and GIXOS data fitting; isotherm measurements; TRXF spectra; Influence of protonation degree on modeled GIXOS curves and on the profiles of counterions and other chemical components; 3D-RDFs; PMFs of MC30 interactions.

We acknowledge DESY (Hamburg, Germany), a member of the Helmholtz Association HGF, for the provision of experimental facilities. Parts of this research were carried out at PETRA III and we would like to thank Chen Shen for assistance in using P08. Beamtime was allocated for proposal I-11009842. We thank Gerald Brezesinski and Tetiana Mukhina for help with the synchrotron experiments. This research was supported by the Bundesministerium für Bildung und Forschung (BMBF) within the Röntgen-Ångström cluster project Medisoft, Grant Nos. 05K18EZA and 05K18WMA, and by the Emmy Noether program of the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation), Grant No. 315221747. We acknowledge the GOETHE HLR for super-computing access. The authors gratefully acknowledge the scientific support and HPC resources provided by the Erlangen National High Performance Computing Center (NHR@FAU) of the Friedrich-Alexander-Universität Erlangen-Nürnberg (FAU) under the NHR project b119ee. NHR funding is provided by federal and Bavarian state authorities. NHR@FAU hardware is partially funded by the German Research Foundation (DFG), Grant No. 440719683.

The authors have no conflicts to disclose.

Miriam Grava: Formal analysis (equal); Investigation (equal); Methodology (equal); Validation (equal); Writing – original draft (equal); Writing – review & editing (equal). Mohd Ibrahim: Formal analysis (equal); Investigation (equal); Methodology (equal); Validation (equal); Writing – original draft (equal); Writing – review & editing (equal). Akhil Sudarsan: Formal analysis (equal); Investigation (equal); Methodology (equal); Writing – original draft (equal); Writing – review & editing (equal). Julio Pusterla: Investigation (equal); Methodology (equal); Writing – review & editing (equal). Julian Philipp: Investigation (equal); Methodology (equal); Writing – review & editing (equal). Joachim O. Rädler: Conceptualization (equal); Writing – review & editing (equal). Nadine Schwierz: Conceptualization (equal); Formal analysis (equal); Investigation (equal); Methodology (equal); Project administration (equal); Validation (equal); Writing – original draft (equal); Writing – review & editing (equal). Emanuel Schneck: Conceptualization (equal); Formal analysis (equal); Investigation (equal); Methodology (equal); Project administration (equal); Validation (equal); Writing – original draft (equal); Writing – review & editing (equal).

The data that support the findings of this study are available from the corresponding author upon reasonable request.

1.
Z.
Li
,
J.
Carter
,
L.
Santos
,
C.
Webster
,
C. F.
van der Walle
,
P.
Li
,
S. E.
Rogers
, and
J. R.
Lu
, “
Acidification-induced structure evolution of lipid nanoparticles correlates with their in vitro gene transfections
,”
ACS Nano
17
,
979
(
2023
).
2.
X.
Hou
,
T.
Zaks
,
R.
Langer
, and
Y.
Dong
, “
Lipid nanoparticles for mRNA delivery
,”
Nat. Rev. Mater.
6
,
1078
1094
(
2021
).
3.
L.
Miao
,
Y.
Zhang
, and
L.
Huang
, “
mRNA vaccine for cancer immunotherapy
,”
Mol. Cancer
20
,
41
(
2021
).
4.
M. J.
Carrasco
,
S.
Alishetty
,
M.-G.
Alameh
,
H.
Said
,
L.
Wright
,
M.
Paige
,
O.
Soliman
,
D.
Weissman
,
T. E.
Cleveland
,
A.
Grishaev
, and
M. D.
Buschmann
, “
Ionization and structural properties of mRNA lipid nanoparticles influence expression in intramuscular and intravascular administration
,”
Commun. Biol.
4
,
956
(
2021
).
5.
M.
Jayaraman
,
S. M.
Ansell
,
B. L.
Mui
,
Y. K.
Tam
,
J.
Chen
,
X.
Du
,
D.
Butler
,
L.
Eltepu
,
S.
Matsuda
,
J. K.
Narayanannair
,
K. G.
Rajeev
,
I. M.
Hafez
,
A.
Akinc
,
M. A.
Maier
,
M. A.
Tracy
,
P. R.
Cullis
,
T. D.
Madden
,
M.
Manoharan
, and
M. J.
Hope
, “
Maximizing the potency of siRNA lipid nanoparticles for hepatic gene silencing in vivo
,”
Angew. Chem.
124
,
8657
8661
(
2012
).
6.
C.
Tanford
and
J. G.
Kirkwood
, “
Theory of protein titration curves. I. General equations for impenetrable spheres
,”
J. Am. Chem. Soc.
79
,
5333
5339
(
1957
).
7.
B.
Roux
and
T.
Simonson
, “
Implicit solvent models
,”
Biophys. Chem.
78
,
1
20
(
1999
).
8.
R. W.
Zwanzig
, “
High-temperature equation of state by a perturbation method. I. Nonpolar gases
,”
J. Chem. Phys.
22
,
1420
1426
(
1954
).
9.
J. G.
Kirkwood
, “
Statistical mechanics of fluid mixtures
,”
J. Chem. Phys.
3
,
300
313
(
1935
).
10.
A.
Warshel
,
F.
Sussman
, and
G.
King
, “
Free energy of charges in solvated proteins: Microscopic calculations using a reversible charging process
,”
Biochemistry
25
,
8368
8372
(
1986
).
11.
T.
Simonson
,
J.
Carlsson
, and
D. A.
Case
, “
Proton binding to proteins: pKa calculations with explicit and implicit solvent models
,”
J. Am. Chem. Soc.
126
,
4167
4180
(
2004
).
12.
A. M.
Baptista
, “
Comment on ‘Explicit-solvent molecular dynamics simulation at constant pH: Methodology and application to small amines’ [J. Chem. Phys. 114, 9706 (2001)]
,”
J. Chem. Phys.
116
,
7766
7768
(
2002
).
13.
S.
Donnini
,
R. T.
Ullmann
,
G.
Groenhof
, and
H.
Grubmüller
, “
Charge-neutral constant pH molecular dynamics simulations using a parsimonious proton buffer
,”
J. Chem. Theory Comput.
12
,
1040
1051
(
2016
).
14.
Philipp
et al, “
pH-dependent structural transitions in cationic ionizable lipid mesophases are critical for lipid nanoparticle function
,” (
2023
) (submitted).
15.
S. L.
Duncan
and
R. G.
Larson
, “
Comparing experimental and simulated pressure-area isotherms for DPPC
,”
Biophys. J.
94
,
2965
2986
(
2008
).
16.
J. B.
Klauda
,
R. M.
Venable
,
J. A.
Freites
,
J. W.
O’Connor
,
D. J.
Tobias
,
C.
Mondragon-Ramirez
,
I.
Vorobyov
,
A. D.
MacKerell
, Jr.
, and
R. W.
Pastor
, “
Update of the CHARMM all-atom additive force field for lipids: Validation on six lipid types
,”
J. Phys. Chem. B
114
,
7830
7843
(
2010
).
17.
A.
Lamberg
and
O. H. S.
Ollila
, “
Comment on ‘Structural properties of POPC monolayers under lateral compression: Computer simulations analysis
,’”
Langmuir
31
,
886
887
(
2015
).
18.
C.
Tempra
,
O. H. S.
Ollila
, and
M.
Javanainen
, “
Accurate simulations of lipid monolayers require a water model with correct surface tension
,”
J. Chem. Theory Comput.
18
,
1862
1869
(
2022
).
19.
C. A.
Helm
,
H.
Möhwald
,
K.
Kjaer
, and
J.
Als-Nielsen
, “
Phospholipid monolayer density distribution perpendicular to the water surface. A synchrotron x-ray reflectivity study
,”
Europhys. Lett.
4
,
697
(
1987
).
20.
L.
Wiegart
,
B.
Struth
,
M.
Tolan
, and
P.
Terech
, “
Thermodynamic and structural properties of phospholipid Langmuir monolayers on hydrosol surfaces
,”
Langmuir
21
,
7349
7357
(
2005
).
21.
E.
Schneck
,
T.
Schubert
,
O. V.
Konovalov
,
B. E.
Quinn
,
T.
Gutsmann
,
K.
Brandenburg
,
R. G.
Oliveira
,
D. A.
Pink
, and
M.
Tanaka
, “
Quantitative determination of ion distributions in bacterial lipopolysaccharide membranes by grazing-incidence x-ray fluorescence
,”
Proc. Natl. Acad. Sci. U. S. A.
107
,
9147
9151
(
2010
).
22.
M.
Sturm
,
O.
Gutowski
, and
G.
Brezesinski
, “
The effect of pH on the structure and lateral organization of cardiolipin in Langmuir monolayers
,”
ChemPhysChem
23
,
e202200218
(
2022
).
23.
V. L.
Shapovalov
and
G.
Brezesinski
, “
Breakdown of the Gouy–Chapman model for highly charged Langmuir Monolayers: Counterion size effect
,”
J. Phys. Chem. B
110
,
10032
10040
(
2006
).
24.
V. L.
Shapovalov
,
M. E.
Ryskin
,
O. V.
Konovalov
,
A.
Hermelink
, and
G.
Brezesinski
, “
Elemental analysis within the electrical double layer using total reflection x-ray fluorescence technique
,”
J. Phys. Chem. B
111
,
3927
3934
(
2007
).
25.
P. R.
Cullis
and
M. J.
Hope
, “
Lipid nanoparticle systems for enabling gene therapies
,”
Mol. Ther.
25
,
1467
1475
(
2017
).
26.
M.
Kanduč
,
A.
Schlaich
,
A. H.
de Vries
,
J.
Jouhet
,
E.
Maréchal
,
B.
Demé
,
R. R.
Netz
, and
E.
Schneck
, “
Tight cohesion between glycolipid membranes results from balanced water–headgroup interactions
,”
Nat. Commun.
8
,
14899
(
2017
).
27.
J.
Pusterla
,
E.
Scoppola
,
C.
Appel
,
T.
Mukhina
,
C.
Shen
,
G.
Brezesinski
, and
E.
Schneck
, “
Characterization of lipid bilayers adsorbed to functionalized air/water interfaces
,”
Nanoscale
14
,
15048
15059
(
2022
).
28.
J.
Als-Nielsen
and
D.
McMorrow
,
Elements of Modern X-Ray Physics
(
John Wiley & Sons
,
2011
).
29.
W.
Yun
and
J.
Bloch
, “
X-ray near total external fluorescence method: Experiment and analysis
,”
J. Appl. Phys.
68
,
1421
1428
(
1990
).
30.
S.
Mora
,
J.
Daillant
,
D.
Luzet
, and
B.
Struth
, “
X-ray surface scattering investigation of Langmuir films: Phase transitions and elastic properties
,”
Europhys. Lett.
66
,
694
(
2004
).
31.
R. G.
Oliveira
,
E.
Schneck
,
B. E.
Quinn
,
O. V.
Konovalov
,
K.
Brandenburg
,
T.
Gutsmann
,
T.
Gill
,
C. B.
Hanna
,
D. A.
Pink
, and
M.
Tanaka
, “
Crucial roles of charged saccharide moieties in survival of gram negative bacteria against protamine revealed by combination of grazing incidence x-ray structural characterizations and Monte Carlo simulations
,”
Phys. Rev. E
81
,
041901
(
2010
).
32.
M.
Kanduč
,
E.
Schneck
, and
C.
Stubenrauch
, “
Intersurfactant h-bonds between head groups of n-dodecyl-β-d-maltoside at the air–water interface
,”
J. Colloid Interface Sci.
586
,
588
595
(
2021
).
33.
G. H.
Vineyard
, “
Grazing-incidence diffraction and the distorted-wave approximation for the study of surfaces
,”
Phys. Rev. B
26
,
4146
(
1982
).
34.
L. G.
Parratt
, “
Surface studies of solids by total reflection of x-rays
,”
Phys. Rev.
95
,
359
(
1954
).
35.
F.
Abelès
and
F. A.
La
, “
La théorie générale des couches minces
,”
J. Phys. Radium
11
,
307
309
(
1950
).
36.
G.
Brezesinski
and
E.
Schneck
, “
Investigating ions at amphiphilic monolayers with x-ray fluorescence
,”
Langmuir
35
,
8531
8542
(
2019
).
37.
A.
Thompson
,
D.
Vaughan
et al,
X-Ray Data Booklet, Tables 1–3
(
Lawrence Berkeley National Laboratory, University of California
,
CA
,
2001
); https://xdb.lbl.gov/xdb.pdf.
38.
J.
Caspers
,
E.
Goormaghtigh
,
J.
Ferreira
,
R.
Brasseur
,
M.
Vandenbranden
, and
J.-M.
Ruysschaert
, “
Acido–basic properties of lipophilic substances: A surface potential approach
,”
J. Colloid Interface Sci.
91
,
546
551
(
1983
).
39.
N.
Schwierz
,
D.
Horinek
, and
R. R.
Netz
, “
Specific ion binding to carboxylic surface groups and the pH dependence of the Hofmeister series
,”
Langmuir
31
,
215
225
(
2015
).
40.
M. J.
Abraham
,
T.
Murtola
,
R.
Schulz
,
S.
Páll
,
J. C.
Smith
,
B.
Hess
, and
E.
Lindahl
, “
GROMACS: High performance molecular simulations through multi-level parallelism from laptops to supercomputers
,”
SoftwareX
1–2
,
19
25
(
2015
).
41.
C. J.
Knight
and
J. S.
Hub
, “
MemGen: A general web server for the setup of lipid membrane simulation systems
,”
Bioinformatics
31
,
2897
2899
(
2015
); https://academic.oup.com/bioinformatics/article-pdf/31/17/2897/17085151/btv292.pdf.
42.
D. A.
Case
, et al.,
Amber
2021 (University of California, San Francisco.,
2021
).
43.
M.
Ibrahim
,
J.
Gilbert
,
M.
Heinz
,
T.
Nylander
, and
N.
Schwierz
, “
Structural insights on ionizable Dlin-MC3-DMA lipids in DOPC layers by combining accurate atomistic force fields, molecular dynamics simulations and neutron reflectivity
,”
Nanoscale
15
,
11647
11656
(
2023
).
44.
S.
Mamatkulov
and
N.
Schwierz
, “
Force fields for monovalent and divalent metal cations in TIP3P water based on thermodynamic and kinetic properties
,”
J. Chem. Phys.
148
,
074504
(
2018
).
45.
W. L.
Jorgensen
,
J.
Chandrasekhar
,
J. D.
Madura
,
R. W.
Impey
, and
M. L.
Klein
, “
Comparison of simple potential functions for simulating liquid water
,”
J. Chem. Phys.
79
,
926
935
(
1983
).
46.
G.
Bussi
,
D.
Donadio
, and
M.
Parrinello
, “
Canonical sampling through velocity rescaling
,”
J. Chem. Phys.
126
,
014101
(
2007
).
47.
M.
Parrinello
and
A.
Rahman
, “
Polymorphic transitions in single crystals: A new molecular dynamics method
,”
J. Appl. Phys.
52
,
7182
7190
(
1981
).
48.
U.
Essmann
,
L.
Perera
,
M. L.
Berkowitz
,
T.
Darden
,
H.
Lee
, and
L. G.
Pedersen
, “
A smooth particle mesh Ewald method
,”
J. Chem. Phys.
103
,
8577
8593
(
1995
).
49.
B.
Hess
,
H.
Bekker
,
H. J. C.
Berendsen
, and
J. G. E. M.
Fraaije
, “
LINCS: A linear constraint solver for molecular simulations
,”
J. Comput. Chem.
18
,
1463
1472
(
1997
).
50.
N.
Michaud-Agrawal
,
E. J.
Denning
,
T. B.
Woolf
, and
O.
Beckstein
, “
MDAnalysis: A toolkit for the analysis of molecular dynamics simulations
,”
J. Comput. Chem.
32
,
2319
2327
(
2011
).
51.
W.
Humphrey
,
A.
Dalke
, and
K.
Schulten
, “
VMD: Visual molecular dynamics
,”
J. Mol. Graphics
14
,
33
38
(
1996
).
52.
C.
Vega
and
E.
de Miguel
, “
Surface tension of the most popular models of water by using the test-area simulation method
,”
J. Chem. Phys.
126
,
154707
(
2007
).
53.
P.
Müller
,
D. J.
Bonthuis
,
R.
Miller
, and
E.
Schneck
, “
Ionic surfactants at air/water and oil/water interfaces: A comparison based on molecular dynamics simulations
,”
J. Phys. Chem. B
125
,
406
415
(
2021
).
54.
A. J.
Nelson
and
S. W.
Prescott
, “
refnx: neutron and x-ray reflectometry analysis in python
,”
J. Appl. Crystallogr.
52
,
193
200
(
2019
).
55.
D.
Marsh
, “
Lateral pressure in membranes
,”
Biochim. Biophys. Acta, Rev. Biomembr.
1286
,
183
223
(
1996
).
56.
J.
Bloch
,
W.
Yun
,
X.
Yang
,
M.
Ramanathan
,
P.
Montano
, and
C.
Capasso
, “
Adsorption of counterions to a stearate monolayer spread at the water-air interface: A synchrotron x-ray study
,”
Phys. Rev. Lett.
61
,
2941
(
1988
).
57.
J.
Daillant
,
L.
Bosio
,
J.
Benattar
, and
C.
Blot
, “
Interaction of cations with a fatty acid monolayer. A grazing incidence x-ray fluorescence and reflectivity study
,”
Langmuir
7
,
611
614
(
1991
).
58.
J. F.
Nagle
and
S.
Tristram-Nagle
, “
Structure of lipid bilayers
,”
Biochim. Biophys. Acta, Rev. Biomembr.
1469
,
159
195
(
2000
).
59.
C.
Stefaniu
,
V. M.
Latza
,
O.
Gutowski
,
P.
Fontaine
,
G.
Brezesinski
, and
E.
Schneck
, “
Headgroup-ordered monolayers of uncharged glycolipids exhibit selective interactions with ions
,”
J. Phys. Chem. Lett.
10
,
1684
1690
(
2019
).
60.
V. M.
Kaganer
,
H.
Möhwald
, and
P.
Dutta
, “
Structure and phase transitions in Langmuir monolayers
,”
Rev. Mod. Phys.
71
,
779
(
1999
).
61.
R.
Zantl
,
L.
Baicu
,
F.
Artzner
,
I.
Sprenger
,
G.
Rapp
, and
J. O.
Rädler
, “
Thermotropic phase behavior of cationic Lipid–DNA complexes compared to binary lipid mixtures
,”
J. Phys. Chem. B
103
,
10300
10310
(
1999
).

Supplementary Material