Bayesian determination of the effect of a deep eutectic solvent on the structure of lipid monolayers

In this work, we present the first example of the self-assembly of phospholipid monolayers at the interface between air and an ionic solvent. Deep eutectic solvents are a novel class of environmentally friendly non-aqueous room temperature liquids with tunable properties, that have wide ranging potential applications and are capable of promoting the self-assembly of surfactant molecules. We use a chemically-consistent Bayesian modelling of X-ray and neutron reflectometry measurements to show that these monolayers broadly behave as they do on water. This method allows for the monolayer structure to be determined, alongside the molecular volumes of the individual monolayer components without the need for water-specific constraints to be introduced. Furthermore, using this method we are able to better understand the correlations present between parameters in the analytical model. This example of a non-aqueous phospholipid monolayer has important implications for the potential uses of these solvents and for our understanding of how biomolecules behave in the absence of water.

have been widely studied as simplistic models for biological membranes. As such, they have been used to gain insight into many biological processes that are technologically and medically relevant. For example, investigations at the air/salt-water interface have identified the importance that interactions between charged phospholipid heads and ions present in solution have on the structure, monomer packing and stability of the monolayer [19,20]. However, the native environment for lipids in-vivo is far from a simple aqueous solution. In fact, it has been suggested [2,4] that DES might form within the crowded cellular environment and could assist in solubilizing biological species in an intermediate environment between that of the hydrophobic phospholipid tails and highly polar water-rich regions, thereby assisting survival under extreme conditions such as freezing temperatures or drought where the water content of cells is restricted. Despite the broad interest in these systems, the presence of stable phospholipid monolayers at the interface between air and an ionic solvent media has not been previously reported, to the best of the authors' knowledge.
The study of lipid monolayers on a non-aqueous solvent like DES is relevant to some potential new scientific areas or technological applications. For example, it may be possible to study the interactions of proteins with model membranes in a zero or low water environment, or even at low temperatures where aqueous based model mem-branes are not practical. To enable such experiments we must first understand how lipid monolayers are affected by a non-aqueous environment. Published studies covering such non-aqueous systems are, however, rather scarce mirroring the availability of solvents that show a solvophobic effect comparable to water. In the absence of water, formamide has been reported to support the formation of lipid monolayers [21,22]. Monolayers have also been observed at the air-mercury interface [23,24], although in this case, the similarity to water is somewhat limited. We have recently shown that self-assembly of surfactant monolayers is possible for DES with behaviour that is much more comparable to that seen in water [11,25].
In this work, we have investigated the structure of phospholipid monolayers at the air-DES interface. To do this we have used the combination of both X-ray and neutron reflectometry (XRR & NR), techniques that allow detailed structures to be determined with near-atomic resolution. In order to interpret our results, we have used a novel analysis method which allows us to sufficiently constrain our model while using Bayesian inference [26,27] to properly quantify the errors and correlations inherent in such models of reflectometry data. Specifically, we have constrained the model to ensure that it is chemically self-consistent across multiple measurements at various surface pressures. Four different phospholipids; 1,2-dipalmitoyl-sn-glycero-3-phosphocholine (DPPC), 1,2-dimyristoyl-sn-glycero-3-phosphocholine (DMPC), 1,2-dilauroyl-sn-glycero-3phosphocholine (DLPC) and 1,2-dimyristoyl-sn-glycero-3-phospho-(1'-rac-glycerol) (sodium salt) (DMPG), were studied at the interface between a 1:2 mixture of choline chloride:glycerol and air. This has allowed a comparison of two chemically distinct phospholipid head components, the zwitterionic phosphocholine (PC) and the anionic phosphoglycerol (PG) to be probed in this ionic solvent ( Figure 1). In addition, the effect of the tail chain length in the assembly of the Langmuir monolayer was explored by comparing molecules with progressively increasing tail length (12, 14 and 16 carbon atoms in each of the two hydrocarbon chains in each molecule). Our study is inevitably limited in extent due to the limited availability of X-ray and neutron beamtime together with the prohibitive expense of deuterated chemicals. However, we believe that this work both hints at interesting differences in behaviour for phospholipids monolayers on DES and demonstrates how our approach to fitting can provide a quantitative understanding of this behaviour.
Recent developments in computational resources and software have enabled powerful methodologies and algorithms to be harnessed by those from non-expert backgrounds. This has benefitted significantly from opensource software projects such as the Python language [28] and the Jupyter notebooks framework [29]. In the area of NR and XRR, the landscape of data-analysis software is diverse, with a range of software packages available from a variety of sources; refnx [30,31], MOTOFIT [32], Rascal [33] Aurore [34], Refl1D [35], and GenX [36].
The use of a Python library, such as refnx, enables the implementation of custom models that contain chemically-relevant information as well as the application of probability distribution function (PDF) sampling techniques. The Python library emcee [37] allows refnx to access the Goodman & Weare Affine Invariant Markov chain Monte Carlo (MCMC) Ensemble method [38]. This allows the sampling of the high-dimensionality parameter space, relevant in reflectometry analysis, in a Bayesian fashion, where the new samples are generated with consideration of those sampled previously [39]. Bayesian inference gives an understanding of the PDF for the fitted parameters and therefore estimations of their inverse uncertainties and inter-parameter correlations.

A. Materials
Choline chloride (99 % %, Sigma-Aldrich) and glycerol (99 %, Sigma-Aldrich), and d 9 -choline chloride (99 %, 98 % D, CK Isotopes) and d 8 -glycerol (99 %, 98 % D, CK Isotopes) were purchased and used without further purification. The DES was prepared by mixing the precursors at a 1:2 ratio of choline chloride:glycerol, and heating at 80 • C until a homogeneous, transparent liquid formed [1]. The solvent was equilibrated overnight at 40 • C and subsequently stored under a dry atmosphere. Due to the limited availability of the deuterated precursors, a fully protonated subphase (hDES) and a partially deuterated subphase (hdDES) were prepared and used during the neutron reflectometry (NR) experiment. The partially deuterated subphase was prepared using the following mixtures of precursors: 1 mol of 0.38 mol fraction of hcholine chloride/0.62 mol fraction of d-choline chloride; and 2 mol of 0.56 mol fraction of h-glycerol/0.44 mol frac-tion of d-glycerol. The deuterated solvent was prepared following the procedure discussed above.
The water content of the DES was determined before and after each experiment by Karl-Fischer titration (Mettler Toledo DL32 Karl-Fischer Coulometer, Aqualine Electrolyte A, Aqualine Catholyte CG A) in order to ensure water presence was kept to a minimum. Those measurements showed that the water content of the solvent was kept below 0.3 wt/% during all the experimental procedures presented here, which we assume to be negligible and have to little impact on the characteristics of the DES [3,4].
In the XRR experiment, sample preparation was performed in situ using the standard method for the spreading of insoluble monolayers on water: a certain amount of the phospholipid solution was spread onto the liquid surface in order to provide a given surface concentration. After the evaporation of the chloroform, it is assumed that the resulting system is a solvent subphase with a monolayer of phospholipid at the interface. Surface concentration was modified by closing and opening the PTFE barriers of a Langmuir trough. In order to minimise the volumes used in the NR experiment (to keep the cost of deuterated compounds to a manageable level), it was not possible to use a Langmuir trough. Instead, small Delrin adsorption troughs were used that did not have controllable barriers. So, although the surface coverage was nominally the same as used in the X-ray studies, the lack of precise control over the surface pressure meant that it was not appropriate to co-refine XRR and NR contrasts together.

B. Methods
XRR measurements were taken on I07 at Diamond Light Source, at 12.5 keV photon energy using the double-crystal deflector [40]. The reflected intensity was measured in a momentum transfer range from 0.018Å −1 to 0.7Å −1 . The data were normalised with respect to the incident beam and the background was measured from off-specular reflection and subsequently subtracted. Samples were equilibrated for at least one hour and preserved under an helium atmosphere to minimise the adsorption of water by the subphase. XRR data were The NR experiments were performed on FIGARO at the Institut Laue-Langevin using the time-of-flight method [41]. Data at two incident angles of 0.62 • and 3.8 • were measured to provide a momentum transfer range from 0.005Å −1 to 0.18Å −1 . Two surface pressures for each system and contrast was measured (DMPC: 20 and 25 mN m −1 , DPPC: 15 and 20 mN m −1 ). Similar to the X-ray procedure, samples were given enough time to equilibrate (at least two hours), kept under an inert atmosphere, and all measurements were conducted at 22 • C.

C. Data analysis
The use of XRR and NR to analyse the structure of phospholipids on the surface of water is well documented [19,20,[42][43][44][45][46]. The models used in the rationalisation of XRR and NR data have varied significantly in numbers of layers present, use of interfacial roughness, and the parameterisation of the physical constraints applied. Frequently, these physical constraints include the volumes of the phospholipid head and tail components, using values taken from other techniques, such as those shown in Table I. Additionally, a recent evaluation of the applicability of different models for surfactant and phospholipid monolayers from the NR perspective has been published [47], that suggests possible oversights in the modelling of NR data.
In Table I, there appears to be a general consensus that the component volume for the phosphocholine (PC) head is around 320Å 3 to 360Å 3 , while the phosphatidylglycerol (PG) head is in the range 289Å 3 to 291Å 3 . However, it is not clear that the head component volumes from the literature, that are derived from water-based measurements, will be appropriate for this work, which involves a non-aqueous solvent. The charged nature of the zwitterionic or anionic lipid heads may have different interactions with the polar, but neutral water as compared to the charged DES components [25], which will affect the phospholipid head component volume. Therefore, herein we apply a chemically-consistent model that allows for the co-refinement of reflectometry measurements at different surface pressure and makes no assumption of the component volume for the lipid head, V h , or tail, V t . Instead, these parameters were allowed to vary for each lipid while being constrained to be self-consistent over different surface pressures in the same phase; Liquid-Condensed (LC) for DPPC and Liquid-Expanded (LE) for DMPC, DMPG, and DLPC. Furthermore, it is known that, on water, increased surface pressure and the associated LE-LC phase transitions lead to a compression of the lipid tail volume [54,55], and this compaction has not necessarily been accounted for in the literature [47]. This model avoids this issue by making no assumption about the molecular volumes and only considers surface pressures that we believe to be in the same phase.
Our chemically-consistent model has been implemented in the Python library refnx [30,31]. This software allows for the inclusion of a custom model to be defined, from which parameters feed into the Abelès reflectivity model (a model that is widely used to calculate reflectivity [26,27]). This custom model, along with a series of Jupyter notebooks showing, in full, the analysis performed, can be found in the ESI and is available under a CC BY-SA 4.0 licence.
This model consists of two layers to define the lipid monolayer; the head layer at the interface with the solvent and the tail layer at the interface with the air. The head components have a calculated scattering length, b h , (found as a summation of the X-ray or neutron atomic scattering lengths), and a component volume, V h . These head components make up a layer with a given thickness, d h , and roughness, σ h , within which some volume fraction of solvent can intercalate, φ h . The tail layer is defined in the same fashion, except that the thickness, d t , is limited such that it may be no greater than the maximum extended length of the lipid tail (the Tanford length, t t [56]), which is given in Table II, and that no solvent may intercalate into the layer (e.g. φ t = 0). The scattering length density (SLD) of the tail and head layers used in the Abelès model can therefore be found as follows, where, SLD s is the scattering length density of the subphase (DES), and i indicates either the tail or head layer.
To ensure that the number density of head components and pairs of tail components is the same, the following constraint was included in the model [57], Based on the work of Campbell et al. [47], a single value for the interfacial roughness was fitted for all of the interfaces, including the subphase (i.e. σ h = σ t = σ s ), as there is only a single lipid molecule type in each monolayer. Therefore, any capillary wave roughness at the air-DES interface is carried equally through the layers. The interfacial roughness was constrained to be greater than 3.3Å in agreement with previous work [10].
In order to justify the use of a single tail volume across many surface pressures, it was necessary to ensure that the lipids remain in the same phase. On water, this can be demonstrated with a Langmuir isotherm. However, while we have confidence that the individual surface pressures measured were reliable, we were unable to collect consistent Langmuir isotherm measurements, due to the high viscosity of the DES. Instead, we have used grazing incidence X-ray diffraction to confirm the phases of DMPC and DPPC at 30 mN m −1 . DPPC was found to be in the LC phase and DMPC in the LE phase at room temperature for the surface pressures measured (see Section S3 in the ESI). We assume that DMPG and DLPC are also in the LE phase since there is no reason to believe that the phase behaviour in these systems differs significantly from DMPC at the same temperature.
In the first of two steps, this custom model was used to co-refine the component volume of the lipid head component, V h , the volume of the tail component, V t , and the head thickness, d h across XRR measurements at four different surface concentrations. In keeping with the work of Campbell et al. [47], a single value for the head thickness was fitted for each lipid across all surface pressures, as the thickness of the head layer was considered to be dependent on molecular dimensions only, and has been shown to vary little with surface pressure [20]. The following parameters were allowed to vary; d t , and σ t,h,s , independently across the surface pressures, while others, shown in Table II, were held constant at the values given. For each co-refinement of four XRR measurements, there were, in total, eleven degrees of freedom in the fitting process. Throughout all of the analyses, the reflectometry scale factor was allowed to vary freely, while the background was constrained to the intensity of either the largest or second-largest q-value.
In the second step, the head and tail component volumes, and head layer thickness determined from XRR were fixed for the refinement of the custom model against the NR measurements. This approach means that the number of variable parameters to fit the NR data can be reduced to two, namely the thickness of the tail layer, d t , and the interfacial roughness, σ t,h,s , for the co-refinement of two datasets. Table II also gives the details of the scattering lengths and SLDs used as invariant parameters for the NR fitting.
In both cases, the refinement of the custom model to the experimental data involved the transformation of the reflectometry calculated from the model and the data into Rq 4 such that the contribution of the Fresnel decay was removed, before using the differential evolution method available to refnx from the scipy library [58], to find the parameters that gave the best fit to the data. The parameter space was then probed using the MCMC method available through emcee [37], which allowed for an estimate of the probability distribution function (PDF) associated with each parameter. In the MCMC sampling, 200 walkers were used over 1000 iterations, following equilibration of 200 iterations. The use of MCMC sampling allowed for a Bayesian inference of the PDF for each of the variables and their respective interactions and the Shapiro test to be used to assess if each PDF was normally distributed. Parameters that were shown to be normally distributed are given with symmetric confidence intervals, while those that failed the Shapiro test are given with asymmetric confidence intervals (95 % confidence intervals in both cases). However, it is important to note that these are not true confidence intervals, and account only for the uncertainty present in the data, i.e. they do not account of systematic uncertainty in the measurement that is underrepresented, or unrepresented, in the experimental dataset.

III. RESULTS & DISCUSSIONS
The chemically-consistent model was co-refined across the four surface pressure XRR measurements for each lipid. Although lipid monolayers on water often reach values above 50 mN m −1 , signifying a decrease in surface tension to ca. 20 mN m −1 upon addition of lipid, it was not possible to reach those values for the monolayers on DES [19]. The surface tension of pure choline chloride:glycerol has been previously reported to be 63.5(5) mN m −1 [10]. This means that a surface pressure of ca. 40 mN m −1 iin choline chloride:glycerol will also result in a reduction to ca. 20 mN m −1 in the surface tension. Interestingly, the absolute values of surface tension in both water and choline chloride:glycerol after addition of lipid are similar. Therefore, surface pressures of up to 40 mN m −1 were measured for the monolayers on DES, when possible. The resulting XRR profiles and associated SLD profiles are shown in Figure 2. Table  III gives details of all varied parameters for each lipid at 30 mN m −1 , as well as the details of φ h which was determined from Eqn. 2 (the same details for the other surface pressures are available in Section S1 of the ESI).  Following the initial structural determination of the monolayer by XRR, NR was used to confirm the structure and show the applicability of the chemically-consistent model for DPPC and DMPC. The resulting NR profiles and associated SLD profiles, at a surface pressure of 20 mNm −1 are given in Figure 3 (the other NR and SLD profiles can be found in Figure S1 of the ESI). Table  IV gives details of the varied parameters at each surface pressure as well as φ h as determined from NR.

A. Effect of compression on monolayer thickness
From Tables III and IV, we can see that, as expected and as found in previous work [19,59], the thickness of the tail layer increases as the number of carbon atoms in the tail chain increases. Furthermore, the thickness of the tail layers in these monolayers appears to agree well with values found for water-analogues; 13.72 ± 0.01 A at 30 mN/m in DES compared with d t = 15.8Å at 30 mN m −1 in water for DMPC [43], and 16.91 ± 0.01 A at 30 mN m −1 in DES compared with d t = 16.7Å at 40 mN m −1 in water for DPPC [45].
The variation of the tail layer thickness in the models with surface pressure is given for each lipid in Figure  4(e). As is commonly observed for lipid monolayers at the air-water interface, an increase in surface pressure results in an increase in the thickness of the tail layer. The variation of this parameter differs for each of the lipids. For DLPC and DMPC a gradual increase in the thickness is observed at the surface pressures measured here. For DMPC and DPPC however, after a gradual increase in thickness between 15 mN m −1 and 25 mN m −1 , a plateau is reached. This is indicative of the formation of a compacted monolayer where further compression does not produce further changes in the thickness of the tail layer. This behaviour has been seen for DMPC and DPPC monolayers, where the formation of the condensed phases is reached at high surface concentrations [19].

B. Effect of compression on solvent concentration
In Figure 4(f), it is clear that for all four lipids, as the surface pressure is increased there is a corresponding decrease in the percentage solvent present in the lipid head layer. This can be rationalised by considering that when the surface pressure is increased, the free volume available to the solvent between the lipid head components reduces forcing the solvent out of the lipid head layer and into the bulk. A similar effect has been observed when increasing the surface pressure from 11 mN m −1 to 31 mN m −1 for a DMPC/DMPG monolayer at the airwater interface [42].

C. Effect of compression on the lipid tail component volumes
It can be seen by comparing Tables I and III that the volumes of the lipid tails are significantly lower in the current measurements than found previously, by other techniques. It is unlikely that this is a result of the DES subphase, due to the hydrophobic nature of the lipid tails. However, this reduction has been shown previously [47], where it was rationalised by the compaction of the monolayer at elevated surface pressure. In that  In this work, a single tail component volume was fitted to each lipid for all four surface pressures that were measured. This is based on the assumption, that at all four surface pressures, the lipids adopt the same phase and therefore any variation in the structure with surface pressure would manifest only as a change in the tail thickness. It is clear when comparing Tables I and III that some of the tail component volumes are also reduced in the current XRR measurements compared to those determined previously. The reduction was found to be between 8 % to 12 % for DPPC, DMPC and DLPC when compared with literature sources at 24 • C to 30 • C, this is in good agreement with the maximum compression percentage of 15 % noted by Small and coworkers [55]. DMPG shows a small increase in the tail volume relative to the literature value quoted at a lower temperature. Notably, this value is similar to that found in this work for DMPC, which has the same tail structure and suggests that our results are at least self-consistent.   Table I. Interestingly, the component volume for the PG head is similar to that for the PC head with a value of 329.95 +0.32 −0.33Å 3 . The PG head component volume in water, from either DMPG using differential vibrating tube densimetry [52] or POPG using molecular dynamics simulations [53], is noticeably smaller. This indicates that there may be some effect arising from the solvation in choline chloride:glycerol causing an apparent increase in the PG component volume when compared with water. However, this has only been shown for a single PG-lipid at the airDES interface.
The major difference between the two head groups of the lipids is that the PG is present as a sodium salt, whereas the PC is zwitterionic. When in solution the anionic PG head is expected to associate with cations in solution, as it does in water [60] where such interactions depended on a variety of factors including the ionic strength. In the case of a DES, the environment is inherently ionic and therefore the interaction of an anionic lipid head may be more complex. As well as interacting with the sodium, the head is likely to interact with the choline cations, similar to behaviour reported previously for surfactant micelles [61]. The extent of interaction with each of the cations is unclear, but regardless it seems likely that the solvation of the PG head is improved in the DES relative to water. This better solvation would explain the apparent increase in the volume of the PG head since it would result in a swelling of this group through its strong interactions with the solvent. In the case of PC, the proximity of a local cation within the molecule results in the same folding of the head group seen in water because this interaction is less transient than the equivalent interactions with the solvent.

E. Refinement of neutron reflectometry
The ability to fit the NR data, as shown in Figure 3 indicates that the value found for the head component volume is consistent between the pair of measurements for the same system. It is clear, that again stable monolayers of the lipids are forming at the air-DES interface, and that the component volumes determined from XRR measurements are robust enough to be used in the modelling of NR data. Furthermore, the trends observed with increasing surface pressure in the XRR models, pertaining to the increasing tail thickness and decreasing solvent concentration in the head components are consistent with that found in the NR models.

F. Interparameter correlations
The use of Bayesian inference and MCMC sampling allowed for the probing of the probability distribution function for each parameter individually. However, it also enables the pairwise inter-parameter PDFs to be investigated, an example for DMPC at 30 mN m −1 is shown in Figure 5 (similar plots for the other XRR datasets and the NR datasets are available in Figures S3 -S21 of the ESI). These two-dimensional PDFs give important information about the correlations that are present between the parameters of the chemically-consistent model. The less circular in nature that the 2D PDF is, the greater the correlation that is present between the parameters, e.g. a north-east/south-westerly skewed PDF indicates a positive correlation (where an increase in one parameter correlates an increase in the other) while a north-west/south-easterly skewed PDF indicates a negative correlation (an increase in one correlates with a decrease in the other).
Substantial correlations are present in the parameters fitted to the XRR datasets, indicating important uncertainty that must be considered. In particular, as can be seen in Figure 5, there is a positive correlation between the lipid head thickness, d h , and the solvent concentration in the head layer, φ h . This correlation can be rationalised as a result of the SLD of the solvent and the head layer (which is ∼50 % solvated) being similar, and therefore the boundary between the two is not easy to define. Such correlations are unavoidable without considering many neutron contrasts of the lipid and solvent, due to the solvophilic nature of the lipid heads. Another important correlation is that between the head thickness and the tail thickness, d t , again this is due to the lack of a well-defined boundary between the head and tail layers. This is partially driven by the interfacial roughness that is present between the layers, with the correlation being more pronounced for phospholipids with shorter tails (e.g. there is a greater correlation for DLPC than DPPC).

IV. CONCLUSIONS
For the first time, stable phosphocholine and phosphatidylglycerol lipid monolayers have been observed and characterised on an ionic solvent surface. Until the emergence of ionic liquids and DES, only a limited number of molecular solvents exhibited the ability to promote selfassembly and, to the best of our knowledge, only water and formamide among those had demonstrated the formation of phospholipid monolayers at the air-liquid interface.
A physically and chemically constrained modelling approach and Bayesian analysis method was used to rationalise these measurements showing that the structures are remarkably similar at the air-DES interface to those previously observed at the air-water interface. This has the important implication that DES, therefore, offer the possibility of performing studies of model membranes in the absence of water. Such applications may include fundamental investigations of phospholipid monolayers in extreme environments (total or partial absence of water, cryogenic temperatures), protein-membrane interactions and development of new technologies for drug delivery. However, the PG component did show a significant difference; having a larger head component volume than observed for the same system in water. This suggests that the transfer of lipids to a DES is not just a simple substitution of the subphase. In this specific case, we have proposed an explanation based on the dissociation of the PG head component salt and the subsequent interaction with the DES. Finally, the use of MCMC sampling of the reflectometry model parameter space gives insight into the correlations present in XRR data modelling, that should be considered in work of this nature.
The ability to determine the head component volume was facilitated by access to easy to use, open-source software that allowed for the straightforward use of a custom, chemically-consistent model within the analysis of the XRR and NR measurements. Furthermore, to our knowledge, this work presents the first use of a chemicallyconsistent parameterisation to co-refine XRR measurements at different surface concentrations.