Open Access Article
This Open Access Article is licensed under a Creative Commons Attribution-Non Commercial 3.0 Unported Licence

The effect of hydronium ions on the structure of phospholipid membranes

Evelyne Deplazes *a, David Poger b, Bruce Cornell c and Charles G. Cranfield d
aSchool of Biomedical Sciences, Curtin Health Innovation Research Institute and Curtin Institute for Computation, Curtin University, Perth, WA 6845, Australia. E-mail:
bSchool of Chemistry and Molecular Biosciences, The University of Queensland, Brisbane, QLD 4072, Australia
cSDx Tethered Membranes Pty. Ltd., Roseville, NSW 2069, Australia
dSchool of Life Sciences, University of Technology Sydney, Ultimo, NSW 2007, Australia

Received 4th October 2017 , Accepted 22nd November 2017

First published on 27th November 2017

This work seeks to identify the mechanisms by which hydronium ions (H3O+) modulate the structure of phospholipid bilayers by studying the interactions of H3O+ with phospholipids at the molecular level. For this, we carried out multiple microsecond-long unrestrained molecular dynamics (MD) simulations of a POPC bilayer at different H3O+ concentrations. The results show that H3O+ accumulates at the membrane surface where it displaces water and forms strong and long-lived hydrogen bonds with the phosphate and carbonyl oxygens in phospholipids. This results in a concentration-dependent reduction of the area per lipid and an increase in bilayer thickness. This study provides an important molecular-level insight into the mechanism of how H3O+ modulates the structure of biological membranes and is a critical step towards a better understanding of the effect of low pH on mammalian and bacterial membranes.

1 Introduction

Phospholipids are the major components of cell membranes in all living organisms.1 Phospholipid bilayers are thus routinely used as model systems for mammalian and bacterial membranes to investigate their structural, chemical and physical properties.2 Besides lipid composition, environmental factors such as temperature, salt concentration, ionic strength and pH affect the properties of membranes and lipid bilayers, including their structure, level of fluidity, interfacial tension, conduction and ion permeability. The effect of temperature,3–5 pressure6,7 and salt8–11 on lipid bilayers has been studied extensively over the past decades. In contrast, less is known on the effect of pH on the structure of membranes despite its important role in a range of biological processes, including the stability and integrity of the lysosomal membrane,12 the formation of cholesterol-enriched domains13 and drug partitioning.14,15 Furthermore, extracellular acidosis has been associated with inflammation16 and various pathological states (for example tumour growth,17,18 ischaemic stroke19 and epileptic seizures20). The effect of pH on biological membranes is also relevant to understanding how some organisms have adapted to living in extreme pH environments21 and the design of acidophilic bacteria for bio-mining and other biotechnology applications.22

Most studies on the effect of pH on membranes have focused on the physical, electrical and mechanical properties of phospholipid bilayers rather than the molecular interactions of the hydronium ions with the phospholipids. As summarised by Koynova and Caffrey,23 at a pH < 3 in the bulk aqueous solution, a range of fully hydrated, saturated diacyl and dialkyl glycerophospholipids bilayers, including phosphatidylcholines and phosphatidylcholine analogues show an increased lamellar gel-to-liquid crystalline phase-transition temperature.24–29 Furthermore, using a range of experimental techniques such as micropipette aspiration, fluorescence spectroscopy and electrophoresis, Zhou and Raphael characterised the mechanical and electrical properties of a 2-oleoyl-1-stearoyl-sn-glycero-3-phosphocholine (SOPC) bilayer at pH 2–9.30 Based on measurements of ζ-potential and dipole moments, the authors concluded that pH-induced changes in membrane bending stiffness stemmed from alterations in interfacial rather than intra-membrane electrostatics. In another study, Petelska and Figaszewski proposed a theoretical model to describe the effect of pH on interfacial tension in the vicinity of the isoelectric point of a phosphatidylcholine lipid bilayer based on interfacial tension measurements of a solvent-formed egg phosphatidylcholine bilayer.31 The authors then extended their model to describe the adsorption of the H+ and OH ions to the surface of phosphatidylcholine and phosphatidylserine bilayers for pH ranges within and outside the isoelectric point.32 Whilst these studies provided some insight into the effect of pH on macroscopic properties of phospholipid bilayers, they did not provide structural information at the molecular level of solvent-free bilayers. That is, none of the methods used can report on the molecular interactions between lipids and H3O+ or OH ions, and how these interactions might differ from and compete with hydrogen bonding of lipids with H2O. Indeed, rather than studying the effect of protons and hydronium ions on the phospholipid bilayer, most experimental and theoretical studies have focused on the transfer of protons or hydronium ions across the surface of biological membranes, a process critical to energy metabolism.33–38 In addition, the intricate interactions of water with phospholipids, in particular how interfacial water can cross-link phospholipid head groups via hydrogen bonds, has been investigated experimentally and in simulations.9,39–44

It has been recognised for some time that water can cross-link phospholipid head groups via hydrogen bonding to the oxygens present at the interface.41–43 Glycerophospholipids provide a number of hydrogen-bonding sites. In phosphatidylcholines such as POPC (2-oleoyl-1-palmitoyl-sn-glycero-3-phosphocholine) these include the carbonyl oxygens in the sn-1 (referred to as O1A) and sn-2 chains (O2A) and the two non-ester phosphate oxygens (OP1 and OP2), as illustrated in Fig. 1. Note, the oxygens OP1 and OP2 are chemically equivalent through resonance and both P–O bonds have a 1.5-bond character owing to sharing of the negative charge. Recently, Cranfield et al. proposed a model to account for the pH-induced structural changes in model membranes through the modification of water–lipid hydrogen bonds.45 Specifically, the conductance and capacitance of tethered phospholipid bilayers in bathing solutions at acidic and alkaline pH were determined using AC impedance spectroscopy and DC-ramped amperometry. In addition, the thickness and water content of the bilayers were examined using neutron reflectometry. It was suggested that the variations in membrane conductivity and water penetration into the bilayer observed at low pH was caused by H3O+ ions competing with and disrupting water-bridged intermolecular hydrogen bonds between lipid molecules. Furthermore, it was hypothesised that the observed decrease in area per lipid was due to an enhanced stability of hydrogen bonds of H3O+ with the phosphate and carbonyl oxygens in lipids.

image file: c7cp06776c-f1.tif
Fig. 1 Structure of a POPC lipid molecule showing the hydrogen bonding sites on the carbonyl and phosphate oxygens (A) and chemical moieties (B). The oxygens OP1 and OP2 are chemically equivalent through resonance.

In this study, we examine the interactions of H3O+ ions with a POPC bilayer at a molecular level using multiple microsecond-long unrestrained molecular dynamics (MD) simulations at H3O+ concentrations of 0.4 and 0.04 M. While bulk-phase concentrations of 0.4 M and 0.04 M for the hydronium ion may seem high (equivalent to bulk-phase pH values of about 0.5 and 1.5, respectively), a range of experimental33,35,46,47 and computational48–51 studies have suggested that the H3O+ concentration is significantly higher on the membrane surface than in the bulk water phase, and that the interfacial water layer may act as a barrier for ions which would contribute to retarding surface-to-bulk H3O+ equilibration and separating the surface of a membrane from the bulk. Indeed, calculations of the potential of mean force of the hydrated excess proton near phospholipid bilayers in atomistic simulations have indicated that the interfacial pH could be as low as 3 when the pH of the bulk water region was taken as 7.51 This would imply that in the experiments performed by Cranfield et al.45 at a bulk-phase pH of 5, the local concentration in hydronium ions at the surface of the model membranes could be within a range comparable to that used in the simulations.

In summary, the results presented here show that H3O+ ions accumulate at the water–lipid interface where they displace water. The H3O+ ions form short and long-lived hydrogen bonds with the phosphate and carbonyl oxygens with a particular preference for the sn-2 carbonyl oxygen. The strong hydrogen bond network results in a decrease of the area per lipid and an increase in the bilayer thickness in agreement with previous experimental studies45 and consistent with the increase in the lamellar gel-to-liquid crystalline phase-transition temperature in phosphatidylcholines and phosphatidylcholine analogues.24–29 The combined results from these simulations are supportive of the hypothesis by Cranfield et al. that the strong hydrogen bonding capacity of H3O+ is the origin of the pH-induced changes in the structure of phospholipid membranes.45

2 Experimental

Model systems

Three simulation systems were examined with H3O+ concentrations, denoted [H3O+], of 0, 0.04 M and 0.4 M (Table 1). Each system consisted of a preassembled bilayer of 512 POPC molecules and at least 13[thin space (1/6-em)]000 water molecules; that is, at least 45 water molecules per lipid to ensure a fully hydrated state. The POPC bilayer in the absence of H3O+ ions was taken from previous studies.52,53 The system consisting of a POPC bilayer with 0.04 M and 0.4 M [H3O+] were prepared by randomly adding 16 and 164 H3O+ ions, respectively.
Table 1 Overview of simulation systems listing the number of H3O+ ions added and the resulting H3O+ concentrations. Each simulation system contains a preassembled bilayer of 512 POPC molecules. The simulations for the neutral system were taken from previous studies52,53
System name H3O+ ions [H3O+] Simulations (μs)
Neutral 2 × 1
0.04 M 16 0.04 M 3 × 1
0.4 M 164 0.4 M 3 × 1

Simulation parameters

All simulations were performed using GROMACS version 4.6.7,54 in conjunction with the GROMOS 54A7 force field55 and parameters for POPC developed by Poger et al.53 This force field has been shown to reproduce a range of lipid-packing, chain-ordering and hydration properties of phosphatidylcholine bilayers.52,53,56,57 Parameters for the H3O+ ion compatible with the GROMOS 54A7 force field were obtained from the Automated Topology Builder (ATB; http://,59 Briefly, the partial charges of the oxygen and hydrogen atoms were −0.587e and 0.529e, respectively, giving an overall charge of +1 to the ion. The O–H bond length and the H–O–H angle were set to 0.0983 nm and 109.50°, respectively (see also Fig. S1 in the ESI). The structure and parameter files can be downloaded from the ATB (molid 3859). The ATB has been thoroughly tested and validated using a large set of small molecules to reproduce hydration free enthalpies. Parameters are produced based on QM calculations and charge group partitioning60 optimised for small molecules used in biomolecular simulations. The model used for the hydronium ion in this study, in particular the atomic partial charges and the ion geometry, was similar to other models employed successfully in previous simulation studies.61–67 Such a model does not allow for bond-breaking or creation and therefore ignores any proton hopping as proposed in the Grotthuss mechanism.68,69 Likewise, it is assumed that the excess proton is localised on a single water molecule forming a hydronium cation, rather than being delocalised over multiple water molecules and changing between different species, such as the Eigen (H9O4+) and Zundel cations (H5O2+).68,70 Solvent water was described using the simple point charge (SPC) water model.71 All systems were subject to rectangular periodic boundary conditions. The systems were simulated using a 2 fs time step. Both temperature and pressure were maintained close to their reference values of 303 K and 1 bar, respectively, using the Berendsen weak-coupling method.72 For temperature coupling, a time constant of 0.1 ps was used. For semi-isotropic pressure coupling a time constant of 1 ps alongside an isothermal compressibility of 4.65 × 10−5 bar−1 were used. The LINCS algorithm73 was used to constrain the lengths of all bonds in lipids. Nonbonded interactions were evaluated using a single-range cutoff scheme whereby interactions within a 1.4 nm cutoff were calculated at every step and the pair list was updated every 5 steps. Previous work with the GROMOS force field used a twin-range cutoff scheme but its association with the current implementation of the integrator in GROMACS leads to artefacts in the collective properties of a lipid bilayer.74 Despite some loss in computational efficiency, using a single-range cutoff scheme was found to counterbalance those artefacts.74 To correct for the truncation of electrostatic interactions beyond 1.4 nm a reaction-field correction75 with a relative dielectric constant (εr) of 62 was applied. The systems with 0.04 and 0.4 M [H3O+] were simulated in triplicate and the reference system without H3O+ ions was simulated in duplicate, each time using randomly assigned starting velocities. All simulations were 1 μs long.

Changing the protonation state of the phosphate group in POPC as a result of the high concentration in hydronium ions used in the simulations was not considered insofar as the value of the pKa of the phosphate group in phosphatidylcholines has remained unclear given the various estimates reported including 0.8,42 <1,76 1.0,41 1.3,77 <1.528 and 2.25.78 Even for the isolated phosphate group alone, a pKa value as low as −1.54 has been proposed for dimethyl phosphate.51 In addition, those estimates relied on critical assumptions that intrinsically limited their accuracy. For example, Monceli et al.79 determined a value of pKa of 0.8 assuming that the area per lipid was constant and corresponded to that of a fluid-phase phosphatidylcholine bilayer (0.65 nm2), and that the properties of interfacial water, that is the water layer residing near the membrane surface, was the same as those of bulk water, in particular the dielectric permittivity. However, both of those assumptions have proven arguable since it was shown that variations in pH induced phase transition in lipid bilayers.23,29,30,80

Setup of membrane simulation systems

The systems containing 0.04 M and 0.4 M [H3O+] were generated by adding 16 and 164 H3O+ ions in the water slab, respectively, to the 521 POPC lipid bilayer at neutral pH (Table 1). The positive net charge of the systems containing H3O+ ions was neutralised by adding Cl anions. Additional Na+ and Cl ions were used to obtain a final concentration of 0.15 M NaCl. All the ions (H3O+, Na+ and Cl) were added to the hydrated systems by replacing water molecules selected randomly. Each system was energy-minimised using a steepest-descent algorithm and then simulated for 5 ns with all lipid atoms position-restrained to allow the water and ions to relax, followed by the production simulations.

Setup of ‘bulk’ water simulation systems

To validate the parameters of H3O+ a 40 ns simulation of ‘bulk’ water with 0.04 M [H3O+] and 0.15 M NaCl was carried out. For this, a cubic box of side lengths 6.5 nm containing 9081 water molecules and 6 randomly positioned H3O+ ions was created. The positive charge of the system was neutralised with 6 Cl anions and additional Na+ and Cl ions were used to obtain a final concentration of 0.15 M NaCl. The system was energy-minimised using a steepest-descent algorithm and then simulated for 40 ns saving frames every 2 ps.

Data analysis

Unless stated otherwise trajectories from independent simulations were analysed separately. Analysis was carried out using GROMACS tools54 and the MDAnalysis package.81 The area per lipid was calculated using the lateral dimensions of the box (in the xy plane) divided by the number of lipids per leaflet. The bilayer thickness was derived from the electron density profile of the phosphorus atoms of the system along the bilayer normal (taken as the z axis) and corresponded to the distance between the two maxima of the electron density. Distribution of distances between H3O+ ions were calculated using oxygen–oxygen distances from the last 500 ns of the trajectory. For hydrogen bonding analysis, a 0.35 nm cut-off for the radius and a 30° for the angle were used. Contacts between phosphorus atoms in the POPC molecules and H3O+ ions and water, respectively, were calculated using the centre of mass of the P atom, the H3O+ ions and water, respectively, and a 0.4 nm cut-off. Unless otherwise stated, uncertainties are given as ±1 standard deviation.

3 Results

The results section is divided into four parts: (i) validation of the H3O+ model, (ii) the effect of H3O+ ions on the structure of the lipid bilayer; (iii) the accumulation and location of the H3O+ ions in the lipid bilayer; and (iv) the hydrogen bonding network in the phosphate, glycerol and carbonyl of the lipid bilayer.

Validation of the H3O+ model

Over the past three decades a wide range of force-field parameters for the hydronium ion have been proposed. The geometric parameters (bond lengths and angles) used as well as the partial charges, van der Waals parameters, the presence or absence of polarisability, and the inclusion of virtual sites in these models vary considerably (see ref. 61–67, 82 and 83 and references therein). In the classical model presented here, the partial charges are 0.529e and −0.587e on the hydrogen and oxygen, respectively, and the Lennard-Jones parameters for the oxygen atom are ε = 0.84 kJ mol−1 and σ = 0.29 nm (see also Fig. S1, ESI). These parameters are based on DFT calculations (B3LYP/6-31G*) and charge group assignments that have been optimised for small molecules in biomolecular simulations,60 and implemented on the ATB.58,59 Also, our model of H3O+ is consistent with the force field and SPC water model used in our simulations of the hydrated bilayers.

As a validation of our model of H3O+, a 40 ns simulation of 0.04 M H3O+ in water with 0.15 M NaCl was carried out. The radial distribution function (RDF; Fig. S2, ESI) between the oxygen in H3O+ and the oxygen in water shows a first peak at 0.25 nm. This is in excellent agreement with earlier simulation studies of H3O+ in solution using classical non-polarisable force fields for H3O+.61,62 The RDF also agrees well with simulations using a more complex model of H3O+ that takes into account the interaction of an exchange charge distribution of the charge-transfer complex.83

Effect of H3O+ ions on membrane structure

The area per lipid AL and the bilayer thickness d are important properties to describe the structure and the degree of packing within a bilayer, which are thus indicative of the level of fluidity of a lipid bilayer.56Fig. 2 shows the time-dependent variation of AL and d calculated from the simulations of a POPC lipid bilayer with 0.04 M and 0.4 M H3O+ in comparison to the AL and d from the simulations of a POPC bilayer in the absence of H3O+ (referred to as neutral). As can be seen, all the systems reached equilibrium within about 0.5 μs and AL and d fluctuated around their average values over the last 0.5 μs. In the absence of H3O+ ions, the average values for AL and d over the last 0.5 μs of the two simulations were 0.628 ± 0.007 nm2 and 3.80 ± 0.16 nm, respectively. This is in good agreement with simulation and experimental values previously reported for a fluid-phase POPC bilayer.53,56,84,85 In contrast, AL decreased and d increased at higher H3O+ concentrations, that is AL = 0.572 ± 0.006 nm2 and d = 4.02 ± 0.14 nm at 0.04 M [H3O+], and AL = 0.514 ± 0.003 nm2 and d = 4.32 ± 0.13 nm, at 0.4 M [H3O+] over the last 0.5 μs of the three simulations.
image file: c7cp06776c-f2.tif
Fig. 2 Effect of H3O+ ions on the structure of a POPC lipid bilayer from 1 μs of unrestrained MD simulation. Area per lipid (AL) (A) and membrane thickness (d) (B) as a function of time from simulations of neutral POPC (grey, brown) and in the presence of 0.04 M [H3O+] (blue, magenta, indigo) and 0.4 M [H3O+] (black, light green, dark green). AL and d are shown as moving averages over 100 and 50 data points, respectively.

Accumulation of H3O+ ions at the water–lipid interface

At the start of the simulation, the H3O+ ions were randomly distributed in the bulk solvent. Shortly after the start of the simulation the H3O+ ions migrated to the water–lipid interface as reflected in the sharp increase in the average number of contacts formed between the phosphorus atom (P) and H3O+ ions at the cost of contacts between P and water molecules (Fig. 3). The time-dependent variation of these contacts suggests that the number of H3O+ ions at the water–lipid interface reached a steady state after about 300–400 ns. This is consistent with the convergence of the area per lipid and the bilayer thickness (Fig. 2), which further highlights that the simulations were at equilibrium over the last 0.5 μs of simulation. For the rest of the simulations the H3O+ ions remained at the water–lipid interface.
image file: c7cp06776c-f3.tif
Fig. 3 Number of contacts between phosphorus atoms in the POPC molecules and H3O+ ions and water, respectively, calculated from simulations of a POPC lipid bilayer in the presence of 0.4 M [H3O+] (A) and 0.04 M [H3O+] (B). Data from the three independent simulations is shown (grey, black and brown). The number of contacts are normalised by the number of lipids in the system.

To determine the location of the hydronium ions at the water–lipid interface, the electron density profiles of the phosphate, choline, sn-1 carbonyl and sn-2 carbonyl groups as well as water and H3O+ ions along the bilayer normal were calculated using data from the last 0.5 μs of the simulations. Fig. 4 displays the electron density profiles from simulations of POPC in the absence of H3O+ (A) and in the presence of 0.04 M [H3O+] (B) and 0.4 M [H3O+] (C). For clarity, the density of H3O+ from simulations with 0.04 M [H3O+] has been multiplied by 10. Comparing the density profiles of water in the three graphs shows a receding of the water density from the centre of the bilayer with increasing H3O+ concentration. This effect is absent in the simulations of POPC without H3O+ ions and thus further confirms that water at the membrane surface is being displaced by H3O+ ions. Furthermore, the position of the electron density for H3O+ reveals that H3O+ ions are located below the phosphate group close to the carbonyl groups.

image file: c7cp06776c-f4.tif
Fig. 4 Electron density profiles of a POPC lipid bilayer in the absence of H3O+ (A) and in the presence of 0.04 M [H3O+] (B) and of 0.4 M [H3O+] (C). Each density profile shows the electron density ρ as a function of the distance from the centre of the bilayer z for the phosphate (red), choline (green), sn-1 (magenta) and sn-2 (orange) carbonyl groups as well as water (blue) and H3O+ ions (black) calculated from the last 0.5 μs of the simulations. For clarity, the density of H3O+ from simulations with 0.04 M H3O+ (B) has been multiplied by 10.

Hydrogen bonding of H3O+

As illustrated in Fig. 1 the POPC molecule has several hydrogen-bonding sites through which it can interact with H3O+ and H2O. Fig. 5A shows histograms of hydrogen bonds per H3O+ ions for the different lipid oxygens calculated from the last 200 ns of the 0.04 M [H3O+] and 0.4 M [H3O+] simulations. The carbonyl oxygen in the sn-2 chain O2A formed between 2.0 to 2.5 hydrogen bonds per H3O+ while the non-ester phosphate oxygens OP2 and OP1 formed about 1.0 and 0.5 hydrogen bonds per H3O+, respectively. In contrast, the carbonyl oxygen in the sn-1 chain O1A formed almost no hydrogen bond with the H3O+ ion. These results suggest that H3O+ has a strong preference for the carbonyl oxygen in the sn-2 chain. The lifetime of the hydrogen bonds was estimated by calculating the occupancy, that is, the fraction of time during which the hydrogen bond is present. Fig. 5B lists the occupancy of hydrogen bonds formed by the different lipid oxygens with water and H3O+ ions. For a given oxygen, the occupancy was averaged over the three independent simulations of a POPC lipid bilayer in the presence of 0.04 [H3O+] and 0.4 [H3O+], respectively. As can be seen, the hydrogen bonds for all lipid oxygens exhibited greater occupancies for H3O+ than for water, in particular for O2A. For all oxygen, there is no significant difference between the 0.4 [H3O+] and 0.04 [H3O+] simulations suggesting that the effect is not dependent on the specific H3O+ concentration.
image file: c7cp06776c-f5.tif
Fig. 5 Hydrogen bonds (H-bonds) between lipid oxygens. (A) Histograms of H-bonds per H3O+ ions for the different lipid oxygens calculated from simulations of a POPC lipid bilayer in the presence of 0.04 M [H3O+] and 0.4 M [H3O+]. For each oxygen, the number of H-bonds were averaged over the last 200 ns of the simulation from the three independent trajectories. (B) Occupancy of H-bond formed by the different lipid oxygens with water and H3O+ ions, respectively. For each oxygen, the occupancy was calculated over the full lengths of the trajectory and then averaged over the three independent simulations a POPC lipid bilayer in the presence of 0.04 M [H3O+] and 0.4 M [H3O+], respectively. (C and D) Probability distribution of donor–acceptor distances distance for the hydrogen bonds formed by the lipid oxygens with the H3O+ ion (C) and water (D) calculated from the last 500 ns of simulations of POPC in the presence of 0.4 M [H3O+]. In summary, the analysis indicates that H3O+ forms stronger and longer-lived hydrogen bonds with the lipid oxygen compared to H2O.

Fig. 5C and D show the probability distribution of the donor–acceptor distance for the hydrogen bonds formed by the lipid oxygens with H3O+ (5C) and water (5D). Comparison of the two plots shows that for both H3O+ and water the average hydrogen bond distance is the shortest for O2A, followed by O1A. The bonds with OP1 and OP2 are the longest. However, for all four lipid oxygens the average hydrogen bond distances are slightly shorter for H3O+ compared to water. In addition, the distribution for the hydrogen bonds with H3O+ are narrower.

Bonding arrangement of lipids around the H3O+ ion

Next, we investigated the packing of lipids around H3O+ ions. To determine the most likely arrangement of lipids around H3O+ a list of all possible hydrogen-bonding arrangements for the H3O+ was prepared. For example, a H3O+ ion can be interacting with two O2A oxygens and one OP1 oxygen or one O2A oxygen and two OP1 oxygens or one O2A, one OP1 and one OP2 etc. The occurrence of each possible arrangements over the full-length of the 1 μs simulation was calculated as the frequency of finding all lipid oxygen within 0.2 nm of a given H3O+. The data was averaged over the three independent simulations with 0.04 M [H3O+] and 0.4 M [H3O+], respectively. The occurrence for the four most populated arrangements is shown in Table 2. Comparing the data from simulations with 0.04 M [H3O+] and 0.4 M [H3O+] shows that the standard deviations are much larger for simulations with 0.04 M [H3O+] due to the much lower number of the H3O+ ions (16 vs. 164). Nevertheless, the order of the four most populated arrangements is the same in two simulation systems.
Table 2 Occurrence for the different arrangements of lipid oxygen observed in simulation of POPC in the presence of 0.04 M [H3O+] and 0.4 M [H3O+]. Occurrences were calculated using the full-length of the 1 μs simulation and a 0.2 nm cut-off, and averaged over the three independent simulations of 0.04 M [H3O+] and 0.4 M [H3O+], respectively. In each arrangement, the H3O+ ion interacts with three oxygen from three lipid molecules. Errors are standard deviations
Arrangement Occupancy
0.04 M [H3O+] (%) 0.4 M [H3O+] (%)
2 O2A, 1 OP1 45 ± 11 36 ± 1
3 O2A 24 ± 11 19 ± 2
2 O2A, 1 OP2 7 ± 1 11 ± 1
1 O2A, 2 OP1 4 ± 3 6 ± 0

Fig. 5A and B show examples of the two arrangements with the highest occurrences. Together, these account for almost two thirds of hydrogen bonding arrangements across the 3 μs of the 0.04 M [H3O+] and 0.4 M [H3O+] simulations. In the most populated arrangement (42%), H3O+ interacts with two O2A oxygen and one OP1 oxygen (Fig. 5A). In the second most common arrangement (22%), H3O+ hydrogen bonds with three O2A oxygen (Fig. 5A). In both cases, the ion is surrounded by three lipids.

4 Discussion

Despite the important role of pH in biological processes and the association of acidosis with various pathophysiological states little is known about the effect of the H3O+ ion on the structure of phospholipid membranes. In this study, we aimed to characterise the molecular interactions with the H3O+ ion and phospholipids and its effect of the structure of a phospholipid bilayer. For this we carried out unrestrained MD simulations of POPC lipid bilayers in the presence 0.04 M and 0.4 M [H3O+] and combined this with data from simulations of POPC in the absence of H3O+.6

Analysis of the structure of the bilayers showed that in simulations with 0.04 M and 0.4 M H3O+ the area per lipid AL dropped by about 8% and 18% compared to a fluid-phase POPC bilayer without H3O+ ions (Fig. 2A). At the same time the bilayer thickness d increased by about 6% at 0.04 M H3O+ and 14% at 0.4 M H3O+ (Fig. 2B). This would indicate that the presence of H3O+ ions induced a reduction in the level of the fluidity of the bilayer. These observations are consistent with recent neutron reflectometry experiments on phospholipids that showed an increase in the outer leaflet tail thickness of about 4% as the pH was decreased from 9 to 5,45 as well as the gel-like characteristics and the reduced elasticity observed in SOPC bilayers at pH 2.30 In addition, a number of studies have shown that for bulk pH < 3, the lamellar gel-to-liquid crystalline phase-transition temperature for a range of phosphatidylcholines and phosphatidylcholine analogues is increased by about 10 °C.23,25,26,28,29,63 The concomitant change in AL and d can be interpreted in terms of the critical packing parameter (CPP). For a system of constant molecular volume, the molecular shape can be determined by the CPP defined as the ratio v/(a0l), where v is the hydrated molecular volume, a0 is the interfacial molecular area and l is the length of the lipid chains. The CPP can reflect a property of a heterogeneous molecular assembly where the individual components which may include lipids, sterols and peptides can combine to provide an averaged CPP. The membrane geometry influenced by the averaged CPP will modulate the physical properties of the membrane such as membrane thickness and conduction.86,87 That the POPC bilayer remained planar (CPP∼1) in all simulations at 0.04 M and 0.4 M [H3O+] suggests that the term a0l approximated as ALd was constant at both [H3O+] concentrations and that the increase in d closely matched the decrease in AL.

The analysis of the location and hydrogen bonding of H3O+ enabled us to identify the molecular interactions that underlie these changes. Results from the electron density profiles indicate that H3O+ sequesters into the water–lipid interface where it displaces water. The H3O+ ion accumulates slightly below the phosphate group close to the carbonyl oxygens (Fig. 4). This agrees with previous simulations studies of a proton excess or hydrated protons at the surface of DMPC,50 DOPC or DOPE bilayers.49,51 In these studies, the free energy profiles showed a strong membrane surface affinity for the protons and the H3O+ was shown to interact with both the phosphate and carbonyl oxygen. These observations from simulations are in line with the results from measurements of the ζ-potential by Zhou and Raphael30 which indicated that H3O+ ions alter the interfacial electrostatics by acting as counterions at the water–lipid interface. The same study also reported that there are no changes in the intra-membrane electric properties, as measured by the dipole moment. This suggest that H3O+ ions do not enter the hydrophobic core of the membrane, which is also in agreement with the data from our simulations.

The hydrogen bond analysis suggests that the H3O+ ions show a preference for the carbonyl oxygen in the sn-2 chain (O2A) followed by the non-ester phosphate oxygen (OP1). There are little to no hydrogen bonds formed by H3O+ and the other non-ester phosphate oxygen (OP2) and the sn-1 carbonyl oxygen (O1A). The three hydrogen bonding sites on H3O+ combined with this preference means that the majority of the time H3O+ is surrounded by three lipid molecules forming hydrogen bonds with either three O2A or two O2A and one OP1 oxygens, as depicted in Fig. 6. This is consistent with the observations reported previously in which the pH dependence of POPC which contains ester carbonyl groups resulted in a 245-fold change in conduction between pH 5 and pH 9 of the bulk solution compared to an 8-fold change in conduction for an ether-lipid containing no carbonyls.45 Previous simulations of hydrated protons with DOPC bilayers suggested that in the deep interface region, closer to the hydrophobic core, a distorted Zundel cation (H5O2+) can bridge between the phosphate and carbonyl oxygen.51 Note that the results from this study cannot be directly compared to data from this study as in our simulations a hydronium cation was used where the excess proton is localised on a single water molecule rather than a Zundel cation.

image file: c7cp06776c-f6.tif
Fig. 6 The two most common hydrogen bonding networks formed by H3O+ ions and POPC lipids observed in simulations of a POPC lipid bilayer in the presence of 0.4 M [H3O+] and 0.04 M [H3O+]. (A) In the most populated arrangement a H3O+ ion forms hydrogen bonds with the sn-2 carbonyl oxygens (O2A) from two lipids and the non-ester oxygen (OP1) from a third lipid. (B) The second most common arrangement in which H3O+ ion forms hydrogen bonds with the sn-2 carbonyl oxygens (O2A) from three lipids.

Combined, these results clearly indicate that H3O+ ions form a distinct hydrogen bonding network with the polar lipid oxygen that is different from the one formed by water molecules. Even at low pH the concentration of water (55 M) and Na+ ions (0.15 M) exceeds the concentration of H3O+ by at least 3 orders of magnitude, yet it is the interactions of H3O+ with the lipid head groups that alter the structure of the membrane. The small size and positive charge of H3O+ results in a high-field strength providing it with the ability to form strong non-covalent interactions with the polar oxygen in the lipid molecules. In our simulations, this is evident in the shorter and longer-lived hydrogen bonds of H3O+ when compared to the ones by H2O. The lifetimes of the hydrogen bonds and the relative occupancies for the different lipid arrangements is likely distorted by the absence of ability for proton transfer. To obtain accurate lifetimes for hydrogen bonds and relative stabilities of the different bonding arrangement would require the use of a model that allows for proton hopping and charge transfer delocalization, which is beyond the scope of this study. Nevertheless, on a qualitative level the longer-lived H-bonds observed in our simulations are in agreement with significantly reduced proton transfer rates at membrane surface reported in previous simulations.50 The increased stability of hydrogen bonds formed by the H3O+ ion can also be related to the diffusion of protons and hydronium ions across membrane surfaces. For example, based on the energetics and mobility of protons near a solvated DMPC bilayer Wolf et al.50 observed anomalous diffusion “characterized by a short subdiffusive regime (1 ns) and a subsequent superdiffusive regime” where the subdiffusive regime results from the H3O+ ion being tightly bound to lipids. Similarly, Yamashita and Voth51 reported that the “lateral diffusion coefficients are 1 order of magnitude smaller in the interface region than in bulk water for all the lipids.”

Overall, the simulations are consistent with the model of lipid bilayers and biological membranes in a low-pH environment that has emerged from a range of experimental studies. Specifically, lipid bilayers and membrane undergo a reduction in the decrease in area per lipid that stems from an increase in the stability of the hydrogen bonds between lipids and hydronium ions H3O+. At high H3O+ concentration (that is acidic pH), H3O+ ions accumulate at the membrane interface where they form strong hydrogen bond networks with the carbonyl and phosphate oxygen. This results in a condensation of the lipids onto the H3O+ ion resulting in a reduced area per lipid AL.

5 Conclusions

To the best of our knowledge this is the first study that reports, at the molecular level, the interactions between H3O+ and phospholipids and its effect on the structure of phospholipid bilayers. Our results indicate that the H3O+ ions partition into the water–lipid interface where they displace water and form strong and long-lived hydrogen bonds with the phospholipid phosphate and carbonyl oxygens. These strong interactions result in a reduced area per lipid and an increased bilayer thickness. In order to quantitatively predict this effect, more concentrations would need to be examined but this was beyond the scope of the present study. Additional studies will also be required to investigate this effect in more complex membranes.

Conflicts of interest

The authors have no conflicts of interest to declare.


This work was supported by resources provided of the Pawsey Supercomputing Centre with funding from the Australian Government and the Government of Western Australia, as well as with the assistance of resources provided at the NCI National Facility systems at the Australian National University through the National Computational Merit Allocation Scheme supported by the Australian Government. This work was supported in part by Australian National Health and Medical Research Council Early Career Fellowship to ED.


  1. G. van Meer, D. R. Voelker and G. W. Feigenson, Nat. Rev. Mol. Cell Biol., 2008, 9, 112–124 CrossRef CAS PubMed .
  2. G. Pabst, N. Kučerka, N. Mu-Ping and J. Katsaras, Liposomes, Lipid Bilayers and Model Membranes: From Basic Research to Application, CRC Press, 2016 Search PubMed .
  3. R. Basu, S. De, D. Ghosh and P. Nandy, Physica A, 2001, 292, 146–152 CrossRef CAS .
  4. A. Blicher, K. Wodzinska, M. Fidorra, M. Winterhalter and T. Heimburg, Biophys. J., 2009, 96, 4581–4591 CrossRef CAS PubMed .
  5. S. L. Veatch and S. L. Keller, Biochim. Biophys. Acta, Mol. Cell Res., 2005, 1746, 172–185 CrossRef CAS PubMed .
  6. R. Chen, D. Poger and A. E. Mark, J. Phys. Chem. B, 2011, 115, 1038–1044 CrossRef CAS PubMed .
  7. R. Winter, Curr. Opin. Colloid Interface Sci., 2001, 6, 303–312 CrossRef CAS .
  8. A. Aroti, E. Leontidis, M. Dubois and T. Zemb, Biophys. J., 2007, 93, 1580–1590 CrossRef CAS PubMed .
  9. M. L. Berkowitz and R. Vácha, Acc. Chem. Res., 2012, 45, 74–82 CrossRef CAS PubMed .
  10. H. I. Petrache, T. Zemb, L. Belloni and V. A. Parsegian, Proc. Natl. Acad. Sci. U. S. A., 2006, 103, 7982–7987 CrossRef CAS PubMed .
  11. J. Song, J. Franck, P. Pincus, M. W. Kim and S. Han, J. Am. Chem. Soc., 2014, 136, 2642–2649 CrossRef CAS PubMed .
  12. A. C. Carreira, R. F. M. de Almeida and L. C. Silva, Sci. Rep., 2017, 7, 3949 CrossRef PubMed .
  13. D. A. Redfern and A. Gericke, J. Lipid Res., 2005, 46, 504–515 CrossRef CAS PubMed .
  14. S. D. Krämer, A. Braun, C. Jakits-Deiser and H. Wunderli-Allenspach, Pharm. Res., 1998, 15, 739–744 CrossRef .
  15. K.-J. Schaper, H. Zhang and O. A. Raevsky, Quant. Struct.-Act. Relat., 2001, 20, 46–54 CrossRef CAS .
  16. J. A. Kellum, M. Song and J. Li, Crit. Care, 2004, 8, 331 CrossRef PubMed .
  17. L. E. Gerweck and K. Seetharaman, Cancer Res., 1996, 56, 1194–1198 CAS .
  18. Y. Kato, S. Ozawa, C. Miyamoto, Y. Maehata, A. Suzuki, T. Maeda and Y. Baba, Cancer Cell Int., 2013, 13, 89 CrossRef CAS PubMed .
  19. Y. Huang and J. O. McNamara, Cell, 2004, 118, 665–666 CrossRef CAS PubMed .
  20. A. E. Ziemann, M. K. Schnizler, G. W. Albert, M. A. Severson, M. A. Howard Iii, M. J. Welsh and J. A. Wemmie, Nat. Neurosci., 2008, 11, 816–822 CrossRef CAS PubMed .
  21. M. de Rosa, A. Gambacorta, B. Nicolaus and W. D. Grant, J. Microbiol., 1983, 129, 2333–2337 CrossRef CAS .
  22. H. N. Khaleque, J. P. Ramsay, R. J. Murphy, A. H. Kaksonen, N. J. Boxall and E. L. Watkin, Genome Announc., 2017, 5, e01469–e01416 Search PubMed .
  23. R. Koynova and M. Caffrey, Biochim. Biophys. Acta, Biomembr., 1998, 1376, 91–145 CrossRef CAS .
  24. J. M. Boggs, G. Rangaraj and K. M. Koshy, Chem. Phys. Lipids, 1986, 40, 23–34 CrossRef CAS .
  25. G. Cevc, J. Biochem., 1987, 26, 6305–6310 CrossRef CAS .
  26. G. Cevc, J. Phys., 1989, 50, 1117–1134 CAS .
  27. G. Cevc, A. Watts and D. Marsh, J. Biochem., 1981, 20, 4955–4965 CrossRef CAS .
  28. H. Eibl and P. Woolley, Biophys. Chem., 1979, 10, 261–271 CrossRef CAS PubMed .
  29. H. Träuble and H. Eibl, Proc. Natl. Acad. Sci. U. S. A., 1974, 71, 214–219 CrossRef .
  30. Y. Zhou and R. M. Raphael, Biophys. J., 2007, 92, 2451–2462 CrossRef CAS PubMed .
  31. A. D. Petelska and Z. A. Figaszewski, Biophys. J., 2000, 78, 812–817 CrossRef CAS PubMed .
  32. A. D. Petelska and Z. A. Figaszewski, Biochim. Biophys. Acta, Biomembr., 2002, 1561, 135–146 CrossRef CAS .
  33. M. Brändén, T. Sandén, P. Brzezinski and J. Widengren, Proc. Natl. Acad. Sci. U. S. A., 2006, 103, 19766–19770 CrossRef PubMed .
  34. B. Gupta, O. Haas and G. G. Scherer, J. Appl. Polym. Sci., 1994, 54, 469–476 CrossRef CAS .
  35. J. Heberle, J. Riesle, G. Thiedemann, D. Oesterhelt and N. A. Dencher, Nature, 1994, 370, 379–382 CrossRef CAS PubMed .
  36. S. J. Marrink, F. Jähnig and H. J. Berendsen, Biophys. J., 1996, 71, 632–647 CrossRef CAS PubMed .
  37. J. F. Nagle and R. A. Dilley, J. Bioenerg. Biomembr., 1986, 18, 55–64 CrossRef CAS PubMed .
  38. H. L. Tepper and G. A. Voth, Biophys. J., 2005, 88, 3095–3108 CrossRef CAS PubMed .
  39. J. Fitter, R. E. Lechner and N. A. Dencher, J. Phys. Chem. B, 1999, 103, 8036–8050 CrossRef CAS .
  40. F. Foglia, M. J. Lawrence, C. D. Lorenz and S. E. McLain, J. Chem. Phys., 2010, 133, 145103 CrossRef PubMed .
  41. M. Pasenkiewicz-Gierula, Y. Takaoka, H. Miyagawa, K. Kitamura and A. Kusumi, J. Phys. Chem. A, 1997, 101, 3677–3691 CrossRef CAS .
  42. M. Pasenkiewicz-Gierula, Y. Takaoka, H. Miyagawa, K. Kitamura and A. Kusumi, Biophys. J., 1999, 76, 1228–1240 CrossRef CAS PubMed .
  43. R. H. Pearson and I. Pascher, Nature, 1979, 281, 499–501 CrossRef CAS PubMed .
  44. W. Zhao, D. E. Moilanen, E. E. Fenn and M. D. Fayer, J. Am. Chem. Soc., 2008, 130, 13927–13937 CrossRef CAS PubMed .
  45. C. G. Cranfield, T. Berry, S. A. Holt, K. R. Hossain, A. P. Le Brun, S. Carne, H. Al Khamici, H. Coster, S. M. Valenzuela and B. Cornell, Langmuir, 2016, 32, 10725–10734 CrossRef CAS PubMed .
  46. O. A. Gopta, D. A. Cherepanov, W. Junge and A. Y. Mulkidjanian, Proc. Natl. Acad. Sci. U. S. A., 1999, 96, 13159–13164 CrossRef CAS .
  47. C. J. Slevin and P. R. Unwin, J. Am. Chem. Soc., 2000, 122, 2597–2602 CrossRef CAS .
  48. A. Mashaghi, P. Partovi-Azar, T. Jadidi, M. Anvari, S. P. Jand, N. Nafari, M. R. R. Tabar, P. Maass, H. J. Bakker and M. Bonn, J. Phys. Chem. B, 2013, 117, 510–514 CrossRef CAS PubMed .
  49. A. M. Smondyrev and G. A. Voth, Biophys. J., 2003, 85, 864–875 CrossRef .
  50. Maarten G. Wolf, H. Grubmüller and G. Groenhof, Biophys. J., 2014, 107, 76–87 CrossRef CAS PubMed .
  51. T. Yamashita and G. A. Voth, J. Phys. Chem. B, 2010, 114, 592–603 CrossRef CAS PubMed .
  52. D. Poger and A. E. Mark, J. Chem. Theory Comput., 2012, 8, 4807–4817 CrossRef CAS PubMed .
  53. D. Poger, W. F. van Gunsteren and A. E. Mark, J. Comput. Chem., 2010, 31, 1117–1125 CrossRef CAS PubMed .
  54. B. Hess, C. Kutzner, D. van Der Spoel and E. Lindahl, J. Chem. Theory Comput., 2008, 4, 435–447 CrossRef CAS PubMed .
  55. N. Schmid, A. P. Eichenberger, A. Choutko, S. Riniker, M. Winger, A. E. Mark and W. F. van Gunsteren, Eur. Biophys. J., 2011, 40, 843–856 CrossRef CAS PubMed .
  56. D. Poger, B. Caron and A. E. Mark, Biochim. Biophys. Acta, Biomembr., 2016, 1858, 1556–1565 CrossRef CAS PubMed .
  57. D. Poger and A. E. Mark, J. Chem. Theory Comput., 2010, 6, 325–336 CrossRef CAS PubMed .
  58. K. B. Koziara, M. Stroet, A. K. Malde and A. E. Mark, J. Comput.-Aided Mol. Des., 2014, 28, 221–233 CrossRef CAS PubMed .
  59. A. K. Malde, L. Zuo, M. Breeze, M. Stroet, D. Poger, P. C. Nair, C. Oostenbrink and A. E. Mark, J. Chem. Theory Comput., 2011, 7, 4026–4037 CrossRef CAS PubMed .
  60. S. Canzar, M. El-Kebir, R. Pool, K. Elbassioni, A. K. Malde, A. E. Mark, D. P. Geerke, L. Stougie and G. W. Klau, J. Comput. Biol., 2013, 20, 188–198 CrossRef CAS PubMed .
  61. D. J. Bonthuis, S. I. Mamatkulov and R. R. Netz, J. Chem. Phys., 2016, 144, 104503 CrossRef PubMed .
  62. A. A. Chialvo, P. T. Cummings and J. M. Simonson, J. Chem. Phys., 2000, 113, 8093–8100 CrossRef CAS .
  63. L. X. Dang, J. Chem. Phys., 2003, 119, 6351–6353 CrossRef CAS .
  64. B. J. Gertner and J. T. Hynes, Faraday Discuss., 1998, 110, 301–322 RSC .
  65. I. Kusaka, Z. G. Wang and J. H. Seinfeld, J. Chem. Phys., 1998, 108, 6829–6848 CrossRef CAS .
  66. S. Urata, J. Irisawa, A. Takada, W. Shinoda, S. Tsuzuki and M. Mikami, J. Phys. Chem. B, 2005, 109, 4269–4278 CrossRef CAS PubMed .
  67. R. Vácha, V. Buch, A. Milet, J. P. Devlin and P. Jungwirth, Phys. Chem. Chem. Phys., 2007, 9, 4736–4747 RSC .
  68. N. Agmon, Chem. Phys. Lett., 1995, 244, 456–462 CrossRef CAS .
  69. C. J. T. Grotthuss, Ann. Chim. Phys., 1806, 58, 54–74 Search PubMed .
  70. J. M. J. Swanson, C. M. Maupin, H. Chen, M. K. Petersen, J. Xu, Y. Wu and G. A. Voth, J. Phys. Chem. B, 2007, 111, 4300–4314 CrossRef CAS PubMed .
  71. H. J. C. Berendsen, J. P. M. Postma, W. F. van Gunsteren and J. Hermans, in Inter molecular Forces, ed. B. Pullman, Springer, Netherlands, 1981, ch. 21, vol. 14, pp. 331–342 Search PubMed .
  72. H. J. C. Berendsen, J. P. M. Postma, W. F. van Gunsteren, A. Dinola and J. R. Haak, J. Chem. Phys., 1984, 81, 3684–3690 CrossRef CAS .
  73. B. Hess, H. Bekker, H. J. C. Berendsen and J. G. E. M. Fraaije, J. Comput. Chem., 1997, 18, 1463–1472 CrossRef CAS .
  74. S. Reißer, D. Poger, M. Stroet and A. E. Mark, J. Chem. Theory Comput., 2017, 13, 2367–2372 CrossRef PubMed .
  75. I. G. Tironi, R. Sperb, P. E. Smith and W. F. van Gunsteren, J. Chem. Phys., 1995, 102, 5451–5459 CrossRef CAS .
  76. E. London and G. W. Feigenson, J. Lipid Res., 1979, 20, 408–412 CAS .
  77. M. S. Fernández and E. Calderón, Berichte der Bunsengesellschaft für physikalische Chemie, 1991, 95, 1669–1674 CrossRef .
  78. M. Gutman and E. Nachliel, Biochim. Biophys. Acta, Bioenerg., 1990, 1015, 391–414 CrossRef CAS .
  79. M. R. Moncelli, L. Becucci and R. Guidelli, Biophys. J., 1994, 66, 1969–1980 CrossRef CAS PubMed .
  80. D. P. Siegel, J. L. Burns, M. H. Chestnut and Y. Talmon, Biophys. J., 1989, 56, 161–169 CrossRef CAS PubMed .
  81. N. Michaud-Agrawal, E. J. Denning, T. B. Woolf and O. Beckstein, J. Comput. Chem., 2011, 32, 2319–2327 CrossRef CAS PubMed .
  82. D. E. Sagnella and G. A. Voth, Biophys. J., 1996, 70, 2043–2051 CrossRef CAS PubMed .
  83. U. W. Schmitt and G. A. Voth, J. Chem. Phys., 1999, 111, 9361–9381 CrossRef CAS .
  84. N. Kučerka, M.-P. Nieh and J. Katsaras, Biochim. Biophys. Acta, Biomembr., 2011, 1808, 2761–2771 CrossRef PubMed .
  85. N. Kučerka, S. Tristram-Nagle and J. F. Nagle, J. Membr. Biol., 2006, 208, 193–202 CrossRef PubMed .
  86. C. G. Cranfield, S. T. Henriques, B. Martinac, P. Duckworth, D. J. Craik and B. Cornell, Langmuir, 2017, 33, 6630–6637 CrossRef CAS PubMed .
  87. S. Nizalapur, K. K. K. Ho, O. Kimyon, E. Yee, T. Berry, M. Manefield, C. G. Cranfield, M. Willcox, D. S. Black and N. Kumar, Org. Biomol. Chem., 2016, 14, 3623–3637 CAS .


Electronic supplementary information (ESI) available. See DOI: 10.1039/c7cp06776c

This journal is © the Owner Societies 2018