In silico model to predict dermal absorption of chemicals in finite dose conditions

Abstract

The development of in silico approaches that can estimate the dermal absorption of chemicals exposed in practical conditions is highly anticipated. In the present study, an in silico model to estimate both the dermal absorption rate and dermal permeation profile was developed for the application of chemicals in finite dose conditions. Forty-three chemicals with molecular weights in the range 116–362 and logK_o/w in the range 1.1–4.5 were used to develop an in silico model. A gradient boosting tree approach was applied to estimate permeation parameters for diffusion and partition coefficients of the chemicals in skin using physicochemical parameters of the chemicals such as molecular weight, lipophilicity, and the highest and lowest occupied molecular orbitals as the descriptor. In addition, 11 chemicals with different molecular weights and lipophilicities were applied on excised human skin in a finite dose condition, and dermal absorption profiles were obtained. Consideration of donor-solvent evaporation time, saturated concentrations of the chemicals, and donor-solvent coverage area on the skin surface, in addition to estimated skin permeation parameters of the chemicals, showed comparatively good dermal absorption profiles, although some cases of underestimation of dermal absorption were identified. It will be necessary to verify the accuracy of this model through experiments using more chemicals. However, the obtained results suggested that the established model may be valid to estimate the dermal absorption of chemicals in practical conditions.

INTRODUCTION

Recently, next-generation risk assessment has been developed as an exposure-led, hypothesis-driven approach that integrates non-animal experiments, such as in silico and chemical methods, to evaluate the safety of cosmetic ingredients (Beal et al., 2022; Paul Friedman et al., 2020; Zobl et al., 2024). This assessment considers the absorption, distribution, metabolism, and excretion (ADME) of exposed chemicals in humans. Chemical transport across the skin is a multi-step process (Guy and Hadgraft, 1985; Sugibayashi et al., 2010). Following the application of a formulation containing chemicals, the first step involves the partitioning of the chemicals from the formulation into the stratum corneum (SC) layer (K_SC), which is then followed by diffusion of the chemical in the SC layer (D_SC). A chemical that has passed through the SC is then distributed into the viable epidermis and dermis (VED) (K_VED), and the chemical in the VED diffuses in the layer (D_VED). The permeability coefficient of chemicals through the skin (K_p) was obtained from the steady state flux of dermally exposed chemicals applied in infinite dose conditions. The reciprocal of K_p (R_p) indicates the resistance to drug permeation through the skin (Ghanem et al., 1992). Therefore K_p can be further expressed as the sum of the resistance to permeation through the SC (1/K_{p_SC}) and the resistance through the VED (1/K_{p_VED}); 1/K_p=1/K_{p_SC}+1/K_{p_VED}. In addition, K_{p_SC} and K_p__VED can be determined as K_{p_SC} = D_SC∙K_SC /L_SC and K_{p_VED} = D_VED∙K_SC/L_VED, respectively, where L_SC and L_VED are the thickness of the SC and VED, respectively. Many approaches (Potts and Guy, 1992; Barratt, 1995; Moss and Cronin, 2002; Lien and Gao, 1995) have been reported to estimate K_p values using physicochemical properties of chemicals such as molecular weight (MW) and n-octanol-water coefficient (K_o/w). The maximum permeation flux was calculated as the product of K_p and the concentration of chemicals, and the steady-state blood concentration was determined using the maximum permeation flux, application skin area, and the total clearance from the body of chemicals (Bouzom et al., 2012; Hartmanshenn et al., 2016).

In addition to the prediction of steady-state blood concentration using K_p values, steady-state skin concentration was also determined using skin permeation parameters of D_SC, D_VED, K_SC, and K_VED derived from Fick’s second law of diffusion (Sugibayashi et al., 2010). Arce et al. (2020) analyzed the in vitro skin permeation profiles in infinite dose conditions to obtain skin permeation parameters, and they reported that the skin concentration was well predicted with skin permeation parameters when rhododendrol solution was applied in a finite dose condition.

The Cosmetics Europe ADME Task Force compared 6 in silico models to estimate systemically absorbed chemicals after dermal application in finite dose conditions (Grégoire et al., 2021). In this comparison, the skin permeation rate against the applied amount of chemical absorbed over 24 hr was evaluated. These estimations were highly variable for the safety assessment of dermally exposed chemicals. However, prediction of the permeation amount-time profile of dermally exposed chemicals can be more quantitative and may be useful in predicting the blood concentration-time profile. Prediction of skin permeation profiles requires prediction of skin permeation parameters, such as D_SC, D_VED, and K_SC, not K_p values. Several in silico models to predict skin permeation parameters, such as D_SC, D_VED, and K_SC, have been published. Wang et al. (2007) published prediction models of D_SC and K_SC of the solute. Lian et al. (2008) and Mitragotri et al. (2011) reported the D of a solute in the lipid matrix in the SC, D_lip, and the partition of the solute to lipid bilayers in the SC, K_lip. In addition, Kretsos et al. (2008) yielded prediction models of K_VED and D_VED in mammalian dermises, such as human, rat, guinea pig, and mouse. However, a few models have been developed using the same chemical data set, and these reported models have been developed with different experimental conditions.

Recently, Ellison et al. (2020) determined D_SC, D_VED, K_SC, and K_VED values of 50 compounds in human skin with the same experimental conditions. Although 50 chemicals are not enough to cover a large number of chemicals, these skin permeation parameters are helpful in developing an in silico penetration model based on Fick`s second law of diffusion. In their report, they also described that empirical in silico models are generally based on data sets obtained from several different publication sources, thus resulting in large heterogeneity. However, they investigated the relationship between single physicochemical parameters of the chemicals and the K or D measured in their work, in spite of a weak relationship between them.

A gradient boosting tree (GBT) is a machine learning method that combines a decision tree and boosting. The tree is fitted based on the residuals of the previous layer, allowing each layer to fit and modify the data of the previous layer. Boosting improves the accuracy of the prediction function and minimizes the sum of the prediction errors by repeatedly applying the prediction function in series and combining the output of each function with weights. The GBT prediction algorithm has seen growth in its adoption as a scientific method because it is an effective, robust, non-overfitting tool for prediction, and it has gained attention for multiple applications. GBT analysis has been used for the estimation of skin permeation parameters such as K_p (Ita and Prinze, 2024). In addition, Yuan et al. (2023) and Abdallah et al. (2024) reported the usefulness of GBT to predict parameters.

In the present study, GBT was used to predict skin permeation parameters, not the percentage of dermal absorption. This approach may be available for the estimation of skin permeation when chemicals are applied through different thicknesses of the SC. Furthermore, when the percentage of dermal absorption was over- or underestimated, verification of the reason is possible based on the estimated skin permeation parameters. An in silico model for predicting skin permeation parameters was developed with the database reported by Ellison et al., and skin permeation profiles of chemicals in finite dose conditions were estimated using the developed in silico model.

MATERIALS AND METHODS

Materials

Benzoic acid (BA), ethyl p-hydroxybenzoate (EP), methyl p-aminobenzoate (M-PABA), ethyl p-aminobenzoate (E-PABA), and diisopropyl fluorophosphate (DFP) as an ester inhibitor were purchased from FUJIFILM Wako Pure Chemical Co. (Osaka, Japan). Isosorbide dinitrate (ISDN), was purchased from Chugokukayaku Co., Inc. (Hiroshima, Japan). 4-Buthylresorcinol (BR), methyl p-hydroxybenzoate (MP), propyl p-hydroxybenzoate (PP), butyl p-hydroxybenzoate (BP), propyl p-aminobenzoate (P-PABA), butyl p-aminobenzoate (B-PABA) were gifts from Toko Pharmaceutical Industries Co., Ltd. (Tokyo, Japan). Other reagents were of commercially available special grade or for high-performance liquid chromatography (HPLC) and were used without further purification. Table 1 summarizes the chemical names, their abbreviations, and physicochemical properties of their MW and calculated logarithm of K_o/w (clogP). The clogP values were obtained from Scigress (Fujitsu, Kawasaki, Japan), which is a computer simulation software. Excised human skin (obtained from the abdomens of White female donors, aged 44–59) were purchased from K.A.C. Corporation (Kyoto, Japan). The skin samples were stored at -30°C until the in vitro skin permeation experiment.

Table 1. Chemicals used in this experiment for in vitro skin permeation test.

Compound	Abbreviation	MW	logKo/w
Isosorbide dinitrate	ISDN	236.14	1.3
Methyl p-aminobenzoate	M-PABA	151.16	1.4
Ethyl p-aminobenzoate	E-PABA	165.19	1.86
Benzoic acid	BA	122.12	1.9
Methylparaben	MP	152.15	1.96
4-Butylresorcinol	BR	166.22	2.4
Propyl p-aminobenzoate	P-PABA	179.22	2.43
Ethylparaben	EP	166.17	2.47
Butyl p-aminobenzoate	B-PABA	193.24	2.87
Propylparaben	PP	180.2	3.04
Butylparaben	BP	194.23	3.57

Preparation of applied chemical solutions

Applied chemical solutions of BA, BR, and FP were prepared in phosphate-buffered saline (PBS, pH 7.4) and those of MP, EP, PP, BP, M-PABA, E-PABA, E-PABA, and B-PABA were prepared with PBS containing DFP at a concentration of 0.54 µmol/mL. The applied concentration of chemicals was set to 60% of the saturated concentration for each chemical. The pH and concentration of the applied solution was selected according to the procedure of Ellison et al. (2020), because we reference the reported in vitro skin permeation profiles through an excised human skin model in finite dose conditions. The saturated concentration of the chemical solution was obtained by adding an excess amount of chemicals to purified water or PBS containing DFP and stirring for 24 hr at 32°C. The resulting suspension was then filtered using a syringe filter (DISMIC13HP 0.20 µm, ADVANTEC, Tokyo, Japan). The filtrate was diluted with purified water or PBS containing DFP, and the concentrations were measured using HPLC to investigate the saturated concentration. The HPLC conditions were described in the Detection of chemical concentrations section.

In vitro skin permeation test of 11 chemicals with excised human skin

Before in vitro skin permeation experiments, frozen excised human skin samples were thawed at room temperature. Then, the skin samples were set in Franz-type diffusion cells with an effective permeation area of 1.77 cm² to determine the cumulative amount of chemical permeation through the skin. PBS (1.0 mL) was applied on the SC side (donor side) for 1 hr to hydrate the skin, then skin impedance was measured using a low frequency AC impedance meter (AschJapan Co. Ltd., Tokyo, Japan, applied voltage: 5 V, frequency: 12 Hz, measuring range: 0.5 to 100 kΩ) to check the integrity of the skin barrier. Excised human skin with an impedance above 10 kΩ was used in the experiment. The applied saline was removed, and the chemical solution (17.7 µL; equating to 10 µL/cm²) and 6.0 mL of PBS were applied on the donor and receiver sides, respectively. In the case of skin permeation experiments with ester compounds, DFP in PBS was selected as the receiver fluid.

During skin permeation tests, the receiver solution was stirred using a magnetic stirrer. An aliquot (500 µL) was withdrawn from the receiver chamber and the same volume of PBS containing DFP or PBS alone was added to the chamber to keep the volume constant. The chemical concentration in the receiver fluid was determined by HPLC.

Detection of chemical concentrations

Chemical concentrations in the receiver solution were determined using an HPLC system (Prominence, Shimadzu, Kyoto, Japan) equipped with a UV detector (SPA-20A; Shimadzu). The collected receiver solution sample was mixed with acetonitrile in a 2:1 volume ratio and centrifuged at 21,500 × g at 4°C for 5 min to remove proteins and contaminants. The obtained supernatant was injected into the HPLC system. The HPLC system consisted of a pump (LC-20AD; Shimadzu), UV detector (SPD-20A; Shimadzu), system controller (SCL-10AVp; Shimadzu), and auto injector (SIL-20A; Shimadzu). The Inertsil^® ODS (4.6 × 150 mm, 3.5 μm; GL Sciences Inc., Tokyo, Japan) column was maintained at 40°C during elution of the mobile phase at 1.0 mL/min rate, 0.1% phosphoric acid: acetonitrile = 55: 45 for BA, BP, B-PABA, PP, and P-PABA, 0.1% phosphoric acid: acetonitrile = 70: 30 for EP, E-PABA, MP, and M-PABA, purified water: acetonitrile =35: 65 for BR, purified water: acetonitrile =55: 45 for ISDN. The flow rate was maintained at 1.0 mL/min. Peaks were detected at 254 nm for BA, at 260 nm for BP, EP, MP, and PP, at 280 nm for B-PABA, E-PABA, M-PABA, and P-PABA, at 210 nm for BR, and at 220 nm for ISDN.

Measurement of applied skin area

The formulation distribution on the human skin surface mounted in a Franz-type diffusion cell was observed using a microscope (photo resolution: 2592 × 1944 pixels) that can be connected to a computer with a USB cable, and the area was analyzed from the obtained image using ImageJ.

Evaporation rate

The evaporation rate of the donor-solvent from the skin surface was measured using weight change of the applied solution. Excised human skin specimens were mounted in Franz diffusion cells as described above. Purified water was applied on the skin surface at the volume of 10 µL/cm² using a micropipette to spread over the entire effective diffusion area (1.77 cm²), then purified water was absorbed using a paper towel (1 cm × 1 cm) from the skin surface 2 min after the application. The weight of the absorbed towel was immediately measured, and the amount of purified water was calculated from the increased paper towel weight. This procedure was repeated every 2 min up to 20 min. This experiment was done in triplicate.

Estimation of skin permeation parameters using a regression model

GBT was performed using JMP Pro (version 17, SAS Institute Inc., Cary, NC, USA). The default parameters of GBT were used. Skin permeation parameters, D_SC, D_VED, K_SC, and K_VED, of 43 chemicals were cited from the published reported by Ellison et al. (2020) that they selected chemicals with 1.2< logarithmic K_o/w (logK_o/w) < 4.5 and MW < 500 Da. In addition, chemical parameters, such as MW, melting point (MP), clogP, dipole moment, heat of formation, electron affinity (EA), ionization energy (IP), number of hydrogen donors and acceptors, absolute hardness, and chemical potential, the highest occupied molecular orbital (HOMO), and lowest unoccupied molecular orbital (LUMO) energies were calculated using Fujitsu SCIGRESS software (version 3.7, Kanagawa, Japan). The LUMO-HOMO gap was also used as the descriptor. Supplemental Table 1 shows molecular descriptors acquired in the present study. The validation of this prediction model was performed with 80% of samples selected at random as the training set, and the remaining 20% of samples was assigned as the test set for the external validation.

The GBT platform produces an additive decision tree model that is based on many smaller decision trees that are constructed in layers. The tree in each layer consists of a small number of splits, typically five or fewer. Each layer is fit using the recursive fitting methodology. The data are recursively divided based on the relationship between the predicted value and the response value, and a decision tree is created. The division algorithm searches for all divisions of the predicted value that best predict the response. For a given tree, the predicted value for an observation in a leaf is the mean of all observations in that leaf (Hastie et al., 2009). The contribution of the explanatory variable was calculated using the error sum of squares (SS) as an indicator of how much it is reduced by the branching. The smaller the SS of each group by branching means the better the branching was conducted, and the explanatory variable used for the branching has a higher contribution (Sall, 2002). The SS of the selected candidates is calculated as follows:

SS_test = SS_parent-(SS_right-SS_left)

where SS in a node is just S² (n-1).

Calculation of skin permeation of chemicals with finite dose

Skin permeation of chemicals was calculated using a two-layered diffusion model composed of SC and VED layers. The two-layered model with Fick’s second law of diffusion was described in our previous report (Arce et al., 2020), using the equation:

Initial and boundary conditions were:

where L_SC and L_VED are the thicknesses of the SC and VED, respectively. C_v is the chemical concentration in the applied formulation. C_v was changed based on the evaporation rate from the applied formulation and the amount of chemical that migrated into the skin. In addition, the cumulative amount of chemical permeated per unit area, Q, is expressed by Eq. (5).

Equations (1) and (2) were expressed with differential equations as:

Change in chemical concentration in the applied formulation over time was expressed with the successive calculation method.

In this study, permeation parameters calculated from GBT were applied to estimate the skin permeation of chemicals after application in finite dose conditions, and the estimation data were compared with the observed data obtained from in vitro skin permeation experiment with excised human skin. For the estimation of skin permeation profiles, the average thickness (at least 4 different points) of the whole skin measured using a thickness gauge was applied. Additionally, the thickness of the SC layer was set to 12.45 µm referenced from the report by Grégoire et al. (2021).

In addition to estimation using the two-layered diffusion model described above, the skin permeation profile of chemicals after finite dose application was also estimated using a TCAT™ model as a module in GastroPlus™ software Version 9.8.3 (Simulation Plus, Lancaster, CA, USA). Chemical structure-data file (SD file), applied chemical concentration, applied volume, and the saturated concentration were inputted. Other parameters such as the skin thickness (ex. SC thickness 12.45 µm) were referenced from the report by Grégoire et al. (2021) and the initial value of the module of dermal/subcutaneous model in GastroPlus. The input parameters, except for SD file, are shown in Supplemental Table 2.

RESULTS

Fig. 1 shows the relationship between molecular descriptors acquired in the present study. The top right region shows correlation values and the bottom left region shows the plot data between the corresponding row and column descriptors. Several descriptors such as LUMO and HOMO showed good relationships to electron affinity, absolute hardness, and chemical potential. Thus, electron affinity, absolute hardness, and chemical potential were removed, and the skin permeation parameters were estimated using MW, clogP, LUMO, HOMO, MP, dipole moment, heat of formation, and numbers of hydrogen donors and acceptors. GBT was applied to our data set to reveal quantitative relationships between molecular descriptors and skin permeation parameters. The predictive ability of GBT was evaluated using 80% of the data set selected at random and 20% of the data set samples as the test set. Fig. 2 shows a typical example for the relationship between the logarithm of observed and estimated values for each skin permeation parameter. The symbols in black and red show the training and test data sets, respectively. A fairly good estimation of skin permeation parameters was achieved because the average R² values and RMSE in the training set were high and appreciably small, respectively (logK_SC: R² 0.776 ± 0.009, RMSE 0.289 ± 0.009, logD_SC: R² 0.776 ± 0.019, RMSE 0.283 ± 0.012, logK_VED: R² 0.912 ± 0.006, RMSE 0.133 ± 0.005, and logD_VED: R² 0.754 ± 0.021, RMSE 0.458 ± 0.019). In addition, except for logD_VED, estimated skin permeation parameters using the test data sets seemed to achieve relatively good estimations. To identify the molecular descriptors that influence skin permeation parameters, the contribution of each parameter was also confirmed (Fig. 3a). clogP values had the maximum contribution to all of the skin permeation parameters. On the other hand, a relatively high influence of MW was only observed in logD_VED, and its contribution to other skin permeation parameters was not high. Supplemental Table 3 shows the skin permeation and chemical parameters used for GBT.

Fig. 1

Scatterplot matrix and color map of the correlations between different descriptors. Bright red indicates a strong correlation, bright blue indicates a negative correlation, and gray shows no relationship. There is gray in the middle, and 10 shades from bright red to bright blue are available, as shown on the right side of the figure. The data are displayed symmetrically from the top left region to the bottom right region, with the top right showing the color map and the bottom left showing the scatter plot.

Fig. 2

A typical result of prediction of skin permeation parameters of logK_SC (a), logD_SC (b), logK_VED (c), and logD_VED (d). The symbols in black and red show the training and test data sets, respectively. The dotted line shows the relationship between observed and estimated values. A fairly good estimation of skin permeation parameters was achieved because the average R² values and RMSE in the validation set were high and appreciably small, respectively (logK_SC: R² 0.776 ± 0.009, RMSE 0.289 ± 0.009, logD_SC: R² 0.776 ± 0.019, RMSE 0.283 ± 0.012, logK_VED: R² 0.912 ± 0.006, RMSE 0.133 ± 0.005, and logD_VED: R² 0.754 ± 0.021, RMSE 0.458 ± 0.019).

Fig. 3

Contribution index of selected descriptors to estimate skin permeation parameters (a) and a typical example of a decision tree in GBT (b). Bar colors; gray bar: clogP, white bar: MW, slashed bar: LUMO, black bar: HOMO. Mean ± S.D. (n=3)

Fig. 3b shows a typical example for the first layer of the decision tree in the model to predict logK_SC, logD_SC, logK_VED, and logD_VED. For logK_SC, if parameters clogP ≤ 3.3 and LUMO ≥ -1.12 were met, the model added -0.00181 to the equation predicting logK_SC. On the other hand, if the parameter clogP < 3.4 was met, the model added -0.004 to the equation predicting logD_SC.

Fig. 4a shows the change in the applied water from the skin surface in in vitro skin permeation conditions. A gradual decrease in applied water was confirmed, and the decrease rate appeared to be constant. A relationship was confirmed from the amount of retained water on the skin surface with time, and the slope was calculated to be -9.6× 10^-2 g/1.77 cm²/hr. In addition, when 17.7 µL of solution (10 µL/cm²) was applied on excised human skin surface mounted in a Franz-type diffusion cell, no uniform distribution of the applied solution was observed over the entire effective permeation area (1.77 cm²) even at time zero after starting the skin permeation experiment, as shown in Fig. 4b. The applied solution existed as several droplets on the skin surface, and the observed area analyzed using ImageJ (version 1.54h) was 0.58 ± 0.014 cm²/1.77 cm ².

Fig. 4

Change in the weight of the applied water on the skin surface (a) and a typical example of distribution of water droplets just after the application of 10 µL/cm² (b).

Fig. 5 shows in vitro skin permeation profiles of another 11 chemicals through excised human skin after finite dose (10 µL/cm²) application. The range of MW and clogP of the applied chemicals were 116–362 and 1.1–4.0, respectively. In addition, skin permeation profiles were estimated using skin permeation parameters obtained from the constructed GBT model in the present study. Supplemental Table 3 shows the estimated skin permeation parameters for the 11 chemicals. Furthermore, profiles calculated with GastroPlus equipped with the TCAT™ module are also shown. PP, B-PABA, EP, P-PABA, and E-PABA showed similar permeation profiles compared with those predicted in our model. On the other hand, M-PABA showed underestimated dermal absorption and BA, MP, BP, ISDN, and BR had overestimated dermal absorption over 8 hr in the model compared with the observed absorption values. In addition, in the model the slopes of the permeation profiles within 2 hr after the application of BA, MP, ISDN, and BR were greater (higher flux) than those observed.

Fig. 5

Percentage cumulative amount of eleven chemicals that permeated through excised human skin. Symbols: ●: observed data, ○: estimated dermal absorption with our model, □: estimated dermal absorption with GastroPlus. Mean ± S.D. (n = 4)

When skin permeation profiles were predicted using GastroPlus equipped with the TCAT™ module, the predicted profiles did not show the same tendencies as our prediction model. For example, the estimated skin permeation profile of ISDN was much slower compared with our predicted model.

Fig. 6 shows the relationship of the dermal absorption rate between logarithms of the observed and predicted values for the 11 chemicals. The dotted lines in black and red show 1:1 relationships and the ratio of estimated to observed dermal absorption rate within a factor of 2, respectively. The upper left and lower right show the chemicals with overestimated and underestimated model results compared with the observed values, respectively. Most of the chemicals plotted within the dotted red lines. However, for our constructed model, M-PABA was located slightly outside near the boundary line of the factor of 2. On the other hand, when the prediction was conducted with GastroPlus equipped with a transdermal module, underestimation was confirmed for ISDN.

Fig. 6

Relationship between the estimated and observed dermal absorption rate of chemicals. The dotted lines in black and red show a 1:1 relationship between observed and estimated values, and corresponding values to a factor of 2, respectively. (a) Estimated with our model in the range of 1–100% dermal absorption, (b) expanded view of (a) for the range 20–200%, (c) estimated with GastroPlus for the range 1–100% dermal absorption, and (d) expanded view of (c) for the range 20–200%.

DISCUSSION

In silico approaches to estimate skin permeation coefficients of chemicals, K_p, have been developed for many years. This approach may be useful to focus on steady-state flux or the maximum permeation flux obtained from the product of K_p and the applied concentration. Therefore, it can be used not only in the development of chemicals with anticipated pharmaceutical effects after topical application but also for the risk assessment of topically exposed chemicals. Buist et al. (2010) proposed a simple model to predict dermal absorption in finite dose conditions with K_p and K_SC obtained from infinite dose data. However, changes in the applied chemical concentration in the formulation on the skin surface caused by the permeation of chemicals through the SC layer and evaporation of solvent from the formulation occurred in practical use. Kunita et al. (2024) reported a model that predicts dermal absorption after finite dose application by taking into account solvent evaporation from the application site as well as the application of permeation parameters obtained from an artificial skin membrane. However, it may be preferable to construct in silico models to predict the dermal absorption rate of a chemical when applied in finite dose conditions. In the present study, evaporation rate was measured with purified water instead of each applied solution. When water evaporation rate was measured with several chemicals used in the present study at 60% saturated concentration, no significant difference was confirmed compared with purified water. However, because the evaporation rate of solvent may affect solutes and additives in the applied formulation, it is also necessary to construct a model to estimate the evaporation rate of the volatile components in an applied formulation from the skin surface in finite dose conditions.

Trans-epidermal (via the SC) and trans-appendage (via hair follicles and sweat ducts) are both well-known penetration and permeation routes of chemicals through the skin (Kováčik et al., 2020; Hadgraft and Lane, 2016). The lipophilicity of chemicals affects their penetration and permeation routes. The two-layer diffusion model used in the present study, which consisted of homogeneous layers of SC and VED connected in parallel may be preferable when applied to chemicals for which the predominant penetration and permeation routes are via the trans-epidermal route. In our previous report (Hatanaka et al., 2015), the skin concentrations of chemicals at a steady state were estimated well for lipophilic chemicals (clogP > 1) applied in infinite dose conditions, whereas there was underestimation of the skin concentration when hydrophilic chemicals (clogP < 1) were applied due to their main permeation route being via the trans-appendage route. Therefore, in the present study, skin permeation parameters were estimated for chemicals with clogP > 1.

Machine learning to estimate K_p values has been developed by many researchers. Ita and Prinze (2024) investigated the prediction of K_p for 175 compounds using random forest and GBT models. They selected solute descriptors such as excess molar refraction, combined dipolarity/polarizability, overall solute hydrogen bond acidity and basicity, and the MacGowan’s characteristic molecular volume as descriptors, and regressions established a statistically significant association between these descriptors and K_p values. Baba et al. (2015) constructed a model to predict K_p values with the support of vector regression and random forest. Furthermore, molecular dynamics simulations have been applied to predict K_SC (Piasentin et al., 2023) and skin permeability enhanced by penetration enhancers (Gupta et al., 2019). The boosting tree algorithm is generally robust to the overfitting problem (Basant et al., 2016). According to evaluation with a gradient boosting tree with several descriptors, the chemical properties MW, clogP, LUMO, and HOMO were selected in the present study. Potts and Guy's model for K_p was defined with MW and logK_o/w (Potts and Guy, 1992). In addition, Lian et al. (2008) also reported a K_p predicting model with MW and logK_o/w. Furthermore, Cronin’s model is composed of not only MW and clogP but also topological and electronic parameters as descriptors (Cronin et al., 1999). LUMO and HOMO are related to the electronic properties and structure–property relationships of chemicals (Zaugg et al., 2019; Manrique et al., 2024). The descriptors MW and clogP affect drug diffusivity and the amount of drug transferred into the skin. In addition, LUMO and HOMO are related to the reactivity of chemicals (Adindu et al., 2023), which is related to chemical bonding in the skin. Chemical bounding in the skin may be related to chemical diffusivity in the skin. Therefore, the result showed that dermal absorption through skin is governed by molecular size, hydrophobicity, and structural properties, which also appear to be appropriate as descriptors compared with other in silico models. Golbraikh and Tropsha (2002) reported that the predictability of a machine learning model with internal cross-validation offers a limited assessment of predictability. In the present experiment, due to a limitation of the data set of skin permeation parameters, an external validation was not conducted with the coverage across the entire chemical space of the data set. Therefore, further experiments with external validation data sets may be necessary to predict skin penetration parameters because the predictability result of this study may be optimistic.

The wettability of an applied solution on the skin affects the effective penetration or permeation area of chemicals. Estimated dermal absorption for most of the chemicals applied in the present study was found to be higher than the observed data (data not shown) when the applied solution was assumed to be homogeneously spread on the skin surface. This may be due to an overestimation of dermal absorption area compared with the observed condition. There was no significant difference in the contact area among chemical solutions, so the contact area of the applied chemical solution on the skin surface was set to 0.58 cm² to estimate dermal absorption.

Each skin penetration parameter was intentionally changed by a factor of 50 or 1/50 to clarify the effect of each parameter on the skin penetration profile of M-PABA (Supplemental Fig. 1). The trigger of skin penetration of chemicals is distribution or partition of chemicals into the SC layer. The partition is an equilibrium phenomenon, which occurs in a very short period. The amount of chemicals migrated into the surface SC layer can be expressed with K_sc×C_v. Then, the chemicals distributed on the surface of the SC layer diffuses through the SC layer. The chemical distribution process may continue while the solvent remains on the skin surface in the case of finite dose conditions. A decrease in K_sc may affect the amount of chemical distributed, whereas a decrease in D_sc may influence transfer of the chemical into the skin. Thus, the decrease in D_sc provides a slower skin permeation profile of the applied chemical. In addition, a decrease in D_sc influences the amount of chemical that permeates thought the skin due to the slower decrease in the amount of the chemical in the surface SC layer. Decrease in K_ved and D_ved also leads to a decrease in drug migration from the applied solution into the skin because they provide a slower decrease in the amount of chemical in the SC layer. In addition, a decrease in D_ved affects the shape of the skin permeation profile.

For the dermal absorption estimated for M-PABA, underestimation was observed in our model. Here, the estimated skin permeation profile for M-PABA was intentionally fitted to the observed curve by changing the skin permeation parameters. When K_SC and D_VED of approximately 13-fold and 0.1-fold, respectively, were applied, the estimated skin permeation profile fitted well with the observed curve. This result suggested the underestimation and overestimation of partitioning and diffusivity, respectively. Intarakumhaeng et al. (2018) reported an assumed mechanism of solute permeation in finite dose conditions: solute deposition on the skin surface, solute dissolution in skin lipids, and then solute diffusion across the SC. Therefore, chemical dissolution in skin lipids should be taken into consideration as a reason for the underestimation, but these chemicals have moderately lipophilic natures (clogP > 1). Sakata et al. (2014) reported enhanced chemical penetration caused by “solvent drag” or “plugging” effects of the formulation when the formulation was applied as a finite dose application. They selected lipophilic chemicals and oil solvents, but the same tendency might also be found with the combination of hydrophilic chemicals and aqueous solvents. This enhancement of chemical partitioning into the SC layer due to a “drag solvent effect” may be clearly observed in highly soluble chemicals, because their partitioning into the SC layer was generally low due to their hydrophilicity. Furthermore, regarding diffusivity, the data set of chemical diffusivities in the skin used for the development of this in silico model was obtained with infinite dose application. The hydration state by water in the SC layer during the skin permeation study was maintained with infinite dose application. On the other hand, water evaporation from the formulation as well as the SC layer occurred in finite dose conditions. Unlike infinite dose conditions, skin is partially hydrated in finite dose conditions, so it might be necessary to take this into account in the contribution of chemical permeation involving polar pathways in the SC layer. In addition, prolonged hydration affected the barrier function of the SC layer (Tan et al., 2010; Bouwstra et al., 2003), and changes in barrier function of the SC layer especially affected the skin permeation of hydrophilic chemicals. Iikura et al. (2019) reported effects of humidification on the skin permeation of chemicals with different lipophilicities. Skin permeation enhancement effects with hydration were observed when methylparaben (logK_o/w: 1.93) was applied. Yousef et al. (2017) reported that the hydration effect on the permeation of salicylate esters, and the enhancement effect due to increases in diffusivity were observed when glycol salicylate (logK_o/w: 1.63) was applied. Therefore, chemical diffusivity in the SC layer obtained from infinite dose conditions might be higher than that observed in finite dose conditions.

In the present study, a good prediction was defined within a factor of 2. The range of prediction accuracy might be discussed, but other reports have applied twofold of the observed as an index of “a good prediction” (Kunita et al., 2024; Davies et al., 2020), which we also applied. In the present study, chemicals with MW in the range of 116–362 and clogP in the range of 1.1–4.5 were used to develop the in silico model. Therefore, prediction accuracy will be lower when estimation is performed with chemicals outside of the chemical space. When the dermal absorption of thioglycolic acid (MW: 92.1, clogP: 0.1), which is outside of the data set range, was estimated with our model, the value was almost 100% (data not shown), whereas the observed value was about 10% over 24 hr (Nitu Verma et al., 2023).

Predicted skin permeation profiles with GastroPlus incorporating the TCAT™ model were slower than those obtained from our model. When skin permeation parameters obtained from our model were compared with those obtained from GastroPlus incorporating the TCAT™ model, about 10-fold higher K_sc and about 100-fold lower D_sc values in GastroPlus incorporating the TCAT™ model were mainly confirmed as the difference with our model. These differences are related to a slower permeation profile of the applied chemicals as shown in in Supplemental Fig. 1. Our model may predict the skin permeation profile at pH 7.4. This is the limitation of our model, whereas the GastroPlus TCAT™ model is able generate predictions in various applied pH conditions.

Further studies should be performed with an expanded learning data set of skin permeation parameters. This will lead to the development of a better model that will enable the prediction of the amount of chemicals that penetrate the skin in finite dose conditions. In addition, for evaluation of the in silico model, chemicals with relatively better dermal absorption exceeding 35% were used. Thus, prediction should be evaluated with chemicals that exhibit low dermal absorption rates.

Conclusion

In the present study, an in silico model to predict skin permeation parameters, such as K_SC, D_SC, K_VED, and D_VED, was developed using a data set obtained from excised human skin, and the dermal absorption of chemicals was estimated when they were applied in finite dose conditions. Considering the evaporation rate and the area of spread over the skin surface, comparatively good estimates of dermal absorption, except for chemicals with a clogP value close to 1 were observed. A learning data set with in vitro experiments from human skin samples was obtained from infinite dose conditions, so change in diffusivity and partitioning of chemicals in the skin might occur in finite dose conditions. Understanding of the differences in skin permeation parameters between finite and infinite dose conditions may help in the further development of in silico models for dermal absorption. Furthermore, in the future, it may be useful to establish methods to improve prediction accuracy by combining multiple in silico models. Finally, we hope that this research result will lead to further understanding and use of research based on the 3Rs (reduction, refinement, replacement) approach.

ACKNOWLEDGMENT

This study was supported in part by a grant of Long-range Research Initiative (LRI) by Japan Chemical Industry Association (JCIA) (project #: 21-3-01).

Conflict of interest

The authors declare that there is no conflict of interest.

REFERENCES

•  Abdallah,  R.M.,  Hasan,  H.E. and  Hammad,  A. (2024): Predictive modeling of skin permeability for molecules: investigating FDA-approved drug permeability with various AI algorithms. PLOS Digit. Health, 3, e0000483.
•  Adindu,  E.A.,  Godfrey,  O.C.,  Agwupuye,  E.I.,  Ekpong,  B.O.,  Agurokpon,  D.C.,  Ogbodo,  S.E.,  Benjamin,  I. and  Louis,  H. (2023): Structural analysis, reactivity descriptors (homo-LUMO, ELF, NBO), effect of polar (DMSO, ETOH, H₂O) solvation, and libido-enhancing potential of resveratrol by molecular docking. Chemical Physics Impact., 7, 100296.
•  Arce,  F. Jr.,  Asano,  N.,  See,  G.L.,  Oshizaka,  T.,  Itakura,  S.,  Todo,  H. and  Sugibayashi,  K. (2020): Prediction of skin permeation and concentration of rhododendrol applied as finite dose from complex cosmetic vehicles. Int. J. Pharm., 578, 119186.
•  Baba,  H.,  Takahara,  J.,  Yamashita,  F. and  Hashida,  M. (2015): Modeling and prediction of solvent effect on human skin permeability using support vector regression and random forest. Pharm. Res., 32, 3604-3617.
•  Barratt,  M.D. (1995): Quantitative structure-activity relationships for skin permeability. Toxicol. In Vitro, 9, 27-37.
•  Basant,  N.,  Gupta,  S. and  Singh,  K.P. (2016): Predicting human intestinal absorption of diverse chemicals using ensemble learning based QSAR modeling approaches. Comput. Biol. Chem., 61, 178-196.
•  Beal,  M.A.,  Gagne,  M.,  Kulkarni,  S.A.,  Patlewicz,  G.,  Thomas,  R.S. and  Barton-Maclaren,  T.S. (2022): Implementing in vitro bioactivity data to modernize priority setting of chemical inventories. Altern. Anim. Exp., 39, 123-139.
•  Bouwstra,  J.A.,  de Graaff,  A.,  Gooris,  G.S.,  Nijsse,  J.,  Wiechers,  J.W. and  van Aelst,  A.C. (2003): Water distribution and related morphology in human stratum corneum at different hydration levels. J. Invest. Dermatol., 120, 750-758.
•  Bouzom,  F.,  Ball,  K.,  Perdaems,  N. and  Walther,  B. (2012): Physiologically based pharmacokinetic (PBPK) modelling tools: how to fit with our needs? Biopharm. Drug Dispos., 33, 55-71.
•  Buist,  H.E.,  van Burgsteden,  J.A.,  Freidig,  A.P. and  Maas de Sandt,  J.J. (2010): New in vitro dermal ab and the prediction of dermal absorption under for risk assessment purposes. Regul. Toxicol. Pharmacol., 57, 200-209.
•  Cronin,  M.T.,  Dearden,  J.C.,  Moss,  G.P. and  Murray-Dickson,  G. (1999): Investigation of the mechanism of flux across human skin in vitro by quantitative structure-permeability relationships. Eur. J. Pharm. Sci., 7, 325-330.
•  Davies,  M.,  Jones,  R.D.,  Grime,  K.,  Jansson-Löfmark,  R.,  Fretland,  A.J.,  Winiwarter,  S.,  Morgan,  P. and  McGinnity,  D.F. (2020): Improving the accuracy of predicted human pharmacokinetics: lessons learned from the AstraZeneca drug pipeline over two decades. Trends Pharmacol. Sci., 41, 390-408.
•  Ellison,  C.A.,  Tankersley,  K.O.,  Obringer,  C.M.,  Carr,  G.J.,  Manwaring,  J.,  Rothe,  H.,  Duplan,  H.,  Géniès,  C.,  Grégoire,  S.,  Hewitt,  N.J.,  Jamin,  C.J.,  Klaric,  M.,  Lange,  D.,  Rolaki,  A. and  Schepky,  A. (2020): Partition coefficient and diffusion coefficient determinations of 50 compounds in human intact skin, isolated skin layers and isolated stratum corneum lipids. Toxicol. In Vitro, 69, 104990.
•  Ghanem,  A.H.,  Mahmoud,  H.,  Higuchi,  W.I.,  Liu,  P. and  Good,  W.R. (1992): The effects of ethanol on the transport of lipophilic and polar permeants across hairless mouse skin: Methods/validation of a novel approach. Int. J. Pharm., 78, 137-156.
•  Grégoire,  S.,  Sorrell,  I.,  Lange,  D.,  Najjar,  A.,  Schepky,  A., et al. (2021): Cosmetics Europe evaluation of 6 in silico skin penetration models. Comput. Toxicol., 19.
•  Golbraikh,  A. and  Tropsha,  A. (2002): Beware of q2! J. Mol. Graph. Model., 20, 269-276.
•  Gupta,  R.,  Dwadasi,  B.S.,  Rai,  B. and  Mitragotri,  S. (2019): Effect of chemical permeation enhancers on skin permeability: in silico screening using molecular dynamics simulations. Sci. Rep., 9, 1456.
•  Guy,  R.H. and  Hadgraft,  J. (1985): Transdermal drug delivery: a simplified pharmacokinetic approach. Int. J. Pharm., 24, 267-274.
•  Hadgraft,  J. and  Lane,  M.E. (2016): Advanced topical formulations (ATF). Int. J. Pharm., 514, 52-57.
•  Hartmanshenn,  C.,  Scherholz,  M. and  Androulakis,  I.P. (2016): Physiologically-based pharmacokinetic models: approaches for enabling personalized medicine. J. Pharmacokinet. Pharmacodyn., 43, 481-504.
• Hastie, T.J., Tibshirani, R.J. and Friedman, J.H. (2009). The Elements of Statistical Learning: Data Mining, Inference, and Prediction. 2nd ed. New York: Springer-Verlag.
•  Hatanaka,  T.,  Yoshida,  S.,  Kadhum,  W.R.,  Todo,  H. and  Sugibayashi,  K. (2015): In silico estimation of skin concentration following the dermal exposure to chemicals. Pharm. Res., 32, 3965-3974.
•  Iikura,  H.,  Uchida,  K.,  Ogawa-Fuse,  C.,  Bito,  K.,  Naitou,  S.,  Hosokawa,  M.,  Uchida,  T.,  Sano,  T.,  Todo,  H. and  Sugibayashi,  K. (2019): Effects of temperature and humidity on the skin permeation of hydrophilic and hydrophobic drugs. AAPS PharmSciTech, 20, 264.
•  Intarakumhaeng,  R.,  Wanasathop,  A. and  Li,  S.K. (2018): Effects of solvents on skin absorption of nonvolatile lipophilic and polar solutes under finite dose conditions. Int. J. Pharm., 536, 405-413.
•  Ita,  K. and  Prinze,  J. (2024): Machine learning for skin permeability prediction: random forest and XG boost regression. J. Drug Target., 32, 57-65.
•  Kretsos,  K.,  Miller,  M.A.,  Zamora-Estrada,  G. and  Kasting,  G.B. (2008): Partitioning, diffusivity and clearance of skin permeants in mammalian dermis. Int. J. Pharm., 346, 64-79.
•  Kováčik,  A.,  Kopečná,  M. and  Vávrová,  K. (2020): Permeation enhancers in transdermal drug delivery: benefits and limitations. Expert Opin. Drug Deliv., 17, 145-155.
•  Kunita,  R.,  Nishijima,  T.,  Todo,  H. and  Miyazawa,  M. (2024): Integrating mathematical approaches (IMAS): novel methodology for predicting dermal absorption rates of chemicals under finite dose conditions. J. Toxicol. Sci., 49, 219-230.
•  Lian,  G.,  Chen,  L. and  Han,  L. (2008): An evaluation of mathematical models for predicting skin permeability. J. Pharm. Sci., 97, 584-598.
•  Lien,  E.J. and  Gao,  H. (1995): QSAR analysis of skin permeability of various drugs in man as compared to in vivo and in vitro studies in rodents. Pharm. Res., 12, 583-587.
•  Manrique,  P.D.,  Leus,  I.V.,  López,  C.A.,  Mehla,  J.,  Malloci,  G.,  Gervasoni,  S.,  Vargiu,  A.V.,  Kinthada,  R.K.,  Herndon,  L.,  Hengartner,  N.W.,  Walker,  J.K.,  Rybenkov,  V.V.,  Ruggerone,  P.,  Zgurskaya,  H.I. and  Gnanakaran,  S. (2024): Predicting permeation of compounds across the outer membrane of P. aeruginosa using molecular descriptors. Commun. Chem., 7, 84.
•  Mitragotri,  S.,  Anissimov,  Y.G.,  Bunge,  A.L.,  Frasch,  H.F.,  Guy,  R.H.,  Hadgraft,  J.,  Kasting,  G.B.,  Lane,  M.E. and  Roberts,  M.S. (2011): Mathematical models of skin permeability: an overview. Int. J. Pharm., 418, 115-129.
•  Moss,  G.P. and  Cronin,  M.T. (2002): Quantitative structure-permeability relationships for percutaneous absorption: re-analysis of steroid data. Int. J. Pharm., 238, 105-109.
•  Nitu Verma,  K.G.,  Rakesh,  G. and  Beena,  R. (2023): Multiscale modeling of molecule transport through skin’s deeper layers. Comput. Toxicol., 26.
•  Paul Friedman,  K.,  Gagne,  M.,  Loo,  L.H.,  Karamertzanis,  P.,  Netzeva,  T.,  Sobanski,  T.,  Franzosa,  J.A.,  Richard,  A.M.,  Lougee,  R.R.,  Gissi,  A.,  Lee,  J.J.,  Angrish,  M.,  Dorne,  J.L.,  Foster,  S.,  Raffaele,  K.,  Bahadori,  T.,  Gwinn,  M.R.,  Lambert,  J.,  Whelan,  M.,  Rasenberg,  M.,  Barton-Maclaren,  T. and  Thomas,  R.S. (2020): Utility of in vitro bioactivity as a lower bound estimate of in vivo adverse effect levels and in risk-based prioritization. Toxicol. Sci., 173, 202-225.
•  Piasentin,  N.,  Lian,  G. and  Cai,  Q. (2023): In silico prediction of stratum corneum partition coefficients via cosmomic and molecular dynamics simulations. J. Phys. Chem. B, 127, 2719-2728.
•  Potts,  R.O. and  Guy,  R.H. (1992): Predicting skin permeability. Pharm. Res., 9, 663-669.
•  Sakata,  O.,  Fujii,  M.,  Koizumi,  N.,  Nakade,  M.,  Kameyama,  K. and  Watanabe,  Y. (2014): Effects of the physicochemical properties of oil vehicles on the skin penetration of poor skin-penetrable materials under finite dose conditions in vitro. J Pharm Sci Technol. Jpn., 74, 84-92.
• Sall, J. (2002). “Monte Carlo Calibration of Distributions of Partition Statistics.” SAS Institute Inc., Cary, NC. https://www.jmp.com/content/dam/jmp/documents/en/white-papers/montecarlocal.pdf. (Nov. 27, 2024, accessed)
•  Sugibayashi,  K.,  Todo,  H.,  Oshizaka,  T. and  Owada,  Y. (2010): Mathematical model to predict skin concentration of drugs: toward utilization of silicone membrane to predict skin concentration of drugs as an animal testing alternative. Pharm. Res., 27, 134-142.
•  Tan,  G.,  Xu,  P.,  Lawson,  L.B.,  He,  J.,  Freytag,  L.C.,  Clements,  J.D. and  John,  V.T. (2010): Hydration effects on skin microstructure as probed by high-resolution cryo-scanning electron microscopy and mechanistic implications to enhanced transcutaneous delivery of biomacromolecules. J. Pharm. Sci., 99, 730-740.
•  Wang,  T.F.,  Kasting,  G.B. and  Nitsche,  J.M. (2007): A multiphase microscopic diffusion model for stratum corneum permeability. II. Estimation of physicochemical parameters, and application to a large permeability database. J. Pharm. Sci., 96, 3024-3051.
•  Yuan,  Y.,  Han,  Y.,  Yap,  C.W.,  Kochhar,  J.S.,  Li,  H.,  Xiang,  X. and  Kang,  L. (2023): Prediction of drug permeation through microneedled skin by machine learning. Bioeng. Transl. Med., 8, e10512.
•  Yousef,  S.,  Mohammed,  Y.,  Namjoshi,  S.,  Grice,  J.,  Sakran,  W. and  Roberts,  M. (2017): Mechanistic evaluation of hydration effects on the human epidermal permeation of salicylate esters. AAPS J., 19, 180-190.
•  Zaugg,  K.,  Velasco,  J.,  Robins,  K.A. and  Lee,  D.C. (2019): Understanding the electronic properties of acceptor-acceptor’-acceptor triads. ACS Omega, 4, 5434-5441.
•  Zobl,  W.,  Bitsch,  A.,  Blum,  J.,  Boei,  J.J.,  Capinha,  L.,  Carta,  G.,  Castell,  J.,  Davoli,  E.,  Drake,  C.,  Fisher,  C.P.,  Heldring,  M.M.,  Islam,  B.,  Jennings,  P.,  Leist,  M.,  Pellegrino-Coppola,  D.,  Schimming,  J.P.,  Snijders,  K.E.,  Tolosa,  L.,  van de Water,  B.,  van Vugt-Lussenburg,  B.M.,  Walker,  P.,  Wehr,  M.M.,  Wijaya,  L.S. and  Escher,  S.E. (2024): Protectiveness of NAM-based hazard assessment - which testing scope is required? Altern. Anim. Exp., 41, 302-319.