Models for the Stratum Corneum Lipid Matrix: Effects of Ceramide

Nov 28, 2018 - †Department of Chemical and Biomolecular Engineering and ‡Biophysics Graduate Program, University of Maryland, College Park , Maryl...
0 downloads 0 Views 1MB Size
Subscriber access provided by Gothenburg University Library

B: Biomaterials and Membranes

Models for the Stratum Corneum Lipid Matrix: Effects of Ceramide Concentration, Ceramide Hydroxylation, and Free Fatty Acid Protonation Eric Wang, and Jeffery B. Klauda J. Phys. Chem. B, Just Accepted Manuscript • DOI: 10.1021/acs.jpcb.8b06188 • Publication Date (Web): 28 Nov 2018 Downloaded from http://pubs.acs.org on December 6, 2018

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Models for the Stratum Corneum Lipid Matrix: Effects of Ceramide Concentration, Ceramide Hydroxylation, and Free Fatty Acid Protonation Eric Wang1 and Jeffery B. Klauda1,2* 1Department

of Chemical and Biomolecular Engineering, 2Biophysics Graduate Program University of Maryland, College Park, MD 20742, USA *Corresponding Author: [email protected] Abstract

The primary barrier of the skin is provided by the outer layer known as the stratum corneum, making it the active target for transdermal drug delivery. A repeating structure of the SC known as the short periodicity phase is modeled as a bilayer composed of ceramides (Cer), cholesterol (Chol), and free fatty acids (FFA). This study simulates Cer/Chol/FFA bilayers composed of N-lignoceroylsphingosine, α-lignoceroylphytosphingosine, Chol, deprotonated lignoceric acid, and protonated lignoceric acid using the CHARMM36 force field to study the effects of Cer concentration, Cer hydroxylation, and FFA protonation. Calculated bilayer properties include surface area per lipid, area compressibility, deuterium order parameters, average Chol tilt angle, neutron scattering length density (NSLD) profiles, electron density profiles, bilayer thicknesses, interdigitation, hydrogen bonding, lipid clustering, and two-dimensional radial distribution functions. Based on comparison of NSLD profiles, the effect of Cer hydroxylation is eliminated as a factor causing bilayer thinning relative to experiment suggesting other reasons for this experimental-simualtion mismatch. Cer concentration, Cer hydroxylation, and FFA protonation all cause significant changes, but the degree of change depends on the specific property. Increasing Cer concentration induces transitions in Cer lipids from postured to hunched conformations, although both conformations possess high chain order. FFA protonation tends to strongly influence properties of the FFA, while not perturbing other lipids or the overall bilayer as significantly.

1 ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

1.

Page 2 of 35

Introduction A major goal of transdermal drug delivery is the circumvention of the stratum corneum

(SC), the outer layer of the skin, which presents a strong barrier to all but lipophilic molecules. Within the SC, there are two major components which form a “brick-and-mortar” structure: corneocytes (brick) and the surrounding lipid matrix (mortar),1 which is the major site of transdermal penetration.2 Ceramides (Cer) are known to be essential to the barrier function of the SC, and atopic dermatitis has been implicated in decreased Cer concentrations3 and changed Cer profiles.4-6 Other major constituents of the SC include cholesterol (Chol) and free fatty acid (FFA). Among Cer and FFA, there is extreme heterogeneity, and Cer in particular varies based on hydroxylation and chain length. The lipid matrix has been experimentally studied with a variety of methods including electron microscopy7-8 and infrared spectroscopy.9 Results from electron microscopy demonstrate the appearance of broad-narrow-broad electron lucent bands which suggest a multilayer structure. Additionally, significant work with X-ray diffraction has been performed2, 10-12 which suggest the existence of two repeating units in the lipid matrix: the short periodicity phase (SPP) and the long periodicity phase (LPP). Although the SPP is characterized as a simple bilayer, electron density calculations from X-ray and the appearance of electron lucent bands in electron microscopy suggest that multilayers correspond to the LPP. Molecular dynamics (MD) is a useful tool for the study of biochemical systems. While past simulations have been limited to the pico-nanosecond timescale due to computational cost,13 modern simulations regularly simulate systems up to the microsecond range and beyond. For studying the SC, approaches include coarse grain methods to observe large scale phenomena and superstructure.14-15 However, for studying the atomic interactions that occur in SC lamellae, allatom simulation remains the gold standard. Several studies using the united-atom GROMOS force field exist,16-19 and studies using the all-atom CHARMM36 (C36) force field, which generally agrees better with experiment for sphingolipids, have also been performed in recent years.20-21 Additionally, simulations using a custom parameters with the C36 force field have been performed, with the inclusion of a random walk molecular dynamics (RWMD) method that allows for overcoming kinetic trapping effects.22 RWMD, which is highly useful for biased initial builds, gives results consistent with conventional MD with a randomized initial configuration.

2 ACS Paragon Plus Environment

Page 3 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

The protonation state of FFA is variable in the SC, with the outer acid mantle possessing proportions of protonated and deprotonated FFA and the interior SC containing nearly completely deprotonated FFA.21 Many of the previously mentioned simulations use protonated FFA, but some also examine deprotonated FFA.21, 23 The effect of FFA protonation in phospholipid bilayers has been characterized as a reduction in rigidity and chain order,24 and one can hypothesize that FFA protonation would have similar effects in Cer/Chol/FFA mixtures, but the high order and rigidity of these bilayers could also inhibit strong effects. Similarly, while the hydroxylation of Cer from the sphingosine to the phytosphingosine backbone has been simulated with phosphatidylcholine (PC) bilayers and found significant ordering capable of inducing the gel phase at low concentrations,20 the effect in Cer/Chol/FFA bilayers is likely to be muted. The purpose of this study is to expand on previous work21 examining the structural properties of Cer/Chol/FFA bilayers as models of the SC. Previously, the C36 force field has been shown to have strong agreement with experimental NMR and neutron scattering results21 which supports the use of C36 for further study. The effects of Cer hydroxylation, Cer concentration, and FFA protonation are studied in light of previous results and connection to experiment. Furthermore, the systems of this work are capable of answering whether the hydroxylation of Cer can explain a previously found discrepancy between peak-to-peak distances between simulation and experimental neutron scattering. The bilayers of this work consist of N-lignoceroylsphingosine (CerNS), α-lignoceroylphytosphingosine (CerAP), Chol, deprotonated lignoceric acid (LA), and protonated lignoceric acid (LAP) at physiological skin temperatures of 32ºC. Systems of CerNS and CerAP are separately simulated over a range of Cer concentrations, and a system of LAP is simulated with CerNS in equimolar concentrations. Based on an array of calculated bilayer properties, Cer hydroxylation, Cer concentration, and FFA protonation affect the bilayer with the significance being property-dependent. 2.

Methods

2.1. System setup and simulation protocol Previously,21 we constructed similar ternary systems of CerNS/Chol/LA based on NMR experiments from Stahlberg et al.25 In this work, similar systems are studied with varying Cer concentration, Cer hydroxylation, and LA protonation. Systems of CerNS and CerAP are simulated at XCER of 0.20, 0.34, 0.50, and 0.70, temperatures are set to a physiological value of 3 ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 4 of 35

32ºC, and CerNS at XCER=0.34 is obtained from previous trajectories.21 An additional system containing CerNS at XCER=0.34 with LAP is simulated. For all systems, Chol and LA (or LAP) are maintained in a 1:1 ratio. All systems contained 100 lipids per leaflet with 24 waters per lipid as previously done, which is likely to yield fully hydrated systems.20 There are several differences with experiment. Namely, the hydration is higher than physiological system, which are composed of repeating lipid stacks at low hydration. However, this is useful in simulation to avoid issues with ions being present at low hydration, and lipid force fields are also parametrized with respect to hydrated systems. Therefore, these simplified models are useful for studying effects due to lipid concentration/structure, but they differ from experiment in hydration, lipid diversity, and structural complexity, and these are considered when comparing with experiment. In particular, the difference in hydration prevents association of lipids through the water layer and inhibits transverse hydrogen bonding. Current work on more complex LPP systems takes these factors into consideration. Potassium counter-ions were placed in LA-containing systems to remain electroneutral while the LAP-containing system did not contain ions. CHARMM-GUI Membrane Builder was used to generate triplicate replicas in rectangular boxes and tetragonal cells (X = Y ≠ Z), and these were equilibrated using the six-step process in NAMD.26-29 The agreement between conventional MD with randomized initial configurations and RWMD simulations22 justifies the use of random placement as done with CHARMM-GUI. Productions runs were simulated for 500-1000 ns such that all replicates had at least 200 ns of equilibrated data, which was determined by inspection of plots of surface area per lipid (SA/lip) as a function of time and the corresponding cumulative average. Table 1 contains the system parameters, including the production times which total to 16.5 μs.

4 ACS Paragon Plus Environment

Page 5 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Figure 1. a) Chemical structure of CerNS, CerAP, Chol, LA, and LAP. Atom labels include terminal tail carbons, RDF selection atoms (shown in red), and component surface area representative atoms. b) Snapshot of CerNS/Chol/LA at XCER=0.20 at the end of the simulation. CerNS molecules are shown in cyan with colored headgroups, Chol in blue, LA in red, water as purple dots, and potassium ions as brown beads. 5 ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 6 of 35

Table 1. System parameters for each bilayer system. Number of lipids and water are reported for the whole system. The reported time is the total production time for each replicate. System

# Cer # Chol

# LA

CerNS/Chol/LA 40 80 80 CerNS/Chol/LA* 68 66 66 CerNS/Chol/LA 100 50 50 CerNS/Chol/LA 140 30 30 CerAP/Chol/LA 40 80 80 CerAP/Chol/LA 68 66 66 CerAP/Chol/LA 100 50 50 CerAP/Chol/LA 140 30 30 CerNS/Chol/LAP 68 66 66 * Trajectories taken from Wang and Klauda.21

# Water

Time (ns)

4800 4800 4800 4800 4800 4800 4800 4800 4800

500 350 700 1000 500 500 800 1000 500

The CHARMM36 lipid force field30 with updates for sphingolipids31 was used with the TIP3P water model.32-33 CerAP is not available by default in the CHARMM36 force field and was constructed from a CerNS template using analogous values in the lipid for partial charges. NAMD was used with the constant molecule number, pressure, and temperature (NPT) ensemble. A constant pressure of 1 bar was maintained through the Nosé-Hoover Langevin-piston34-35 and constant temperature of 32 ºC was maintained through Langevin dynamics with a damping coefficient of 1 ps-1. van der Waals interactions were modeled using the Lennard-Jones potential with a force-switching function36 from 10 to 12 Å. Long-range electrostatics were modeled with the particle-mesh Ewald (PME) method using an interpolation order of 6 and a direct space tolerance of 10-6, and electrostatics were calculated every step. Bond lengths involving hydrogen were kept constant by the RATTLE algorithm.37 The simulations time step was 2 fs and coordinates were saved every 5 ps. Visual Molecular Dynamics38 (VMD) was used to create snapshots, and gnuplot39 was used to create plots. 2.2. Analysis Analysis was consistently performed using the last 200 ns of data, which was equilibrated based on SA/lip plots. The analyses presented are overall SA/lip, area compressibility modulus (KA), component SA/lip, apparent deuterium order parameters (SCD,app), Chol tilt, electron density profiles (EDPs), bilayer thickness, interdigitation, hydrogen bonding (H bonding), lipid clustering, and two-dimensional radial distribution functions (2D RDFs).

6 ACS Paragon Plus Environment

Page 7 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

CHARMM40-41 was used to perform most analysis and triplicates were used to calculate averages with standard errors. SCD,app as a function of carbon number was calculated by the formula, 𝑆𝐶𝐷,𝑎𝑝𝑝 =

|〈 𝑐𝑜𝑠 𝜃 ― 〉| 3

2

2

1

(1)

2

where θ is the angle between the C-H vector and the bilayer normal. Note that Eq. 1 is only applicable to calculating deuterium order parameters (SCD) if the lipids rotate around their long axis within the simulation timescale. Since this is not always the case in our simulations, they have been renamed to apparent order parameters and are not directly comparable to experimentally measured values. Chol tilt angle distribution was calculated by determining the angle between the C3 and C17 (Figure 1) carbons against the bilayer normal. The overall SA/lip was calculated as the area of the simulation box divided by the number of lipids per leaflet (always 100 in this study). Component SA/lip values were calculated by determining the X and Y coordinates of representative atoms for each lipid: C2S, C1F, and C4S for Cer lipids, O3 for Chol, and C1 for LA as shown in Figure 1. Quickhull42 constructs a Voronoi tessellation in which each atom is represented by a polygon which is then used to calculate polygonal areas. The sum of all polygonal areas is the total surface area, and summing areas of representative atoms gives individual SA/lip values. The area compressibility modulus (KA) was calculated by the formula,

𝐾𝐴 =

𝑘𝐵𝑇〈𝐴〉

(2)

𝑁𝜎2〈𝐴〉

where kB is Boltzmann’s constant, T is the absolute temperature, 〈𝐴〉 is the average SA of the entire leaflet, N is the number of lipids per leaflet, and 𝜎2〈𝐴〉 is the variance of the area. EDPs were calculated by recentering the bilayer such that the center is at Z=0 Å, calculating atomic densities with a slab thickness of 0.2 Å from -50 Å to 50 Å, and combining individual densities into group and total densities. The overall bilayer thickness (DB), headgroupto-headgroup distance (DHH), and hydrophobic distance (2DC) were calculated from EDPs. DB is the half-maximum distance of the water EDP, DHH is the peak-to-peak distance of the total EDP, and 2DC is the half-maximum distance of the hydrophobic EDP. SIMtoEXP43 was used to obtain a lipid NSLD profile to compare against experiment.44 From the profile, the peak-to-peak distance is calculated as the repeat distance. Interdigitation was calculated by integration of the density overlap profile,19 7 ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

+𝐿

𝜆𝑜𝑣 = ∫ ―𝐿4(

𝜌𝑡(𝑧)𝜌𝑏(𝑧)

𝑑𝑧

𝜌𝑡(𝑧) + 𝜌𝑏(𝑧))2

Page 8 of 35

(3)

where λov is the interdigitation parameter, ρt is the density of the top leaflet, ρb is the density of the bottom leaflet, z is the axis along the bilayer normal, and L is a length that extends past regions of non-zero density (50 Å in this work). Two-dimensional (2D) RDFs were calculate using a ∆r of 0.1 Å. The C2S atom of Cer, the O3 of Chol, and C1 of LA were considered for all RDFs (red atoms in Figure 1a). To calculate the number of intermolecular hydrogen (H) bonds (NHB), an H bond was defined as the donor-acceptor distance being less than 2.4 Å and the donor-acceptor angle being greater than 150°. The fraction of lipids with the considered H bond is reported. Intramolecular H bonds were not considered since they were previously found to be insignificant.20 To calculate lipid clustering, a python scikit-learn software45 was used with a density-based spatial clustering of applications with noise (DBSCAN) algorithm.46 A cut-off distance of 4.5 Å between head groups was used, which was chosen despite a distance of 5.0 Å better agreeing with RDFs because of errors in which clusters exceed the box size at 5.0 Å. Only cluster sizes of at least three lipids were considered. The C2S atom of Cer, O3 atom of Chol, and C1 atom of LA are considered. All reported p-values are based on the t-test for equal means. 3.

Results

3.1. Overall Surface Area per Lipid (SA/lip) and Area Compressibility Modulus (KA) The evolution of the overall SA/lip over time is a strong measure of the equilibration of the system. The time series plots, shown in Figure S1, demonstrate quick equilibration for lower Cer concentrations. XCER=0.20 and 0.34 generally show equilibration within 100 ns, XCER=0.50 equilibrates in 200-400 ns, and XCER=0.70 requires greater than 400 ns. Since Cer contains two chains while LA contains one chain, rotation is restricted and convergence to equilibrium diminishes as Cer concentration increases. As another view of equilibration, a dihedral angle in the upper sphingosine of Cer is plotted in Figure S2, which does not drift over the course of simulation. The dihedral is not stuck in the initial state since it samples small angles as well. A limitation of this study is that equilibration of phase separation, which occurs in experiment, is not possible in simulation. Such simulations would require systems containing up to a thousand lipids and micro-milliseconds of simulation, which is unavailable to current all-atom studies. 8 ACS Paragon Plus Environment

Page 9 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

The average SA/lip (Table 2) demonstrates a strong linear increase with XCER as the SA/lip increases nearly 17-20% for both CerNS and CerAP systems. This is to be expected given the naturally larger headgroup of Cer over Chol and LA. In contrast, the effects of Cer hydroxylation and LA protonation are much more subtle. CerAP has significantly more area than CerNS at XCER=0.70 (p=0.001) where the difference in hydroxylation has the greatest effect, but is generally insignificant at lower concentrations. The increase in SA/lip is somewhat surprising, as the increase in H bonding provided by the hydroxylation of CerAP generally causes a decrease in SA/lip. The SA/lip of CerNS with LAP is slightly lower than with LA (32.6 vs 32.8 Å2) but is nevertheless significant (p=0.01). This can be related to the SCD,app as discussed in a later section. Table 2. Overall SA/lip (always 100 lipids per leaflet including FFA) and area compressibility modulus (KA) for all systems. With respect to XCER, the SA/lip for CerNS and CerAP is linear. The slopes are 0.111 and 0.124, y-intercepts are 29.0 and 28.6, and R2 are 0.999 and 0.995 for CerNS and CerAP, respectively. System

XCER

SA/lip (Å2)

KA (N/m)

CerNS/Chol/LA CerNS/Chol/LA CerNS/Chol/LA CerNS/Chol/LA CerAP/Chol/LA CerAP/Chol/LA CerAP/Chol/LA CerAP/Chol/LA CerNS/Chol/LAP

0.20 0.34 0.50 0.70 0.20 0.34 0.50 0.70 0.34

31.18±0.02 32.81±0.01 34.50±0.08 36.79±0.03 31.18±0.08 32.62±0.07 34.53±0.04 37.41±0.06 32.63±0.04

3.9±0.1 3.4±0.1 3.78±0.06 3.26±0.03 3.9±0.2 4.5±0.2 4.2±0.4 3.3±0.1 3.7±0.1

The area compressibility modulus (KA) is a simple yet useful measure of a bilayer’s rigidity with fluid bilayers generally being around 0.3 N/m. Both Cer and Chol are known rigidifying agents and can increase KA by an order of magnitude over phospholipids. In Table 2, this is recapitulated as KA is consistently in the 3-5 N/m range. The relationship of KA with Cer concentration is not clear as KA both increases and decreases over XCER. At XCER of 0.34 and 0.50, CerAP causes a significant increase compared to CerNS, which is likely due to the greater Hbonding potential in CerAP. However, XCER=0.70 is not affected, indicating that hydroxylation is only significant within a range of Cer concentrations and does not necessarily become more significant at higher concentrations. The protonation of LA does not significantly affect KA (p=0.2) which contradicts previous simulation results that protonation of FFA leads to lower rigidity.24 9 ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 10 of 35

The critical difference lies in the membrane composition as the previous simulation analyzed FFA in a fluid phospholipid bilayer, while the SC systems of this work are solid-ordered. Therefore, factors that would easily perturb phospholipid bilayers such as FFA protonation become insufficient in SC models. This is well-exemplified by the inability of ethanol to penetrate into SC models at skin temperature using the C36 force field (Wang & Klauda, unpublished results). 3.2. Component Surface Areas The overall SA/lip averages between all lipids, but there is naturally a difference in SA/lip among different headgroups. The contribution of each lipid type is reported in Table 3. The SA/lip of Cer increases as XCER decreases, particularly in the XCER=0.20 case which yields values near 48 Å2. High XCER of 0.70 yield lower SA/lip around 43 Å2, which agrees with the pure Cer value of 43 Å2. The decrease from XCER=0.20 to 0.70 yields a 10% decrease. Chol and LA also undergo drops in SA/lip as XCER increases, although LA drops much more mildly. Chol drops from around 27 Å2 to around 19 Å2, a 30% decrease, while LA drops from 27.3 Å2 to 25 Å2, a 9% decrease. Based on visual snapshots, hexagonal packing seems to be present in both high and low Cer concentrations (Figure S3). Although close, the SA/lip of LA is not exactly half the SA/lip of Cer which can be attributed to lack of H-bonding pairs that limit packing. All SA/lip decrease with XCER, but the overall SA/lip still increases because Cer naturally possesses a larger SA/lip. The effect of Cer hydroxylation is generally negligible, but LA protonation induces a significant decrease in Cer (p=0.03). The effect on Chol is insignificant (p=0.06) and the effect on LA even more so (p=0.1). These SA/lip values are lower than Cer/Chol monolayers studied in experiment which found the SA/lip of Chol to be near 38 Å2.47 Differences contributing to this discrepancy are the lack of FFA and the use of monolayers. Experimentally, long-chain FFA increases packing density, and monolayers depend on the surface pressure. The authors note that their monolayer systems cannot be used to accurately draw conclusions about physiological membranes (lipid rafts in their case). Compared to a simulation study by Hofsass et al.48, the Chol SA/lip are in good agreement at low XCER, and their simulations contained phosphatidylcholine instead of Cer and FFA. Table 3. SA/lip (Å2) of specific lipid types in all systems. System

XCER

Cer (Å2)

Chol (Å2)

LA (Å2)

10 ACS Paragon Plus Environment

Page 11 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

CerNS/Chol/LA CerNS/Chol/LA CerNS/Chol/LA CerNS/Chol/LA CerAP/Chol/LA CerAP/Chol/LA CerAP/Chol/LA CerAP/Chol/LA CerNS/Chol/LAP

0.20 0.34 0.50 0.70 0.20 0.34 0.50 0.70 0.34

48.3±0.1 47.3±0.3 45.4±0.2 43.4±0.2 47.8±0.5 46.1±0.2 44.7±0.1 43.8±0.1 46.3±0.4

26.6±0.2 24.6±0.2 21.8±0.2 19.0±0.2 26.7±0.1 24.6±0.3 22.5±0.4 19.0±0.2 26.2±0.6

27.3±0.2 26.3±0.1 25.4±0.1 25.0±0.7 27.33±0.03 26.8±0.2 26.1±0.3 25.9±0.2 25.0±0.4

3.3. Apparent Deuterium Order Parameters (SCD,app) Under the assumption of lacking lamellar inversion (lipid flip-flop), deuterium order parameters (SCD) provide a quantitative measure of lipid chain order. Values near 0.4-0.5 generally reflect high order while values near and below 0.3 reflect disordered states. Gel and solid phases produce SCD that are difficult to interpret, but values are still calculable from simulation and are useful for illustrating order. Since the lipids do not rotate significantly about their long axis, the values obtained are not the same that experimental investigations would yield, so the nomenclature is adjusted to “apparent order parameters” (SCD,app). SCD,app as function of tail carbon number is provided in Figure 2, and angle distributions for select angles are given in Figure S4. In the fatty acid tails of Cer and LA, a characteristic plateau and drop-off is seen, which is due to free volume towards the bilayer center which allows for more rotational flexibility. In contrast, the sphingosine chain of Cer demonstrates a sharp drop in the double bond region before returning to the plateau region. CerAP does not demonstrate a drop due to the lack of a double bond and shows a SCD,app that is more similar to that of fatty acid chains. However, the plateau region of sphingosine still exists between the C8 and C14 carbons as CerNS does.

11 ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 12 of 35

Figure 2. SCD,app as a function of chain carbon number for a) fatty acid chain of Cer, b) sphingosine chain of Cer, and c) fatty acid chain of LA. Squares refer to CerNS systems while crosses refer to CerAP systems. Red refers to XCER=0.20, blue refers to XCER=0.34, green refers to XCER=0.50, orange refers to XCER=0.70, and black refers to LAP at XCER=0.34. Both Cer concentration and hydroxylation have significant effects on the SCD,app of Cer and LA, mostly in carbons towards the end of the lipid chain. For fatty acid chains, the plateau region 12 ACS Paragon Plus Environment

Page 13 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

between the C5 and C10 carbons is relatively consistent between all systems with values near 0.45 indicating the presence of solid-ordered phases. In general, increasing the Cer concentration causes increases in SCD,app for both CerNS and CerAP, and increasing the hydroxylation from CerNS to CerAP also causes greater order. Greater SCD,app is expected for CerAP, which possesses naturally greater H bonding from the phytosphingosine backbone, and greater SCD,app with increasing Cer concentration can be attributed to decreasing Chol concentration since Chol induces disorder because of its short length.49 This increase can also be related to the decrease in Cer component SA/lip at higher Cer concentrations. Additionally, the greater concentration of Cer increases H bonding and possibly allows for easier packing of the lipid chains into all trans conformations. The protonation of LA causes a modest increase in order for the fatty acid chain of Cer and no change in the sphingosine chain, but the order of LAP itself is significantly impacted. While LA generally does not increase until the C4 carbon due to the charged carboxylate, LAP increases in a similar fashion to the fatty acid tail of Cer. Furthermore, past the C8 carbon, LAP has lower order than LA which reflects previous results demonstrating lower order due to protonation in PC bilayers.24 This also indicates that the plateau, which normally exists between the C5 and C10 carbons, is significantly shorter in LAP. As seen in section 3.5, this is correlated with LAP being shifted toward the center relative to LA, which may be due to Born repulsion of like-charged atoms in LA. The overall SA/lip decreases due to protonation, and bilayer contraction is commonly attributed to greater chain order. The decrease in order for LA is incongruous, but this is resolved by considering the increase in Cer chain order. Since Cer chain order increases to a lesser extent than LA decreases, the overall SA/lip is dominantly affected by Cer, which is sensible due to its naturally larger headgroup and dual-chain structure. As previously mentioned, the component SA/lip for Cer increases in low XCER mixtures. While SCD,app decreases with XCER, another factor can be found by examining the fatty acid tail conformation in low XCER mixtures compared to high XCER mixtures (Figure 3). At XCER=0.20, the head-proximal end of the fatty acid tail juts out while the mid-chain carbons remain parallel to the bilayer normal, producing a noticeable kink. This conformation, termed here as the “posturing” conformation, is in stark contrast to the “hunched” conformation at XCER=0.70 in which the fatty acid tail is in close proximity with the sphingosine tail. The probability distributions as functions of chain distance reveal an increase in high-distance frequency (near 8 Å) for XCER=0.20, although 13 ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 14 of 35

both conformations exist and the posturing conformation is less common overall. As illustrated in Figure S5, the posturing conformation is aided by the positioning of Chol between the lipid chains, but the Chol is not necessarily engaged in an H-bond with the amide.

Figure 3. Snapshot of CerNS lipids at a) XCER=0.20 and b) XCER=0.70. c) Probability density of the two-dimensional lateral distance between lipid chains (defined using C16 carbons). 3.4. Tilt Angle of Chol Although the systems of this work do not possess significant tilt that is common in gel phases, Chol itself still possesses a slight propensity for tilt. The average tilt angle of the sterol ring is reported in Table 4. In general, Cer concentration does not have a significant effect on Chol tilt, which generally has values around 9º. Similarly, the protonation of LA does not change Chol tilt, which is not too surprising because the most significant effects would naturally occur for LAspecific properties. There is, however, an anomalous case in CerAP at XCER=0.70 which displays a higher tilt than in other systems (11.7º). Whereas other CerAP systems display slightly lower tilt than CerNS, the XCER=0.70 system displays the opposite effect. This is quite unexpected given the decrease in Chol SA/lip and increase in chain order, both of which would point towards decreased tilt. This discrepancy is most likely due transitioning into the gel phase (Figure S6). Whereas lower 14 ACS Paragon Plus Environment

Page 15 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Cer concentrations exhibit very little chain tilt, CerAP at 70% does possess a slight tilt characteristic of the gel phase. However, the tilt is not pronounced, so the system exists in an intermediate between the gel and solid ordered phase. Table 4. Average tilt angle for Chol in for all systems. System

XCER

Chol (º)

CerNS/Chol/LA CerNS/Chol/LA CerNS/Chol/LA CerNS/Chol/LA CerAP/Chol/LA CerAP/Chol/LA CerAP/Chol/LA CerAP/Chol/LA CerNS/Chol/LAP

0.20 0.34 0.50 0.70 0.20 0.34 0.50 0.70 0.34

8.9±0.1 8.9±0.1 9.3±0.1 9.2±0.2 8.7±0.2 8.4±0.1 8.5±0.2 11.7±0.5 8.8±0.1

3.5. Electron Density Profiles (EDPs), Bilayer Thickness, and Interdigitation Previously, CerNS at XCER=0.34 was compared to experimental neutron scattering length density (NSLD) profiles from Groen et al.44 and found relatively strong agreement in the shape of the NSLD profile, particularly in the bilayer center and presence of secondary peaks.21 However, a discrepancy in the peak-to-peak distance was found, which is 4.9 nm in C36 and 5.36 nm in experiment. In CerAP, this discrepancy remains true, and discussion follows in section 4.

Figure 4. Comparison of NSLD profiles obtained from C36 for CerNS/Chol/LA (CerNS) and CerAP/Chol/LA (CerAP) at XCER=0.34 and experimental neutron scattering of the SPP. A more detailed analysis is allowed through the inspection of EDP profiles (total and component) and the calculation of bilayer thicknesses. Based on overall EDPs (Figure 5), there is 15 ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 16 of 35

a strong dependence on Cer concentration. In particular, the secondary peaks at ±15 Å are most pronounced at low Cer concentration, which is due to the greater concentration of Chol, specifically due to the fused ring structure. As Cer concentration increases, the secondary peak amplitude decreases, but the primary peaks do not change as significantly. However, there is a shift outward which is indicative of general bilayer thickening, which is supported by thickness calculations (Table 5). Additionally, the effect of protonating LA causes lower primary peak density, likely due to the lack of significant charge on LA. Aside from the primary peak, much of the total EDP profile is similar to a deprotonated LA model. In the bilayer center, the appearance of a slight bump is indicative of significant interdigitation, which is provided by the relatively long fatty acid chain. As with secondary peak density, there is a profound dependence on XCER. This effect is subtle as snapshots (Figure S7) of the fluid center show only slight changes in interdigitation. However, the rise in methylene density (Figure S8) combined with the inward shift of LA density suggests that the deeper insertion of LA into the bilayer contributes to the greater interdigitation at low XCER.

Figure 5. Total electron density profiles for a) CerNS systems and b) CerAP systems. LAP refers to CerNS/Chol/LAP at XCER=0.34. 16 ACS Paragon Plus Environment

Page 17 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

The component EDPs for specific groups are shown in Figures S7-9. EDPs for Cer (Figure S8) show a strong dependence on Cer concentration. At XCER=0.20, the methylene and methyl densities display strong secondary peaks, while these attenuate and ultimately disappear at XCER=0.70. These secondary peaks are correlated with a lack of orthogonal movement due to greater rigidity provided by Chol. The EDP is not a dynamic property, so this lack of orthogonal movement is not to be thought of as a slower sampling of position but rather a change in the positional distribution. Additionally, Cer concentration induces a slight shifting outward of the sphingosine peaks which is correlated with a thickening of the bilayer. With this outward shift, the water density half-maximum adjusts from being inside of the sphingosine peak to outside, which is reflective of a decrease in water penetration. The hydroxylation of Cer induced by changing CerNS to CerAP largely causes no noticeable change in EDP shape. At XCER=0.70, there is a noticeable outward shift of the water EDP relative to the sphingosine EDP, reflecting lower water penetration. As expected, Chol EDPs exhibit greater densities at lower Cer concentrations (Figure S9), roughly proportionately to the change in Chol concentration. The most significant peak is the methine peak, with a maximum at ± 15 Å that explains the previously reported secondary peaks in total EDP. These, along with carbon peaks, show slight secondary peaks that are relatively consistent across Cer concentrations. CerAP systems also show similar curves to CerNS systems. In Figure S10, a comparison between CerNS at XCER=0.34 with LA and with LAP is shown. Overall, the curves are very similar, showing the same secondary peaks and peak positions irrespective of LA protonation. However, there is greater overlap of the Chol hydroxyl with carboxyl, which reflects a wider distribution of the carboxyl as well as a shift into the bilayer. From the EDPs, various qualitative trends about the bilayer thickness have been determined, i.e., DHH, DB, and 2DC, which vary in definition. Increasing Cer concentration causes increase in all three bilayer thicknesses, although they do not increase by the same amount. For example, DB increases ~12% from XCER=0.20 to 0.70 while DHH increases only 4% and 2DC increases 6%. 2DC is generally around 90% of DHH, which is higher than in standard phospholipid bilayers due to the lack of a phosphate group in Cer which decreases DHH. CerAP induces a thinner bilayer relative to CerNS based on DHH and 2DC. DB, by comparison, is either unaffected or actually increases due to CerAP. One could use DB-DHH as a simple measure of water penetration 17 ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 18 of 35

in which smaller values reflect greater penetration, with negative values being possible but rare. Universally, increasing Cer concentration and changing CerNS to CerAP increases DB-DHH, with CerNS at XCER=0.20 demonstrating a negative value. The effect of Cer concentration is attributed both to localization of Cer farther out of the bilayer than Chol and greater chain order which restricts water penetration. The combination of restricted water penetration into the headgroup region and a thickening of the hydrophobic region makes increasing Cer concentration a strong factor for barrier function. Increasing hydroxylation in CerAP restricts water penetration as well, which is simply attributed to greater H-bonding capacity. The protonation of LA causes significant thickening in DB (p=0.0001) and 2DC (p= 0.02) while DHH is not significantly affected (p=0.07). The reduced water penetrations is likely due to the lack of charge on LAP which reduces interaction with polar water molecules while the increase in 2DC is likely related to the small decrease in overall and component SA/lip. Table 5. Thicknesses (Å) of each system for different regions. DHH is the head-to-head distance, DB is the overall bilayer thickness, and 2DC is the hydrophobic thickness. Error bars are based on standard errors calculated from triplicates. System

XCER

DHH

DB

2DC

CerNS/Chol/LA CerNS/Chol/LA CerNS/Chol/LA CerNS/Chol/LA CerAP/Chol/LA CerAP/Chol/LA CerAP/Chol/LA CerAP/Chol/LA CerNS/Chol/LAP

0.20 0.34 0.50 0.70 0.20 0.34 0.50 0.70 0.34

48.9±0.2 49.1±0.1 49.9±0.2 51.1±0.5 47.9±0.2 48.6±0.1 49.6±0.1 50.2±0.1 49.7±0.2

48.07±0.02 50.27±0.03 52.7±0.1 54.72±0.02 48.4±0.1 51.2±0.2 53.3±0.2 54.55±0.01 51.4±0.1

43.62±0.03 44.4±0.1 45.33±0.03 46.4±0.1 43.16±0.04 43.98±0.04 44.87±0.01 45.6±0.1 44.7±0.1

In Table 6, the integration parameters, which can be thought of as the average length of each lipid that is overlapping/interdigitating, of all systems are reported. As with component SA/lip, there is a natural difference between lipids in interdigitation. Cer and LA, which both have long fatty acid chains 24 carbons long, demonstrate the most interdigitation (5-11 Å) while Chol demonstrates the least (1-2 Å). Unlike most other properties, Cer concentration does not have nearly as significant an impact on interdigitation. In CerNS, there is some slight decrease overall as Cer concentration increases, but the difference is not always statistically significant and shows some increases as well. While Cer concentration is not significant, Cer hydroxylation has a strong 18 ACS Paragon Plus Environment

Page 19 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

tendency to decrease interdigitation in Cer, while not strongly affecting other lipids. Perhaps the most interesting effect is the protonation of LA, which causes an increase in interdigitation of LA. As previously mentioned, the carboxyl of LA is shifted inward toward the hydroxyl of Chol, which increases the interdigitation of LA chains. The increase in LAP also affects Cer which interdigitates less in response. Chol does not significantly change (p=0.2), although the natural interdigitation of Chol is small to begin with.

19 ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 20 of 35

Table 6. Interdigitation (Å) of individual lipids for each system. System

XCER

Cer (Å)

CerNS/Chol/LA 0.20 7.66±0.07 CerNS/Chol/LA 0.34 7.46±0.04 CerNS/Chol/LA 0.50 7.07±0.01 CerNS/Chol/LA 0.70 7.4±0.1 CerAP/Chol/LA 0.20 5.66±0.1 CerAP/Chol/LA 0.34 5.52±0.03 CerAP/Chol/LA 0.50 5.4±0.2 CerAP/Chol/LA 0.70 6.0±0.4 CerNS/Chol/LAP 0.34 6.4±0.1 3.6. Hydrogen Bonding (H bonding)

Chol (Å)

LA (Å)

1.14±0.06 1.01±0.08 1.3±0.3 1.2±0.1 1.04±0.07 1.1±0.1 1.2±0.07 2.2±0.4 0.88±0.01

7.8±0.1 7.21±0.04 6.81±0.2 6.7±0.3 7.63±0.08 7.0±0.2 6.7±0.1 6.13±0.08 10.59±0.04

As sphingolipids, Cer have a high capacity for H bonding, which is markedly increased by the greater hydroxylation of the phytosphingosine backbone in CerAP. Guo et al.50 found that the hydroxyl orientation also changes from parallel to perpendicular to the bilayer normal, which also contributes to increased bonding. The overall H-bonding frequency is reported in Table 7 for specific lipid types. As expected, the most notable effect is due to hydroxylation. While CerNS frequency ranges between 0.2-0.4, CerAP ranges between 0.7-0.8. This also suggests that the effect of Cer concentration, which generally serves to increase H bonding, is not as significant in CerAP as it is in CerNS. In CerNS, however, the major shift occurs in the transition between XCER=0.20 and 0.34, while further increasing XCER causes little effect. Chol and LA are also affected by the hydroxylation of Cer as both increase in H bonding. In Chol, the effect is more subtle as H bonding only increases around 14%, while LA increases around 30% at XCER=0.70. Chol and LA also change with Cer concentration, both tending to increase in H bonding. It is interesting how all lipids increase in H bonding with increasing Cer concentration since the concentrations of Chol and LA decrease in response. Likely, this is related to the dual role of Cer as a bond donor and acceptor, which allows for H-bond network formation. The protonation of LA to LAP generally does not significantly affect H bonding except for Chol in which the H bonding drops from 0.195 to 0.13. Since the tilt angle is not affected and LAP moves closer to Chol according to the component EDPs, this is somewhat unexpected. As seen in section 3.7, the Chol-LA RDFs shift right due to LA protonation, indicating greater distance which decreases H bonding. Also, the protonation of LA reduces the ability of the protonable oxygen to serve as an H-bond acceptor, which is indicated by the sharp decrease in H bonding (Table S1). The EDP shift is driven not by 20 ACS Paragon Plus Environment

Page 21 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Chol, but by increased bonding with the amide of Cer through the carbonyl of LAP (Table S2). The ability of LAP to serve as a bond donor to Cer also contributes, although it is weaker than bonds in which LAP is an acceptor (Table S1). Table 7. Hydrogen bond (defined as when the distance between donor and acceptor is less than 2.4 Å and the angle is greater than 150º) count for overall lipids in each system. Values are reported as the fraction of lipids with the considered hydrogen bond. Lipid

System

XCER

NHB/lip

Cer

CerNS/Chol/LA CerNS/Chol/LA CerNS/Chol/LA CerNS/Chol/LA CerAP/Chol/LA CerAP/Chol/LA CerAP/Chol/LA CerAP/Chol/LA CerNS/Chol/LAP CerNS/Chol/LA CerNS/Chol/LA CerNS/Chol/LA CerNS/Chol/LA CerAP/Chol/LA CerAP/Chol/LA CerAP/Chol/LA CerAP/Chol/LA CerNS/Chol/LAP CerNS/Chol/LA CerNS/Chol/LA CerNS/Chol/LA CerNS/Chol/LA CerAP/Chol/LA CerAP/Chol/LA CerAP/Chol/LA CerAP/Chol/LA CerNS/Chol/LAP

0.20 0.34 0.50 0.70 0.20 0.34 0.50 0.70 0.34 0.20 0.34 0.50 0.70 0.20 0.34 0.50 0.70 0.34 0.20 0.34 0.50 0.70 0.20 0.34 0.50 0.70 0.34

0.204±0.001 0.37±0.03 0.35±0.03 0.388±0.007 0.7±0.03 0.75±0.04 0.819±0.008 0.81±0.02 0.37±0.02 0.164±0.001 0.195±0.008 0.23±0.02 0.29±0.03 0.183±0.008 0.22±0.02 0.263±0.005 0.33±0.03 0.13±0.01 0.18±0.01 0.22±0.01 0.28±0.02 0.43±0.05 0.254±0.009 0.34±0.04 0.55±0.01 0.57±0.05 0.19±0.01

Chol

LA

The frequencies of individual bond pairs, reported in Tables S1-2, provide greater mechanistic insight into the H bonding of the SC. In Table S1, bond pairs involving the Chol hydroxyl donor reveal negligible H bonding between Chol and either oxygen of the LA carboxyl. 21 ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 22 of 35

However, H bonding between Chol and Cer can be significant, particularly between the Chol hydroxyl and Cer carbonyl. For this same bond, Cer and Chol display inverse relationships with respect to Cer concentration. Adjusting the hydroxylation from CerNS to CerAP causes decreases in Chol-LA bonding, which is somewhat unexpected given the increase in overall H bonding. The answer lies in Chol-Cer bonding, which is quite significant at high Cer concentrations. In particular, the OH-O4 bond reaches a maximum of 0.157 at XCER=0.70. In comparison, the OHO2 bond is negligible. In Table S2, the bond pairs involving Cer as a donor are reported. For CerNS, the amide-carbonyl Cer-Cer bond is most significant, which is common in sphingolipids. However, this bond attenuates significantly in CerAP, up to a four-fold decrease at XCER=0.20. Since Cer H bonding increases in CerAP, this must be due to a shift of bonding towards hydroxyl bond pairs. The O4 hydroxyl is the most significant but not dominant over the other hydroxyls. The O3 hydroxyl, which is present in CerNS, undergoes an increase in H bonding in CerAP systems. Cer-LA bonding generally occurs through the amide donor over hydroxyl donors but is less significant than Cer-Cer bonding. The protonation of LA causes a two-fold increase in H bonding of the carbonyl oxygen but a ten-fold decrease in the bonding of the acidic oxygen due to the loss of charge caused by protonation. LAP-Cer bonding, in which LAP is the bond donor, is not reported in Table S2 but possesses a frequency of 0.036±0.005 for the carboxyl-carbonyl pair. This is similar to values for Cer-LA bonding and is ultimately not very significant to overall H bonding. 3.7. Lipid Clustering and 2D Radial Distribution Functions (2D RDFs) Cluster formation is a bilayer property that is commonly driven by headgroup interactions. In Table 8, the overall number of clustered lipid (NLT) and number of clusters (NCT) are reported. NLT and NCT generally increase with increasing Cer concentration, but the effect is not very pronounced and can decrease between XCER=0.50 and 0.70. The hydroxylation of Cer similarly increases NLT and NCT to a small degree, which agrees with an increase in H bonding. However, the protonation of LA to LAP causes no significant change in either. The effect of Cer concentration becomes clearer when the lipid ratios are considered, as is visually demonstrated in Figure 6. Naturally, lower Cer concentrations will see lower proportions of Cer participating in clusters. The cluster ratio (RC) is reported in Table S3 and demonstrates the expected increase in Cer ratio with XCER. RC for Chol and LA decrease in response. Chol and LA, despite being in

22 ACS Paragon Plus Environment

Page 23 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

equimolar ratios, display some differences in RC as Chol decreases less distinctly with XCER than LA. For example, RC for Chol in CerNS decreases 50% from XCER=0.20 to 0.70 while RC for LA decreases 70%. The protonation of LA does not significantly affect the RC of Cer, but does increase the RC of LA, possibly through the ability of the carboxyl to act as both an H-bond donor and acceptor. RC must also be considered relative to the bulk membrane concentration (RN) to better evaluate clustering tendencies. For Cer, RN is simply XCER. The values for RC-RN in Table S3 indicate that Cer is largely overrepresented in clusters while LA and Chol are underrepresented relative to the membrane concentration. Therefore, Cer possesses a greater clustering tendency than Chol and LA do. However, this does change at high Cer concentrations as clustering interactions are unable to match the increase in membrane concentration. Table 8. Average number of clustered lipids (NLT) and number of clusters (NCT) for overall clustering with a cutoff distance of 4.5 Å. Cluster sizes of 3 lipids and greater were considered. System

XCER

NLT

NCT

CerNS/Chol/LA CerNS/Chol/LA CerNS/Chol/LA CerNS/Chol/LA CerAP/Chol/LA CerAP/Chol/LA CerAP/Chol/LA CerAP/Chol/LA CerNS/Chol/LAP

0.20 0.34 0.50 0.70 0.20 0.34 0.50 0.70 0.34

45.7±0.2 50.1±0.7 52.7±0.2 50.6±0.9 48.7±0.4 52.9±0.7 57.7±0.8 58±2 49.3±0.2

11.6±0.1 12.3±0.2 12.7±0.1 12.20±0.03 11.82±0.04 12.4±0.1 13.11±0.04 13.1±0.4 12.0±0.1

23 ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 24 of 35

Figure 6. Visual representation of lipid clustering for top leaflet in a) CerNS/Chol/LA at XCER=0.20 and b) CerNS/Chol/LA at XCER=0.70. Yellow circles represent Cer lipids, blue circles represent Chol, and red circles represent LA. A cutoff distance of 4.5 Å is used. Clusters are taken from the final snapshot of the simulation rather than from a time average. In order to distinguish the interactions between differing lipids that contribute to clustering, RDFs can be calculated for various atom pairs. In Figure 7, RDFs for Cer-Cer and Chol-Chol clustering are reported. Cer-Cer RDFs in CerNS do not generally show a strong dependence on Cer concentration. CerAP differs since the RDF peaks are larger and show a clear trend of decreasing amplitude as Cer concentrations increases. The distance of the primary peak also tends to shift toward zero, indicating a slight contraction of Cer-Cer interactions at high Cer concentrations, possibly facilitated by increases in H bonding. The protonation of LA causes a significant decrease in the Cer-Cer RDF from 1.6 to 1.2, but does not cause a significant change in Chol-Chol RDFs. Due to the lowered concentration of Chol, Chol-Chol RDFs at XCER=0.70 display lower peaks and larger error intervals. Generally, there is a loss of the steady secondary peaks that are indicative of equidistant Chol shells, which have been reported in previous studies of high Chol systems.21,

51-52

LA-LA RDFs (Figure S11) see a significant increase in peak

amplitude from 1.2 to 1.6 due to the protonation of LA. The XCER=0.70 case displays large error 24 ACS Paragon Plus Environment

Page 25 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

intervals to the lowered LA concentration. Cer-Chol RDFs (Figure S12) similarly show a lack of Cer concentration dependence in CerNS but a strong dependence in CerAP. The relationship is inverse to Cer-Cer RDFs because the RDFs generally increase with Cer concentration rather than decrease. However, these interactions are weak and do not exceed 1.2. Similarly, Cer-LA RDFs (Figure S13) attain maximum values near 1.2 and do not show significant changes with Cer concentration and hydroxylation. LA protonation shows a strong smoothing of the RDF such that the short-range (< 5 Å) interactions are more frequent. Chol-LA RDFs (Figure S14) are much stronger with values reaching up to 1.6. Here, the effect of LA protonation appears to cause a slight shift away from zero while not affecting peak amplitude, indicating an expansion between Chol and LA.

Figure 7. Comparison of 2D RDFs for all systems. (a,b) represent Cer-Cer (C2S-C2S) RDFs and (c,d) represent Chol-Chol (O3-O3) RDFs. (a,c) represent systems involving CerNS and (b,d) represent systems involving CerAP. Lighter shades of each color represent ±standard error. 4.

Discussion and Conclusion Simulations of the SC vary in many aspects: lipid type, lipid compositions, FFA

protonation, number of layers, hydration, and force field to name a few. Although the majority of simulations use protonated FFA, Macdermaid et al.23 have simulated FFA dual layers with protonated and deprotonated FFA using C36 and found that the presence of deprotonated FFA allows for the sequestering of cations around points of contact between lamellae for systems large 25 ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 26 of 35

enough to see appreciable undulation of the bilayer. Cer and Chol were not included in these systems, likely due to the unavailability of the C36 force field for sphingolipids, and bilayer properties were not calculated, so comparison with this work is restricted. However, it does provide an explanation for the purpose of deprotonated FFA that complements this work. Our results indicate that while protonated FFA does induce significant changes in certain membrane properties, these are strongest in FFA-specific properties, and the effect on other lipids or the overall bilayer tend to be more muted. Thus, the role of FFA protonation is drawn into question: is it structurally important, or is it simply a consequence of a larger mechanism (i.e. the increased acidity necessary to kill microbes and activate acidic enzymes)? Currently, it appears that both aspects are important. FFA protonation regulates binding between lamellae and potentially affects permeability and the stability of low hydration multilayers. At the same time, sphingomyelinase enzymes in the SC which generate Cer are also activated by low pH, leading to protonated FFA.53 As mentioned in section 3.5, the peak-to-peak distance in the NSLD profiles is significantly different between simulation and experiment. Several factors were conjectured to cause this discrepancy including temperature, hydration differences, and hydroxylation. Among these, hydroxylation was suspected to have the greatest effect. Now, the effect of hydroxylation is probed by the replacement of CerNS with CerAP and is shown in Figure 4. However, the original hypothesis is not accurate as there is little difference in the NSLD profile between CerNS and CerAP. Since hydroxylation actually tends to produce thinner bilayers, further increasing the hydroxylation will not resolve the discrepancy either. At ±15 Å, the appearance of secondary peaks agrees well with experiment, suggesting that headgroup interactions could be the cause. The inclusion of acetate buffer in experiments could bind to the lipid headgroups, causing the widening of the headgroup peak. Additionally, at ±8 Å, the experimental profile possesses a minimum which corresponds to the terminal ends of the fatty acid and sphingosine chain meeting due to interdigitation. In the simulated profile, there is instead a shoulder which is potentially due to increased disorder of the bilayer center since the fatty acid tails do not protrude straight through the center. For comparison of bilayer properties, work from Gupta et al.16 using the GROMOS force field is useful. Generally, the results are in strong agreement with this work. For example, the component SA/lip of CerNS at XCER=0.70 is 0.43 nm2, which is near the pure Cer value of 0.43 26 ACS Paragon Plus Environment

Page 27 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

nm2 for C36 and 0.39 nm2 for GROMOS. The SA/lip for XCER=0.34 is 0.47 nm2 which is also similar to the GROMOS value of 0.45 nm2. The component EDP calculations for Cer and FFA show bumps in the bilayer center consistent with significant interdigitation found in this work. A discrepancy appears in the calculation of KA which is three-fold higher for GROMOS systems. As the KA is not significantly affected by protonation in C36, the possibility of this being a protonation effect is eliminated. The discrepancy is likely a force field effect because GROMOS is known to produce more compact bilayers.50 For lack of experimental KA data, one could use comparison between SCD,app as a general judge of accuracy. For this, C36 has been found to better replicate experimental results for Cer.20 Das et al.19 have also simulated Chol/Cer/FFA bilayers, although using the OPLS force field instead of C36 or GROMOS. For an equimolar mixture, the KA is roughly 5 N/m, which is intermediate between C36 and GROMOS. The interdigitation parameter λov attains a maximum of 30 Å, which is much larger than the maximum of 10 Å found in this work. This is likely due to differentiating the interdigitation of individual lipids in this work, which is useful given the different chain lengths. Das et al.54 have also simulated more complex systems, including a multilayer stack that contains the exceptionally long ceramide 1. In bilayers, they found that the number of H-bonding groups does not affect clustering, which disagrees with our result indicating that clustering increases, although somewhat marginally. Details regarding their cluster analysis are lacking, but presumably are performed through visual inspection, which could explain the inability to discern subtle clustering differences. Wang and Klauda previously simulated mixtures of Cer and PC,20 including both sphingosine and phytosphingosine Cer. Although Cer/PC bilayers are quite different from the Cer/Chol/FFA bilayers studied here, there are several similarities concerning the effect of phytosphingosine Cer. Primarily, the H bonding is found to increase due to greater hydroxylation, which is in agreement with experimental results from Rerek et al.55 Additionally, the O4 hydroxyl was found to be particularly prominent in H bonding, and the phytosphingosine Cer induces slightly greater clustering than sphingosine Cer in solid bilayers. However, a unique effect of phytosphingosine that is not seen here is the ability to induce a gel phase transition even at low concentrations. The effect on bilayer properties over a phase change is profound, and this cannot be seen in Cer/Chol/FFA bilayers since they naturally exist in solid-ordered states. Therefore, while the use of Cer/PC bilayers as SC models is surprisingly helpful for sphingosine Cer, the use of phytosphingosine Cer is precarious. Experimental work from Moore et al.56 agrees with an 27 ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 28 of 35

increase in H bonding with greater hydroxylation, although the phytosphingosine backbone was not specifically studied. NSLD work from both Groen et al.44 and Mojumdar et al.57 have indicated SPP repeat distances near 5.3 nm, which is larger than the 4.9 nm found in this work. Greater hydroxylation does not account for this difference since CerAP actually tends to slightly thin the bilayer compared to CerNS. Previous results from Bouwstra et al.12 for the SPP indicate that the peak-topeak distance (4.75 nm) can be lower than the repeat distance (5.2 nm) due to spacing between the peaks of adjacent SPP units, although peak-to-peak distances from X-ray are not necessarily comparable to values from NSLD. Bouwstra et al. also mention that the 4.75 nm distance was not accurately obtained due to the loss of phase information, suggesting that the repeat distance of 5.2 nm is better used as comparison with other studies. In both Bouwstra et al. and Groen et al., the lipids were prepared in acetate buffer, which potentially increases the repeat distance through headgroup binding and accounts for the difference with simulation. However, verifying this in simulation is not simple because the timescale of binding to SC headgroups is potentially long. At low XCER, Cer chains spread out resulting in a postured conformation, whereas the chains close together at high XCER, resulting in a hunched conformation. The lack of Cer crowding encourages spreading at low XCER, and the higher concentrations of Chol and LA have less impact on encouraging the hunched conformation. The posturing conformation is consistent with the inverted cone shape of Cer, which has been suggested to exist in the extended conformation at the interface between lamellae.58 In this conformation, the tails point in opposite directions with the polar head in the center. One can imagine that the posturing conformation, especially if it can exist in more extreme forms than found here, is a possible intermediate between the hairpin and extended conformations. Although this cannot be verified using a single bilayer simulation, future multilayer simulations possess the ability to determine whether the posturing conformation is recapitulated and whether it exists as an intermediate between the hairpin and extended conformations. In this work, ternary bilayers composed of Cer, Chol, and LA are studied for the effects of Cer concentration, Cer hydroxylation, and FFA protonation. Significant effects are noted for all bilayer properties, although the extent is property-dependent. For example, the component SA/lip of Cer displays dependence on Cer concentration with lower concentrations seeing increases in 28 ACS Paragon Plus Environment

Page 29 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

SA/lip, but the area compressibility and chain interdigitation is much less affected. Similarly, while the protonation of FFA causes significant changes in LA, such as the decrease of chain order and corresponding increase in interdigitation, the influence on other lipids and the overall bilayer is dampened. Cer hydroxylation, as expected, produces the most significant changes in H bonding and headgroup-dominated membrane properties. FFA protonation is eliminated as a possibility for causing differences with area compressibility in GROMOS, and Cer hydroxylation is eliminated as a factor for differences with NSLD thicknesses in experiment, although the presence of buffers in experiment remains a possible factor. For future studies, the protonation of FFA can be used to represent differing regions of the SC, whether in the outer acid mantle or the neutral interior, but the differences are not dramatic. For low hydration systems, the use of protonated FFA is simpler because it removes the necessity for counter ions. In future multilayer simulations, the existence of hairpin, extended, and posturing conformations remain an interesting subject to study, especially since the full story cannot be explored using single bilayers. 5.

Supporting Information Available Tables of H bonding including individual bond pairs and lipid clustering fractions. Figures

of SA/lip vs time, component EDPs, and 2D RDFs. Force field parameters for CerAP. 6.

Acknowledgements This research is supported by the NSF grant MCB-1149187, a grant to the University of

Maryland from the Howard Hughes Medical Institute through the Science Education Program, and a Goldwater Scholarship. The high performance computational resource used for this research is Deepthought2 which is maintained by the Division of Information Technology at the University of Maryland. Production runs used the Extreme Science and Engineering Discovery Environment (XSEDE) allocation by grant number MCB-100139 on Stampede2 and the computational resources at the Maryland Advanced Research Computing Center (MARCC). 7.

References

1.

Harding Clive, R., The Stratum Corneum: Structure and Function in Health and Disease.

Dermatologic Therapy 2004, 17 (s1), 6-15.

29 ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

2.

Page 30 of 35

Mojumdar, E. H.; Gooris, G. S.; Groen, D.; Barlow, D. J.; Lawrence, M. J.; Deme, B.; Bouwstra,

J. A., Stratum Corneum Lipid Matrix: Location of Acyl Ceramide and Cholesterol in the Unit Cell of the Long Periodicity Phase. Biochimica Et Biophysica Acta-Biomembranes 2016, 1858 (8), 1926-1934. 3.

Imokawa, G.; Abe, A.; Jin, K.; Higaki, Y.; Kawashima, M.; Hidano, A., Decreased Level of

Ceramides in Stratum Corneum of Atopic Dermatitis: An Etiologic Factor in Atopic Dry Skin? Journal of Investigative Dermatology 1991, 96 (4), 523-526. 4.

Ishikawa, J.; Narita, H.; Kondo, N.; Hotta, M.; Takagi, Y.; Masukawa, Y.; Kitahara, T.; Takema,

Y.; Koyano, S.; Yamazaki, S., et al., Changes in the Ceramide Profile of Atopic Dermatitis Patients. In J Invest Dermatol, United States, 2010; Vol. 130, pp 2511-4. 5.

Janssens, M.; van Smeden, J.; Gooris, G. S.; Bras, W.; Portale, G.; Caspers, P. J.; Vreeken, R. J.;

Hankemeier, T.; Kezic, S.; Wolterbeek, R., et al., Increase in Short-Chain Ceramides Correlates with an Altered Lipid Organization and Decreased Barrier Function in Atopic Eczema Patients. J Lipid Res 2012, 53 (12), 2755-66. 6.

van Smeden, J.; Janssens, M.; Kaye, E. C.; Caspers, P. J.; Lavrijsen, A. P.; Vreeken, R. J.;

Bouwstra, J. A., The Importance of Free Fatty Acid Chain Length for the Skin Barrier Function in Atopic Eczema Patients. Exp Dermatol 2014, 23 (1), 45-52. 7.

Breathnach, A. S.; Goodman, T.; Stolinski, C.; Gross, M., Freeze-Fracture Replication of Cells of

Stratum Corneum of Human Epidermis. Journal of Anatomy 1973, 114 (JAN), 65-81. 8.

Swartzendruber, D. C.; Wertz, P. W.; Kitko, D. J.; Madison, K. C.; Downing, D. T., Molecular-

Models of the Intercellular Lipid Lamellae in Mammalian Stratum Corneum. Journal of Investigative Dermatology 1989, 92 (2), 251-257. 9.

Mak, V. H. W.; Potts, R. O.; Guy, R. H., Does Hydration Affect Intercellular Lipid Organization

in the Stratum-Corneum. Pharmaceutical Research 1991, 8 (8), 1064-1065. 10.

Bouwstra, J. A.; Cheng, K.; Gooris, G. S.; Weerheim, A.; Ponec, M., The Role of Ceramides 1

and 2 in the Stratum Corneum Lipid Organisation. Biochimica et Biophysica Acta (BBA) - Lipids and Lipid Metabolism 1996, 1300 (3), 177-186. 11.

White, S. H.; Mirejovsky, D.; King, G. I., Structure of Lamellar Lipid Domains and Corneocyte

Envelopes of Murine Stratum-Corneum - An X-ray Diffraction Study. Biochemistry 1988, 27 (10), 37253732. 12.

Bouwstra, J. A.; Gooris, G. S.; Dubbelaar, F. E. R.; Weerheim, A. M.; Ijzerman, A. P.; Ponec, M.,

Role of Ceramide 1 in the Molecular Organization of the Stratum Corneum Lipids. Journal of Lipid Research 1998, 39 (1), 186-196.

30 ACS Paragon Plus Environment

Page 31 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

13.

Venable, R. M.; Zhang, Y. H.; Hardy, B. J.; Pastor, R. W., Molecular-Dynamics Simulations of a

Lipid Bilayer and of Hexadecane - an Investigation of Membrane Fluidity. Science 1993, 262 (5131), 223-226. 14.

Wennberg, C. L.; Narangifard, A.; Lundborg, M.; Norlén, L.; Lindahl, E., Structural Transitions

in Ceramide Cubic Phases during Formation of the Human Skin Barrier. Biophys. J. 2018, 114 (5), 11161127. 15.

Paloncyova, M.; Vavrova, K.; Sovova, Z.; DeVane, R.; Otyepka, M.; Berka, K., Structural

Changes in Ceramide Bilayers Rationalize Increased Permeation through Stratum Corneum Models with Shorter Acyl Tails. Journal of Physical Chemistry B 2015, 119 (30), 9811-9819. 16.

Gupta, R.; Rai, B., Molecular Dynamics Simulation Study of Skin Lipids: Effects of the Molar

Ratio of Individual Components over a Wide Temperature Range. Journal of Physical Chemistry B 2015, 119 (35), 11643-11655. 17.

Gupta, R.; Dwadasi, B. S.; Rai, B., Molecular Dynamics Simulation of Skin Lipids: Effect of

Ceramide Chain Lengths on Bilayer Properties. Journal of Physical Chemistry B 2016, 120 (49), 1253612546. 18.

Gupta, R.; Sridhar, D. B.; Rai, B., Molecular Dynamics Simulation Study of Permeation of

Molecules through Skin Lipid Bilayer. Journal of Physical Chemistry B 2016, 120 (34), 8987-8996. 19.

Das, C.; Noro, M. G.; Olmsted, P. D., Simulation Studies of Stratum Corneum Lipid Mixtures.

Biophysical Journal 2009, 97 (7), 1941-1951. 20.

Wang, E.; Klauda, J. B., Molecular Dynamics Simulations of Ceramide and Ceramide-

Phosphatidylcholine Bilayers. J. Phys. Chem. B 2017, 121 (43), 10091-10104. 21.

Wang, E.; Klauda, J. B., Simulations of Pure Ceramide and Ternary Lipid Mixtures as Simple

Interior Stratum Corneum Models. J. Phys. Chem. B 2018, 122 (10), 2757-2768. 22.

Moore, T. C.; Hartkamp, R.; Iacovella, C. R.; Bunge, A. L.; McCabe, C., Effect of Ceramide Tail

Length on the Structure of Model Stratum Corneum Lipid Bilayers. Biophys J 2018, 114 (1), 113-125. 23.

MacDermaid, C. M.; DeVane, R. H.; Klein, M. L.; Fiorin, G., Dehydration of multilamellar fatty

acid membranes: Towards a computational model of the stratum corneum. The Journal of Chemical Physics 2014, 141 (22), 22D526. 24.

Ramos, A. P.; Lague, P.; Lamoureux, G.; Lafleur, M., Effect of Saturated Very Long-Chain Fatty

Acids on the Organization of Lipid Membranes: A Study Combining H-2 NMR Spectroscopy and Molecular Dynamics Simulations. Journal of Physical Chemistry B 2016, 120 (28), 6951-6960. 25.

Stahlberg, S.; Skolova, B.; Madhu, P. K.; Vogel, A.; Vavrova, K.; Huster, D., Probing the Role of

the Ceramide Acyl Chain Length and Sphingosine Unsaturation in Model Skin Barrier Lipid Mixtures by H-2 Solid-State NMR Spectroscopy. Langmuir 2015, 31 (17), 4906-4915. 31 ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

26.

Page 32 of 35

Jo, S.; Kim, T.; Iyer, V. G.; Im, W., CHARMM-GUI: a Web-based Graphical User Interface for

CHARMM. J Comput Chem 2008, 29 (11), 1859-65. 27.

Jo, S.; Kim, T.; Im, W., Automated Builder and Database of Protein/Membrane Complexes for

Molecular Dynamics Simulations. Plos One 2007, 2 (9), 9. 28.

Wu, E. L.; Cheng, X.; Jo, S.; Rui, H.; Song, K. C.; Davila-Contreras, E. M.; Qi, Y. F.; Lee, J. M.;

Monje-Galvan, V.; Venable, R. M., et al., CHARMM-GUI Membrane Builder Toward Realistic Biological Membrane Simulations. Journal of Computational Chemistry 2014, 35 (27), 1997-2004. 29.

Phillips, J. C.; Braun, R.; Wang, W.; Gumbart, J.; Tajkhorshid, E.; Villa, E.; Chipot, C.; Skeel, R.

D.; Kale, L.; Schulten, K., Scalable molecular dynamics with NAMD. J. Comput. Chem. 2005, 26 (16), 1781-1802. 30.

Klauda, J. B.; Venable, R. M.; Freites, J. A.; O'Connor, J. W.; Mondragon-Ramirez, C.;

Vorobyov, I.; Tobias, D. J.; MacKerell, A. D.; Pastor, R. W., Update of the CHARMM all-atom additive force field for lipids: Validation on six lipid types. J. Phys. Chem. B 2010, 114 (23), 7830-7843. 31.

Venable, R. M.; Sodt, A. J.; Rogaski, B.; Rui, H.; Hatcher, E.; MacKerell, A. D.; Pastor, R. W.;

Klauda, J. B., CHARMM All-Atom Additive Force Field for Sphingomyelin: Elucidation of Hydrogen Bonding and of Positive Curvature. Biophysical Journal 2014, 107 (1), 134-145. 32.

Jorgensen, W. L.; Chandrasekhar, J.; Madura, J. D.; Impey, R. W.; Klein, M. L., Comparison of

Simple Potential Functions for Simulating Liquid Water. J. Chem. Phys. 1983, 79 (2), 926-935. 33.

Durell, S. R.; Brooks, B. R.; Bennaim, A., Solvent-Induced Forces between Two Hydrophilic

Groups. J. Phys. Chem. 1994, 98 (8), 2198-2202. 34.

Feller, S. E.; Zhang, Y. H.; Pastor, R. W.; Brooks, B. R., Constant-Pressure Molecular-Dynamics

Simulation - the Langevin Piston Method. Journal of Chemical Physics 1995, 103 (11), 4613-4621. 35.

Martyna, G. J.; Tobias, D. J.; Klein, M. L., Constant pressure molecular dynamics algorithms. J.

Chem. Phys. 1994, 101 (5), 4177-4189. 36.

Steinbach, P. J.; Brooks, B. R., New spherical-cutoff methods for long-range forces in

macromolecular simulation. J. Comput. Chem. 1994, 15 (7), 667-683. 37.

Andersen, H. C., Rattle - a Velocity Version of the Shake Algorithm for Molecular-Dynamics

Calculations. Journal of Computational Physics 1983, 52 (1), 24-34. 38.

Humphrey, W.; Dalke, A.; Schulten, K., VMD: Visual Molecular Dynamics. Journal of

Molecular Graphics & Modelling 1996, 14 (1), 33-38. 39.

Racine, J., gnuplot 4.0: A Portable Interactive Plotting Utility. Journal of Applied Econometrics

2006, 21 (1), 133-141.

32 ACS Paragon Plus Environment

Page 33 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

40.

Brooks, B. R.; Brooks, C. L.; Mackerell, A. D.; Nilsson, L.; Petrella, R. J.; Roux, B.; Won, Y.;

Archontis, G.; Bartels, C.; Boresch, S., et al., CHARMM: The biomolecular simulation program. J. Comput. Chem. 2009, 30 (10), 1545-1614. 41.

Brooks, B. R.; Bruccoleri, R. E.; Olafson, B. D.; States, D. J.; Swaminathan, S.; Karplus, M.,

CHARMM - a Program for Macromolecular Energy, Minimization, and Dynamics Calculations. J. Comput. Chem. 1983, 4 (2), 187-217. 42.

Barber, C. B.; Dobkin, D. P.; Huhdanpaa, H., The Quickhull algorithm for convex hulls. Acm

Transactions on Mathematical Software 1996, 22 (4), 469-483. 43.

Kučerka, N.; Katsaras, J.; Nagle, J., Comparing Membrane Simulations to Scattering

Experiments: Introducing the SIMtoEXP Software. J. Membr. Biol. 2010, 235 (1), 43-50. 44.

Groen, D.; Gooris, G. S.; Barlow, D. J.; Lawrence, M. J.; van Mechelen, J. B.; Demé, B.;

Bouwstra, J. A., Disposition of Ceramide in Model Lipid Membranes Determined by Neutron Diffraction. Biophysical Journal 2011, 100 (6), 1481-1489. 45.

Pedregosa, F.; Varoquaux, G.; Gramfort, A.; Michel, V.; Thirion, B.; Grisel, O.; Blondel, M.;

Prettenhofer, P.; Weiss, R.; Dubourg, V., et al., Scikit-learn: Machine Learning in Python. Journal of Machine Learning Research 2011, 12, 2825-2830. 46.

Martin Ester, H.-P. K., Jörg Sander, Xiaowei Xu, A Density-based Algorithm for Discovering

Clusters in Large Spatial Databases with Noise. Proceedings of 2nd International Conference on Knowledge Discovery and Data Mining 1996, 96 (34), 226-231. 47.

Scheffer, L.; Solomonov, I.; Weygand, M. J.; Kjaer, K.; Leiserowitz, L.; Addadi, L., Structure of

Cholesterol/Ceramide Monolayer Mixtures: Implications to the Molecular Organization of Lipid Rafts. Biophysical journal 2005, 88 (5), 3381-3391. 48.

Hofsäß, C.; Lindahl, E.; Edholm, O., Molecular Dynamics Simulations of Phospholipid Bilayers

with Cholesterol. Biophysical Journal 2003, 84 (4), 2192-2206. 49.

Schroeter, A.; Stahlberg, S.; Školová, B.; Sonnenberger, S.; Eichner, A.; Huster, D.; Vávrová, K.;

Hauß, T.; Dobner, B.; Neubert, R. H. H., et al., Phase Separation in Ceramide[NP] Containing Lipid Model Membranes: Neutron Diffraction and Solid-State NMR. Soft Matter 2017, 13 (10), 2107-2119. 50.

Guo, S.; Moore, T. C.; Iacovella, C. R.; Strickland, L. A.; McCabe, C., Simulation Study of the

Structure and Phase Behavior of Ceramide Bilayers and the Role of Lipid Headgroup Chemistry. Journal of Chemical Theory and Computation 2013, 9 (11), 5116-5126. 51.

Wang, E.; Klauda, J. B., Examination of Mixtures Containing Sphingomyelin and Cholesterol by

Molecular Dynamics Simulations. J. Phys. Chem. B 2017, 121 (18), 4833-4844.

33 ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

52.

Page 34 of 35

Adams, M.; Wang, E.; Zhuang, X.; Klauda, J. B., Simulations of simple Bovine and Homo

sapiens outer cortex ocular lens membrane models with a majority concentration of cholesterol. Biochim. Biophys. Acta, Biomembr. 2018, 1860 (10), 2134-2144. 53.

Elias, P. M., Stratum Corneum Acidification: How and Why? Experimental dermatology 2015,

24 (3), 179-180. 54.

Das, C.; Olmsted, P. D., The Physics of Stratum Corneum Lipid Membranes. Philosophical

Transactions of the Royal Society a-Mathematical Physical and Engineering Sciences 2016, 374 (2072), 16. 55.

Rerek, M. E.; Chen, H. C.; Markovic, B.; Van Wyck, D.; Garidel, P.; Mendelsohn, R.; Moore, D.

J., Phytosphingosine and Sphingosine Ceramide Headgroup Hydrogen Bonding: Structural Insights through Thermotropic Hydrogen/Deuterium Exchange. Journal of Physical Chemistry B 2001, 105 (38), 9355-9362. 56.

Moore David, J.; Rerek Mark, E.; Mendelsohn, R., Role of Ceramides 2 and 5 in the Structure of

the Stratum Corneum Lipid Barrier. International Journal of Cosmetic Science 2001, 21 (5), 353-368. 57.

Mojumdar, E. H.; Groen, D.; Gooris, G. S.; Barlow, D. J.; Lawrence, M. J.; Deme, B.; Bouwstra,

J. A., Localization of Cholesterol and Fatty Acid in a Model Lipid Membrane: A Neutron Diffraction Approach. Biophysical Journal 2013, 105 (4), 911-918. 58.

Kiselev, M. A., Conformation of Ceramide 6 Molecules and Chain-Flip Transitions in the Lipid

Matrix of the Outermost Layer of Mammalian Skin, the Stratum Corneum. Crystallography Reports 2007, 52 (3), 525-528.

34 ACS Paragon Plus Environment

Page 35 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

TOC Graphic

35 ACS Paragon Plus Environment