An Atomistic First-Principles Density Functional Theory Model for Single Layer Dry \textit{Stratum Corneum}

TL;DR

This study develops a first-principles density functional theory model for single-layer dry stratum corneum, revealing molecular interactions and accurately predicting vibrational spectra and mechanical properties relevant to skin biophysics.

Contribution

The paper introduces a novel atomistic DFT model for dry stratum corneum, linking molecular interactions to physical properties and matching experimental data.

Findings

Correct assignment of FTIR fingerprint bands

Identification of lipid-NMF vibrational coupling

Accurate prediction of shear modulus of dry SC

Abstract

Many questions concerning the biophysical and physiological properties of skin are still open. Skin aging, permeability, dermal absorption, hydration and drug transdermal delivery, are few examples of processes with its underlying mechanisms unveiled. In this work we present a first-principles density functional quantum atomistic model for single layer stratum corneum (SC) in order to contribute to unveil the molecular interactions behind the skin properties at this scale. The molecular structure of SC was modeled by an archetype of its hygroscopic proteic portion inside of the corneocytes, the natural moisturizing factor (NMF), coupled to glycerol molecules which represent the lipid fraction of SC. The vibrational spectra was calculated and compared to Fourier-Transform Infrared Absorption spectroscopy (FTIR) experimental data obtained on animal model of SC. We noticed that bands in the fingerprint region (800-1800 cm$^{-1}$) were correctly assigned. Moreover, our calculations revealed the existence of two coupled vibration between hydroxyl group of lipid and NMF methylene (1120 and 1160 cm$^{-1}$), which are of special interest since they probe the lipid-amino acid coupling. The model was also able to predict the shear modulus of dry SC in excellent agreement with the reported value on literature. Others physical/chemical properties could be calculated exploring the chemical accuracy and molecular resolution of our model. Research in dermatology, cosmetology, and biomedical engineering in the specific topics of drug delivery and/or mechanical properties of skin are examples of fields that would potentially take advantage of our approach.

Tables

Table 1. Table 1: Vibrational modes assignment for some FTIR bands observed in stratum corneum tissue (Fig. 3) based on SCmod:C(v) model . Bands indicated in red are of special interest since they probe the brick-mortar interaction.

Experimental	SCmod	Assignment of some bands
1040	1040	(CH2+CH3) asymmetric torsion of proline ring and cage
1120	1118	[OH stretching+CH2 wagging]GOL
1160	1160	[OH stretching+CH2]GOL+CH2 wagging
1240	1238	[COH-CHα-CH${}^{{\delta }_2}$] proline ring symmetric torsion
1400	1403	CH2 twisting ring proline
1500	1496	C-CH3 symmetrical bend (umbrella mode)
1560	1561	Amide II
1580	1581	Amide IV
1620	1619	Amide II+(CH2+CH3) asymmetric torsion (ring and cage)
1650	1653	Amide II+(CH2+CH3) asymmetric torsion (ring and cage)
1700	1694	Amide II+CH2 twisting

Equations

$$E_B=E_{STmod}-(E+E_{GOL})$$

$$\nu =\frac{1}{2\pi }\sqrt{\frac{S}{l_0}I}$$

Full text

Keywords: stratum corneum; skin; vibrational spectroscopy; FTIR; density functional theory; computational simulation; tissue computer model; drug delivery.

I Introduction

Skin is a complex, integrated, and dynamic organ that has many functions that go far beyond its role as a barrier to the environmentHongbo and Maibach (2004). It protects the body against external chemical and physical factors, takes part in the metabolic processes, plays a resorptive and thermoregulatory function, and it partakes in immunological processesHongbo and Maibach (2004); Boer et al. (2016). The basic structure of the skin comprises three component layers: epidermis, dermis, and subcutis or hypodermisBurkitt et al. (2014). The epidermis is the outer layer of the three layers that make up the skinBurkitt et al. (2014). The epidermis is divided into four layers, starting at the dermal junction with the basal cell layer and eventuating at the outer surface in the stratum corneum (SC)MD and MD (2013).

A large number of research tools is available for the investigation of biophysical parameters of skin, e.g., transepidermal water loss, electrical impedance, Raman spectroscopy, confocal spectroscopy, optical coherence tomography, magnetic resonance imaging, among othersByrne (2010). In spite of large experimental data available, several issues such as permeability and dermal absorption processes, remain unsolved. Moreover, skin aging processes, hydration and the role of water in life as well still lack clarityBall (2008).

The general morphology of SC is described within the “bricks and mortar” model (see Fig.1a and ref. Del Rosso and Levin (2011)). The corneocytes represent the bricks and the intercellular lamellar lipid membrane represents the mortar between the bricksDel Rosso and Levin (2011). Corneocytes are primarily proteic fiber comprised of keratin macrofibrils. They are protected externally by a cornified cell envelope, and are cohesively held together by corneodesmosomesDel Rosso and Levin (2011). Ceramides (40 – 50 %), cholesterol (25 %), and fatty acids (10 – 15 %) comprises the main components of the intercellular lamellar lipid membraneDel Rosso and Levin (2011). SC is composed by   20% of water, of which a fraction is tightly bound to hygroscopic molecules (the natural moisturizing factor, NMF) and lipids in the skin. NMF compounds are present in high concentrations within corneocytes and represent up to 20% to 30% of the dry weight of the SC. Much of the NMF is represented by amino acids (40%) and their derivatives derived from the proteolysis of epidermal filaggrinVerdier-Sévrain and Bonté (2007). Approximately 34% of the amino acids are neutralVerdier-Sévrain and Bonté (2007). The role of the NMF lies in the fact that its constituent chemicals are highly water soluble and hygroscopic. It plasticizes the SC keeping it resilient and preventing cracking and flaking due to mechanical stresses and is a hydrolitic medium for several important enzymatic reactions in the skinDel Rosso and Levin (2011); Verdier-Sévrain and Bonté (2007); Jokura et al. (1995).

Histological changes in the subcutaneous structures which give rise to biomechanical alterations are reported to occur along the aging process. The tissue roughness and elasticity are the main biomechanical parameters employed to determine the aging degreeMaiti et al. (2016). The above-mentioned changes are correlated to collagen site properties since collagen is directly related to structural functions of skin. In epidermis it was reported decreasing synthesis of keratinocytes, filaggrin, tranglutaminase, and collagen type VII due to agingMaiti et al. (2016). The collagen VII is also crucial for the function and stability of the extracellular matrix as it is an anchoring fibril collagen, where it binds type I and type III collagen, providing stability to the interstitial membrane and the basement membrane extracellular matrix structuresKarsdal (2016). The basal layer is responsible in large extension by the skin mechanical propertiesKarsdal (2016). Amao et al.Amano (2009) relates that the basal layer of sun-exposed skin becomes damaged, multilayered, and partly disrupted compared with that of sun-protected skin. Basal layer plays important roles in maintaining a healthy epidermis and dermis, and repeated damage destabilizes the skin, accelerating the aging processAmano (2009).

Skin is the the most accessible organ for drug delivery purposes. The transdermal delivery could for example, reduce first-pass metabolism associated with oral delivery, and is not painful as injectionsSchoellhammer et al. (2014). However, the major limitation for active compounds penetration through the skin is overcoming the SC. It serves as a rate-limiting lipophilic barrier against the uptake of chemical and biological toxins, as well as trans-epidermal water lossFox et al. (2011); Herman and Herman (2015).

A number of factors affects the dermal absorption, including, e.g., type of skin; physical-chemical properties of delivery systems; skin moisture level; external temperature; skin pre-treatment; among others. Thus, numerous in-vitro and in-vivo models are used to examine the penetration of active compounds through the skinHerman and Herman (2015); Godin and Touitou (2007).

In-silico approaches are in principle able to predict adsorption, distribution, metabolism, excretion, and toxicity on the basis of input data describing physical-chemical properties of the compounds to be delivered and on physiological properties of the exposed skinBaba et al. (2015); Selzer et al. (2015); Hatanaka et al. (2015); Verheyen et al. (2017). In-silico models also represent promising strategies to reduce or to avoid the usage of animal studies for drug developmentFrohlich (2016).

Classical molecular dynamics have been used to perform permeability studies at molecular levelMachado et al. (2016); Rocco et al. (2017). However, these approaches can obscure some relevant physical-chemical data from the environment since in the hard-sphere approximation the neast-nearest interactions are not implicitly considered. For example, the explicit structure of water and its interactions with the neighboring are not accounted forSato et al. (2015). At scale of real biology they are in fact only a small part of the overall pictureSaunders and Voth (2013). The main drawbacks of the these computational methods are the extensive computational burdenRen et al. (2010), absence of realistic phase function, and elastic light scattering modelsZhu and Liu (2013).

Atomistic models based on quantum mechanics calculations have the greatest predictive capability for materials properties. Detailed atomistic models for biological tissues are rare in the literature due to their inherent complexity (aperiodicity and large number of atoms). In order to correctly describe soft tissues and be computationally viable the atomistic model should be small enough ($<$ 150 atoms for actual quantum methods and computational processing capabilities) and include water-water and water-molecule effectsSato et al. (2015). From the quantum point of view, the Density Functional Theory (DFT) has been shown to be a promising method. It is considered a mature technique with good experimental reproducibilityCole and Hine (2016). In biological applications, their results are also very accurate. DFT methods have been employed to unravel complex subjects such as biomolecular structure and electronic properties, and the development of medicinal chemistry. Furthermore, studies of nanomaterials at the biomolecular interface have been considered new challenges and opportunities for application of DFT methodsCui et al. (2016).

PreviouslySato et al. (2015); Sato and Martinho (2018) we presented a DFT-based atomistic computational model for dermis tissue (STmod), which was able to successfully explain important experimental structural and general biochemical trends of normal and inflammatory oral tissues. Such model was pioneer in including quantum details for macroscopic tissue modeling. In the present work we report a SC model (SCmod) built from the STmod as starting point. The stability of some proposed variations in SCmod was also considered and compared to FTIR experimental data.

II Methodology

II.1 Computational details

It is computationally inviable represent the full molecular complexity of SC with available DFT methods. Our choice was to represent the SC by its hygroscopic proteic portion inside of the corneocytes, the NMF. The choice of a suitable unit cell, which forms an infinite reticulate representing SC (see Fig. 1b), was the main step for building the model. We considered the proline-rich C0 and D0 unit cells of STmod Sato et al. (2015) as starting point. Proline is a higroscopic neutral amino acid. It is present in NMF and collagen structuresBurkitt et al. (2014). The C structure has orthorhombic symmetry (P$2_{1}2_{1}2_{1}$ space group); unit cell parameters $a=13.5$ Å, $b=11.5$ Å, and $c=10.5$ Å, and internal cage volume (oblate spheroid) $V_{c}=1.857.3$ Å3. Another structural variant is the compact version of C structure, named D model, with unit cell parameters $a=17.4$ Å, $b=13.3$ Å, $c=9.8$ Åand $V_{c}=182.8$ Å3.

The inclusion of one glycerol (GOL) molecule in the unit cell mimetics the lipidic SC structure. We notice that the inclusion of a larger number of GOL molecules become computational cost of the simulation prohibitive. GOL molecule was added at 6 possible positions (indicated in i) to vi) in Fig. 1c). The periodic boundary conditions were established in two-dimensions (see Fig.1b). Thus, our model represents a single layer of SC. Each structure was first optimized using molecular mechanics. The MMFF94s force fieldHalgren (1999) uses the van der Waals term and includes cubic and quartic terms in the bond stretch, and cubic terms in angle bending potential energy. MMFF94s produces relatively flat dynamically averaged structuresHanwell et al. (2012). The models were pre-optimized in a first step using molecular mechanics with the MMFF94s force field implemented in the Avogadro softwareHalgren (1999).

DFTHohenberg and Kohn (1964) was used in order to obtain equilibrium geometries and harmonic frequencies. Calculations were implemented in the CPMD programcpm (2001) using the BLYP functiona Lee et al. (1988) augmented with dispersion corrections for the proper description of van der Waals interactionsLee et al. (1988); von Lilienfeld et al. (2005). For all simulations, the cutoff energy was considered up to 100 Ry. The linear response for the values of polarization and polar tensors of each atom in the system was calculated to evaluate the eigenvectors of each vibrational mode. In order to verify the feasibility of the model, vibrational calculations were performed and were compared to Fourier-Transform Infrared Absorption Spectroscopy (FTIR) experimental data obtained from SC animal model.

II.2 Animal model for skin

Full thickness pig ear skin was obtained from a local slaughter-house, and protocol was approved by Animal Ethical Committee (#007/2011). The SC was isolated from the dermis according to the methodology previously described by Kaushik and Michniak-KohnKaushik et al. (2012) with minor modifications. The skin was excised post-sacrifice from the outer part of pig ears separated from the underlying cartilage with a scalpel and skin samples were kept in a water bath at 37 0C with 5 mM Hepes buffer and 154 mM NaCl, pH 7.4, supplemented with trypsin (0.1%) for at least 4 hours. Subsequently, the skin segments were washed with the buffer solution (without trypsin) for 1 h. The skin samples were dispersed in deionized water, dried on filter paper and stored under silica atmosphere until use. The layers obtained in this method are $\sim$ 100 $\mu$m thick enclosing $\sim$ 80 $\mu$m of epidermis and $\sim$ 20 $\mu$m of SCJacobi et al. (2007).

II.3 FTIR

The spectra were collected in the Varian 610 FT-IR Micro-spectrometer in the reflectance method with spectral resolution of 2 cm*-1*, 200 scan per sample and 800 scans per background. The sample holder was a platinum recovered surface plate. Three spectra per sample were obtained, and each one was baseline corrected manually with Fityk softwarefit (2009) and normalized to 1. Then the average spectrum was calculated. For data comparison with computational calculations, each spectrum was simulated as a convolution of Gaussian lineshape peaks centered on the calculated frequencies using Fityk software. The linewidth was chosen to be 20 cm*-1*.

III Results and discussion

In order to infer about the influence of the GOL position on the unit cell stability we calculated the GOL binding energy ($E_{B}$) according to

[TABLE]

where $E_{STmod}$ is the energy of the SCmod structure (using C or D STmod unit cells), $E$ is the $C_{0}$ or $D_{0}$ energy, and $E_{GOL}$ is the energy of glycerol molecule. Results are shown in Figure 2. All energies appeared to be negative, which indicates that no external work needs be done to introduce one GOL inside the amino acid cage since and the GOL molecule spontaneously tends to bond in these sites. Thus, these structures are viable to accommodate glycerol molecules.

However, since SCmod:C with GOL in the iv and v positions and SCmod:D with GOL in the i, did not present negative frequencies in the vibrational spectra (not shown), we concluded that these models represent configurations closer to the global minimum energy without mechanical instabilities. Thus, we argued that the lowest energy SCmod:C with GOL in the v (labelled as SCmod:C(v) herein) will be the suitable model to mimic SC from the energetic and mechanical point of view. Figure 3 presents the FTIR experimental data for SC tissue in the fingerprint region (800 - 1800 cm*-1*). The peptidic bond Amide vibrations (1500 - 1700 cm*-1*) dominates the spectra as expected for a sample with high protein content. The solid line is the simulated spectra using the SCmod:C(v) model. The agreement between experimental and simulation is excellent below 1500 cm*-1*. For Amides region, the general spectral vibrations are present; however, important contributions arising from confined water have not been not taken into account in SCmod, which explains the observed difference.

Important pieces of information arise from the vibrational assignment indicated on Table I. This table contains some representative bands (Amides II and IV, some side chains and skeletal characteristics vibrations) since the full assignment would appear exhaustive smearing out the main features captured by the model. Those bands related to the GOL vibrations are indicated as red. The 1120 cm*-1* band is the hydroxyl stretching methylene wagging coupled vibration of GOL. The next one, at 1160 cm*-1*, has similar symmetry. It is a coupled vibration between hydroxyl group of GOL and NMF methylene. Both bands are of special interest since they would be useful to probe the lipid-amino acid coupling, a parameter of special relevance for studies in Dermatology, Cosmetics, and Biomedical Engineering. For example, it would be used to evaluate the interactions of substances applied in the skin with the proposes to contribute to organize the structure or improve the broken/damaged structures, chance the moisture level, among others. In particular, spectral changes were identified as arising from alterations in the concentration of the major constituents of NMF, important in maintaining SC hydrationXu et al. (2016); Choe et al. (2017).

It is well known that thickness and roughness of SC changes with gender, age and anatomical siteMaiti et al. (2016). Thus, determination of mechanical parameters of skin correlated to its structure or biochemical composition is of direct interestMaiti et al. (2016). The results of our calculations can also be used, e.g., to evaluate mechanical properties of SC. We will present an estimate of the shear modulus (S) of the SC. Considering the corneocyte-lipid interaction as going from NMF amino acids – mortar lipids modeled as an harmonic potential, the perturbation frequency would be approximated by as natural frequency of a spring-mass-like oscillatorFrench (2003)

[TABLE]

where $l_{0}$ and $I$ are the corneocyte depth and inertia moment, respectively. Using $\nu=1,160$ cm*-1* , and $l_{0}$ and $I$ from ref. Marks and Plewig (1983) we found $S=0.200$ MPa in good agreement with dry SC indentation measurements from ref. Leyva-Mendivil et al. (2015). The actual skin mechanical measurements protocols evolves indentations tests, rheometry, and others tribologycal techniques (see, e.g.Lewis and Masen (2015); Chen et al. (2015); Bostan et al. (2016)) which are unespecific to chemical composition and unresolved microscopically.

As other possible applications of the SCmod proposed model we can cite i) drug delivery studies where previsions about the outcomes of experiments designed to tracking specific molecules crossing SC (e.g., single molecule absorption spectroscopy, vibrational spectroscopy) could be simulated with high accuracy; ii) mechanical properties (like Young’s modulus, shear modulus, Poisson’s ratio, mechanical toughness) previsions of SC as function of moisture content variations or applied load; iii) study of corneocyte interactions by including more layers into the model. It is useful to understand physical/biochemical properties of tissues since they enable study the molecular interactions which would be present; to quantify the relative strength of each interaction; to perform dynamical simulations; and to predict properties such as electronic structure, transport, mechanical, with high degree of accuracySato et al. (2015).

IV Conclusions

We conclude that our proposed SCmod minimalist model for single dry SC is able to capture the main experimental trends of SC as probed by FTIR spectroscopy performed on animal model of SC. In the specific case considered the calculated vibrational spectra enabled assign the main vibrations active observed in SC and reveal specific vibrational probes to the “brick”-”mortar” interactions. It is important to stress that more pieces of information could be explored and analyzed using the large amount of data obtained. A large amount of modes were not described for shake of simplicity. Finally, the high-wavenumber region (up to 4000 cm*-1*) was not considered in the present work because important modes from confined water which are active in this spectral window are silent in our version of dry Scmod.

Research in dermatology, cosmetology, and biomedical engineering in specific fields of drug delivery and/or mechanical properties of skin can be performed taking advantage of chemical accuracy and molecular resolution of our proposed model quantum model. It is important to mention that SCmod could be improved by introducing specific characteristics under study and other physical, chemical and biological properties beyond of the selected examples included in the present paper could be calculated.

Acknowledgements. The authors would like to thank the Brazilian agencies Conselho Nacional de Desenvolvimento Científico e Tecnológico (CNPq - 311146/2015-5) and Fundação de Amparo à Pesquisa do Estado de São Paulo (FAPESP - 2011/19924-2) for the financial support. The authors would also thank the computational resources provided by Centro Nacional de Processamento de Alto Desempenho em São Paulo (CENAPAD-UNICAMP) and Sistema de Computação Petaflópica (Tier 0) (Santos Dumont-LNCC) under Sistema Nacional de Processamento de Alto Desempenho (SINAPAD) of the Ministério da Ciência, Tecnologia e Inovação (MCTI).

V References