In silico prediction of dermal absorption of pesticides – an evaluation of selected models against results from in vitro testing

ABSTRACT

Current guidance for the estimation of dermal absorption (DA) of pesticides recommends the use of default values, read-across of information between formulations and in vitro testing. While QSARs exist to estimate percutaneous absorption, their use is currently not encouraged. Therefore, the potential of publicly available models for DA estimation was investigated based on data from 564 human in vitro DA experiments on pesticides. The classic Potts Guy model, the correction of Cleek Bunge for highly lipophilic chemicals, the mechanistic model of Mitragotri, and the COSMOS model were used to estimate the permeability coefficient k_p. Different approaches were explored to calculate the percentage of external dose absorbed. IH SkinPerm was examined as stand-alone model. The models generally failed to accurately predict experimental values. For 30–40% of the predictions, there was overestimation by one order of magnitude. Three models underpredicted >10% of the cases, the remaining models <5%. DA of hydrophilic substances was typically underpredicted. Overprediction was more prominent for solid preparations and suspensions. The molecular weight, irritation potential and skin thickness did not correlate with the models’ predictivity. Of the models investigated, IH SkinPerm performed best with 38% of the predictions within one order of magnitude and 2% underpredicted cases.

KEYWORDS:

• Dermal absorption
• in silico model
• pesticides
• QSAR
• skin

Introduction

The term dermal absorption is used to describe the transport of chemicals from the outer surface of the skin to the systemic circulation [1]. Dermal absorption can be separated into three processes: penetration, which is the entry of the substance into a skin layer (e.g. the stratum corneum), permeation, which is the passage of a substance from one layer of the skin to another structurally and functionally different layer, and resorption, which is the uptake of the substance from the dermis to the bloodstream [2]. By considering the anatomy of the skin, there are three possible ways by which a chemical can reach the deeper skin layers, and thus be available for systemic uptake: (i) the transcellular route, which assumes that the chemicals move through the skin cells by partitioning between the cell membrane and the extracellular area, e.g. partitioning between the corneocytes and the extra-lipid area in the stratum corneum, (ii) the intercellular route, which assumes that the chemicals move inside the extracellular matrix without entering or exiting the cells, e.g, the lipid matrix of the stratum corneum, and (iii) the appendageal route, which assumes that the chemicals reach the deeper skin layers directly through appendages (e.g. hair follicles) [2]. Active transport is not considered to play a role in dermal absorption, meaning that chemicals principally penetrate the outermost skin structure, the lipophilic stratum corneum, by passive diffusion. Subsequently, they diffuse into the viable epidermis and the dermis, areas dominated by a mainly aqueous environment. Chemicals present in the dermis can then be resorbed into the blood flow, therefore becoming systemically available [2]. Even though the hair follicle density is highly variable depending on species, sex, race and skin region etc., the appendagal route can play a major role by providing a shortcut route to the epidermis.

When it comes to occupational risk assessment for pesticides, the skin typically represents the most relevant site of exposure. Dermal exposure to pesticides occurs mainly during mixing and application of the products as well as through re-entry activities, i.e. when professionals enter an area that has been previously treated with pesticides thereby getting exposed to pesticide residues [3]. To ensure that operators coming in contact with a particular pesticide are safe, the systemic exposure to the active substance has to be lower than or equal to the systemic AOEL (Acceptable Operator Exposure Level) in each application scenario examined [4]. According to the current EFSA guidance, in cases where experimental data are not available, the use of default values dependent on physical and chemical properties of the active substance and/or the degree of dilution is recommended [5]. Although currently not recommended by regulators [5], a probably more accurate way to estimate the extent of dermal absorption without the need for experimental data could be by using computational models. These can range from Quantitative Structure-Activity Relationships (QSARs) to pharmacokinetic models or analytical and numerical solutions of Fick’s diffusion laws. The QSARs for predicting skin permeability are usually derived empirically by multivariate linear regression, where the dependent variable is an experimental index of permeation (e.g. permeability coefficient, k_p) [6,7]. Since the lipid-rich stratum corneum is considered the main barrier, descriptors for lipophilicity are expected to be highly influential. Additional parameters may be the molecular size, which may determine whether a molecule is small enough to permeate through the skin, hydrogen bonding descriptors, which may determine interactions with skin protein structures [8], as well as water solubility, which may give a hint about the thermodynamic stability of a substance in its vehicle.

By treating the skin as a pseudo-homogenous membrane that lies between a donor and acceptor phase, compound permeation into the skin can be described mathematically by Fick’s diffusion equation, as the flux J, i.e. the amount M that crosses the membrane per time t per area A:

(1) $\mathrm{f}\mathrm{l}\mathrm{u}\mathrm{x}=\mathrm{J}=-\frac{\mathrm{d}\mathrm{M}}{\mathrm{d}\mathrm{t}\ast \mathrm{A}}=\mathrm{D}\ast \frac{\mathrm{d}\mathrm{C}\left(\mathrm{x},\mathrm{t}\right)}{\mathrm{d}\mathrm{x}}\left[\mathrm{a}\mathrm{m}\mathrm{o}\mathrm{u}\mathrm{n}\mathrm{t}\ast \mathrm{c}{\mathrm{m}}^{-2}\ast {\mathrm{h}}^{-1}\right]$(1)

where D denotes the diffusion coefficient and C the concentration of the permeant at depth x and time t. When the donor fluid volume is large enough to keep the concentration of a solvent in it constant, infinite dose conditions dominate the donor-membrane-acceptor system. Under these circumstances, the steady state can be reached, and Equation (1) can be transformed to:

(2) $\mathrm{m}=\frac{\mathrm{D}\ast \mathrm{A}\ast \mathrm{t}\ast \left({\mathrm{C}}_{\mathrm{d}}-{\mathrm{C}}_{\mathrm{a}}\right)}{\mathrm{l}}$(2)

${\mathrm{C}}_{\mathrm{d}}$ and ${\mathrm{C}}_{\mathrm{a}}$ are the concentrations on the donor and the acceptor site and $\mathrm{l}$ the path length. Under sink conditions ${\mathrm{C}}_{\mathrm{a}}\approx 0$ and ${\mathrm{C}}_{\mathrm{d}}=\mathrm{P}\mathrm{C}\ast {\mathrm{C}}_{\mathrm{v}}$, where $\mathrm{P}\mathrm{C}$ is the partition coefficient of the compound between the stratum corneum and the vehicle (or membrane and donor) and ${\mathrm{C}}_{\mathrm{v}}$ the compound’s concentration in the vehicle. Division of Equation (2) by A *t gives the steady state flux ${\mathrm{J}}_{\mathrm{s}\mathrm{s}}\left[\mathrm{a}\mathrm{m}\mathrm{o}\mathrm{u}\mathrm{n}\mathrm{t}\ast \mathrm{c}{\mathrm{m}}^{-2}\ast {\mathrm{h}}^{-1}\right]$:

(3) ${\mathrm{J}}_{\mathrm{s}\mathrm{s}}=\frac{\mathrm{D}\ast \mathrm{P}\mathrm{C}\ast {\mathrm{C}}_{\mathrm{v}}}{\mathrm{l}}$(3)

Apart from describing the flux, Equation (3) gives information about another useful parameter, the permeability coefficient k_p. The permeability coefficient describes the potency of a compound to permeate from a specific vehicle through a homogenous barrier, provided that this vehicle does not alter the barrier properties. By transforming Equation (3), k_p can be defined as:

(4) $k_p=\frac{{\mathrm{J}}_{\mathrm{s}\mathrm{s}}}{{\mathrm{C}}_{\mathrm{v}}}=\frac{\mathrm{D}\ast \mathrm{P}\mathrm{C}}{\mathrm{l}}\left[\mathrm{c}\mathrm{m}\ast {\mathrm{h}}^{-1}\right]$(4)

An additional parameter is the lag time t_lag which is a time point associated with the steady-state flux and is defined as:

(5) $t_{lag}=\frac{{\mathrm{l}}^2}{6\mathrm{D}}\left[\mathrm{h}\right]$(5)

In general, three major factors influence dermal uptake: the compound and formulation properties, the exposure conditions and the structure of the skin. The lipophilicity of a compound plays a crucial role in its absorption and is usually expressed as the octanol/water partition coefficient log k_ow [3]. Lipophilicity and permeability coefficients are correlated in terms of a convex parabolic function [9,10]. Due to the lipophilic environment of the stratum corneum, lipid-soluble compounds can easily cross this layer. However, as soon as they reach the hydrophilic epidermis, their transport rate will decrease resulting in a gradual yet temporary accumulation within the skin. This process is called reservoir effect [3]. The substances may be released as soon as the concentration gradient is high enough and allows their further migration through the skin layers and into the blood vessels. The lipophilic nature is not the only property that determines a compound’s tendency to build a skin reservoir. High keratin-binding and slow desorption kinetics also facilitate the formation of a reservoir [11]. For example, it has been shown that, among the most commonly used anti-fungal substances, their potential to bind to the protein structures of the stratum corneum is a better indicator of their ability to penetrate the epidermis than their lipid-binding affinity [12]. The percutaneous permeation is also affected by the molecular size, usually expressed as the molecular weight MW. Even though the molecular volume MV is a more accurate descriptor for size, the use of molecular weight is more common, since its estimation is more straightforward. This approximation is acceptable because for most compounds of interest – mainly hydrocarbons – the molecular weight to molecular volume ratio is nearly constant. It has been suggested that molecular weight and penetration rate are inversely related [13,14].

Based on the considerations described above, several mathematical models were developed to predict the absorption of chemicals across the skin. Here, we evaluate the capacity of well-established models and/or their combinations to estimate the dermal absorption of agrochemicals from complex plant protection products and in-use dilutions thereof (see Table 1 for a summary of the selected models). Four of the models typically predict k_p, while in regulatory risk assessment, the percentage absorbed relative to the external dose is required as a descriptor of absorption. Therefore, we explored two approaches for the conversion of k_p into the percentage of external dose absorbed for a given exposure situation, resulting in a total of eight model combinations (four k_p prediction models combined with two k_p-% external dose absorbed prediction models). In addition, one stand-alone model, which directly computes the percentage of dose absorbed was included. We chose to use models developed mainly for single substances and their solutions rather than complex preparations. This was because information on the identity of co-formulants and the detailed composition of commercial pesticide products is treated as confidential business information and was not available as input parameters for this study.

Table 1. Summary of the applied models. VE: Viable epidermis, SC: Stratum corneum, RF: Receptor fluid.

Model	Calculates	Approach	Diffusion through	Comment
Potts Guy [14]	kp	Regression	Lipid pathway	State of the art model, simple and straightforward
Mitragotri [19]	kp	Mechanistic – analytical	Tortuous lipid pathway	Equation was adapted as proposed by Lian [7] resembles that of Potts Guy
Cleek Bunge [20]	kp	Mechanistic	Lipid pathway	Corrects Potts Guy by additionally accounting for the resistance provided by the hydrophilic VE
COSMOS [21]	kp	Regression & statistical weighting	Lipid pathway
TNO [22]	% absorbed	Simplistic transition from kp to mass absorbed	n/a
Buist [23]	% absorbed in SC and RF separately	Analytical, more elaborate transition from kp to mass absorbed	n/a	Adaptable for finite doses, lag time and distribution of chemical in the stratum corneum are accounted for
IH SkinPerm [24]	% evaporated & % absorbed in SC and VE, separately	Regression – analytical	Lipid & transcellular pathway	Physical and chemical parameters such as density and vapour pressure are accounted for

Materials and methods

In vitro dermal absorption dataset for model implementation

We used a database compiled from original study reports on dermal absorption submitted to the German Federal Institute for Risk Assessment within the regulatory framework for approval of pesticides. Only data from human in vitro dermal absorption experiments meeting the following quality criteria was included for the purpose of this analysis: (1) the sampling time did not exceed 24 hours, (2) the recovery was between 90–110%, (3) the experiments were conducted according to the OECD Guidelines for the testing of chemicals [15] and (4) the facility that performed the experiments complied with the principles of Good Laboratory Practice. The data was also evaluated against the recommendations of the OECD GD 28 [1] and OECD GN 156 [16] as well as the respective EFSA Guidance [17] and experiments with major deviations were not included. The resulting dataset consisted of measurement means from 564 in vitro experiments based on 134 different active substances present in 224 different mixtures or dilutions thereof and has recently been published by EFSA [18]. Data related to commercial or experimental plant protection product formulations except for ten cases, where solutions in water, ethanol or acetonitrile were tested. Parameters recorded for each experiment included information on the dose applied (i.e. concentration and exposed area), exposure, observation and sampling times, skin type, occlusion state, and formulation type and physical state of the applied specimens. The dataset also included information on the molecular weight of the active substances. More information on the design of individual studies (e.g. flow-through vs. static cells, skin type and region, etc.) is available from the dataset [18]. Additional physicochemical parameters that were required for the implementation of the selected models such as lipophilicity, water solubility, and vapour pressure were collected from peer-reviewed regulatory documents that are publicly available, other public access literature and databases. Information on the experimental octanol/water partition coefficient log k_ow and the water solubility was either retrieved from the EFSA conclusions on the peer review of the pesticide active substances, the European Commission pesticide database, the US EPA, FAO, or PubChem. For log k_ow, values measured at pH 5.5 were available only for Captan, Fluroxypyr and Acetochlor. For water solubility, values measured at 20°C were preferred per OECD Testing Guideline 105 [25]; solubility values measured at 25°C were also accepted. If no solubility values were available at pH 7, known values for the pH closest to pH 7 were used, provided that the pH does not significantly influence water solubility of the active substance or that the typical formulation’s pH is the same to that of the measurement. In cases where active substances consist of two isomers, the overall water solubility was calculated from the solubilities of both isomers after accounting for the isomer ratio. Wherever the water solubility is given as S > X mg/l the lowest known value was used (e.g. S = X mg/l). In rare cases where no experimental values were retrievable, the parameters were computed with Episuite^TM (WSKOWWIN for Diniconazole M (R)-Isomer). Vapour pressure values were obtained either from the EFSA conclusions on the peer review of the pesticide risk assessment of each active substance or from the European Commission pesticide database. Values measured at 25°C were preferred; values at 20°C and extrapolated from higher temperatures were also accepted. For Triticonazole and Asulam, only values at 50°C and 45°C were available. In cases where vapour pressure values were not found, they were computed with EPISUITE^TM (Azafenidin, Diniconazole-M (R)-Isomer, Fenamidone).

Applied models

Three types of dermal absorption models were employed: Four permeability coefficient (k_p) predicting models: Potts Guy, Mitragotri, Cleek Bunge, COSMOS; two k_p -% amount absorbed transition formulae: TNO, Buist; and one direct % amount absorbed predicting model: IH SkinPerm [14,19–24]. Since the aim of this study was to examine the suitability of the models to derive the percentage of the external dose absorbed, the k_p values predicted by the permeability coefficient models Potts Guy, Mitragotri, Cleek Bunge, COSMOS were fed into the two different k_p-% absorbed transition formulae. This two-step calculation yielded a total of eight different model-combinations. An overview of all the applied models and their underlying assumptions is provided in Table 1.

The permeability coefficients were estimated according to the models of Potts Guy, Mitragotri, and Cleek Bunge, by using the following equations:

For Potts Guy [14]:

(6) $log{k}_p=-6.3+0.71\ast log{k}_{ow}-0.0061\ast \mathrm{M}\mathrm{W}\left[cm\ast s^{-1}\right]$(6)

For Mitragotri [7,19]:

(7) $log{k}_p=-5.2+0.7\ast log{k}_{ow}-0.072\ast \mathrm{M}{\mathrm{W}}^{2/3}\left[cm\ast s^{-1}\right]$(7)

For Cleek Bunge correction [20], for highly lipophilic chemicals where B ≥ 0.1 and exposure time long enough to achieve steady state:

(8) $k_p=\frac{k_{p\left(sc\right)}}{1+B}$(8)

With $B=\frac{k_{p\left(sc\right)}\ast \sqrt{MW}}{2.6cm/h}$ and $k_{p\left(sc\right)}$, the permeability coefficient for aqueous vehicles calculated with the Potts and Guy equation (Equation 6).

The permeability coefficients from the COSMOS model were calculated automatically after input of the molecular structure of the substance of interest in the web-based application (https://knimewebportal.cosmostox.eu/). The molecular structures of the active substances were obtained from the EFSA conclusion on the peer review of the pesticide risk assessment of each active substance or, if not available there, from PubChem.

For the approximation of the % of the external dose absorbed from the substance-specific permeability coefficients k_p, further parameters including exposure time, area and concentration were taken into account. The following equation was suggested by TNO [22]:

(9) $A_T\left[\mathrm{\%}\right]=\frac{k_p\left[\frac{cm}{h}\right]\ast time\left[h\right]\ast area\left[c{m}^2\right]\ast concentrationinvehicle\left[\frac{mg}{c{m}^3}\right]}{totaldoseapplied\left[mg\right]}\ast 100$(9)

Time was taken as the experimental sample time, namely the time of exposure to the preparation plus the observation time after the exposure. This was made in order to account for the reservoir effect, i.e. the depot of the substance within the upper layers of the skin.

The Buist formula provides a more elaborate transition from the k_p to the percentage of the amount absorbed by adapting infinite dose experimental data such as lag time t_lag and permeability coefficient k_p to finite dose conditions [23]. The resulting set of equations is summarized in the supplementary document. For the implementation of this model, time t was taken as the experimental exposure time. Since this model accounts for the amount retained in the stratum corneum after exposure, there is no need to take time as sample time to account for reservoir effects in this case.

IH SkinPerm version 1.21 was used to estimate the amount of active substance present in the stratum corneum and the viable epidermis as a whole (https://www.aiha.org/get-involved/VolunteerGroups/Documents/IHSkinPerm_V1-21.xlsm). The input requires information on the active substances’ physical and chemical properties including molecular weight, water solubility, octanol/water partition coefficient, vapour pressure, and density, as well as information on the experimental settings. The developers recommend the use of values measured at pH 5.5 (skin pH) for lipophilicity and water solubility. The log k_ow values measured at pH 5.5 were used only for few substances, due to low data availability; for all other substances, log k_ow values at pH 7 were used. For water solubility, any value measured close to pH 5.5 was accepted, provided that there is no significant influence of the pH on the solubility, or that the typical formulation’s pH is the same to that of the measurement. If no value close to pH 5.5 was available, the solubility at pH 7 was used instead. Density values were approximated as 1000 mg/cm³, as after randomized examination of a limited number of study reports, it was concluded that most reported density values do not substantially differ from 1000 mg/cm³. The experimental input required for the simulation includes the absolute amount applied to the skin surface, the area of the skin exposed to the product and the experimental exposure time. The maximum skin adherence was left as default (−1 mg/cm²) for liquids and was adjusted to 1 mg/cm² in the case of solids (according to the developers any value between 0 and 2 mg/cm² is accepted). The thickness of stagnant air is recommended to be set at 1 cm for bare skin; however, it was left as default at 3 cm to simulate (semi-)occlusion.

Data analysis

In accordance with current risk assessment practice, the experimental dermal absorption was taken as the sum of the relative amount present in the receptor fluid, the skin sample and the tape strips – except for the first and second tapes [5]. Comparison to the predicted values was made by first exploring the general concordance for each tested model. A prediction factor was introduced, which is effectively the ratio between the predicted and experimental absorption values (10).

(10) $\mathrm{P}\mathrm{r}\mathrm{e}\mathrm{d}\mathrm{i}\mathrm{c}\mathrm{t}\mathrm{i}\mathrm{o}\mathrm{n}\mathrm{f}\mathrm{a}\mathrm{c}\mathrm{t}\mathrm{o}\mathrm{r}\left(\mathrm{P}\mathrm{F}\right)=\mathrm{p}\mathrm{r}\mathrm{e}\mathrm{d}\mathrm{i}\mathrm{c}\mathrm{t}\mathrm{e}\mathrm{d}\mathrm{a}\mathrm{b}\mathrm{s}\mathrm{o}\mathrm{r}\mathrm{p}\mathrm{t}\mathrm{i}\mathrm{o}\mathrm{n}\left[\mathrm{\%}\right]/\mathrm{e}\mathrm{x}\mathrm{p}\mathrm{e}\mathrm{r}\mathrm{i}\mathrm{m}\mathrm{e}\mathrm{n}\mathrm{t}\mathrm{a}\mathrm{l}\mathrm{a}\mathrm{b}\mathrm{s}\mathrm{o}\mathrm{r}\mathrm{p}\mathrm{t}\mathrm{i}\mathrm{o}\mathrm{n}\left[\mathrm{\%}\right]$(10)

For regression analysis, logarithmic transformation was applied when appropriate. The models were ranked according to the following criteria:

1. Low percentage of underestimated cases (log PF<0, predicted values lower than experimental)

2. Low factor of underestimation

3. High percentage of cases estimated within one order of magnitude (0≤ log PF <1,)

4. Low percentage of highly overestimated cases (log PF ≥3, predicted values more than 100-fold over the experimental)

Factors potentially influencing the predicted to experimental concordance were explored. Data on the skin irritancy of the active substances were collected from the EFSA conclusions on the peer review of the pesticide risk assessment, the European Commission pesticide database or the ECHA database. The irritation/corrosion potential of the tested formulation was inferred according to the ECHA Guidance on the Application of the CLP Criteria [26] as described in the supplement. Preparations were assigned to four categories: non-irritants (90.4% of the sample), slight irritants as indicated by the EFSA peer-reviewed conclusions (4% of the sample), irritants (CLP category 2) (4.6% of the sample), corrosives (CLP category 1C) (0.4% of the sample). Cases, for which no data were found, were excluded from the analysis (0.6% of the sample). For analysis of the influence of formulation type, the following five categories were formed and considered for further analysis: suspensions (37.6% of the sample), emulsions (23.8% of the sample), dispersions (15.1% of the sample), solutions (8.7% of the sample), and ‘other’ (14.8% of the sample, cases with no data or not assignable). The physical state was defined by the physical state of the applied preparation in each experiment, as indicated by the study reports. When a solid product was moistened with artificial sweat to mimic exposure conditions it was marked as ‘paste’ and not as ‘solid.’ Following categories were formed: liquids (94.5% of the sample), pastes (2.7% of the sample), and solids (2.3% of the sample). Cases, where no data were found, were excluded from the analysis (0.5% of the sample). Information on the skin thickness/type was provided in the database. Cases, where the entries were marked as ‘scissors’ or ‘full thickness’, were treated as ‘dermatomed’ after confirming from the original study reports that the skin thickness was lower than 1 mm. The resulting ‘dermatomed’ category comprised cases with skin thickness varying from 200 μm to 800 μm.

Non-parametric testing was performed to compare the distributions of the logPF for the above-described categories within every model. Kruskal Wallis and Jonckheere-Terpstra non-parametric tests were performed pairwise with the alpha values generally set to 0.05. Bonferroni correction for multiple pairwise comparisons was applied where necessary, in order to reduce type I error (rejection of the null hypothesis).

Software

Predictions were performed with Microsoft Excel, except for the COSMOS model which is available as a web-based application (https://knimewebportal.cosmostox.eu/). For the implementation of IH SkinPerm, the xlsm version 1.21 was used (https://www.aiha.org/get-involved/VolunteerGroups/Documents/IHSkinPerm_V1-21.xlsm). The results were plotted and analysed with IBM SPSS Statistics version 23. The figures were edited with Inkscape version 0.92.3. and GIMP version 2.10.10.

Results

General performance of the models

Figure 1 provides an overview of the log predicted absorption of each model plotted against the respective log experimental absorption. Poor relationships are apparent for all models. To investigate potential relationships between the predicted and experimental absorption values more systematically, we performed a set of regression analyses. With the exception of predictions for IH SkinPerm, none of those identified a notable relationship between experimental and computed values (r² ranging from 0 to 0.011). For IH SkinPerm, the relationship between predicted and experimental values could be described by the linear equation y = αx + β, where α = 0.231 (p < 0.01) and β = 1.628 (p < 0.01), with r² = 0.237. An only marginally better fit was obtained for the cubic equation y = αx³+ βx²+ γx + δ (α = −0.076 (p < 0.05), β = 0.051 (p = 0.065), γ = 0.328 (p < 0.01), δ = 1.601 (p < 0.01)), with r² = 0.253.

Figure 1. Relation of the experimental (logExperimentalAbsorption) to the predicted absorption (logPredictedAbsorption). The experimental absorption was taken as the sum of the relative amount present in the receptor fluid, the skin sample and the tape strips – except for the first and second ones. The predicted absorption represents the absorption as calculated by each model. The diagonal line represents the function y = x.

Model ranking

Since the models generally failed to predict the experimental values, the logarithm of the ratio between the predicted and experimental values, the log prediction factor (log PF), was used to rank the models according to the criteria outlined in the materials and methods section. This analysis aimed at identifying those models that, although unable to yield accurate results, can still provide estimations that will almost always lie either close to or above the experimental values. As shown in Table 2, predictions with COSMOS & TNO underestimated the experimental absorption in 1.4%, those with IH SkinPerm in 2% of the cases. However, the degree of underestimation did not exceed more than one order of magnitude for IH SkinPerm (minimum value for log PF is −0.6), whereas for COSMOS & TNO this was by up to almost 3 orders of magnitude (minimum value for log PF −2.7). When applying the rest of the criteria, IH SkinPerm was identified as the best performing model in the present study. Substantial overestimation of experimental data was 15% by this model. This percentage is small compared to the other models. Finally, for IH SkinPerm, the cases that were best estimated (0≤ log PF<1) were the most amongst the examined combinations, with a value of 38%.

Table 2. Best performing models ranked according to the four criteria.

	Criterion 1		Criterion 2	Criterion 3		Criterion 4
	Low percentage of underestimated cases (log PF<0)		Low factor of under-estimation	High percentage of cases estimated within one order of magnitude (0≤ log PF<1)		Low percentage of highly overestimated cases (log PF≥3)
	Absolute number	Percentage	Minimum log PF	Absolute number	Percentage	Absolute number	Percentage
IH Skin Perm	11	2	−0.6	216	38	84	15
COSMOS & TNO	8	1	−2.7	187	33	145	26
COSMOS & Buist	14	2	−2.7	188	33	145	26
Potts Guy & Buist	22	4	−2.8	206	37	127	23
Cleek Bunge & Buist	21	4	−2.7	205	36	127	23
Mitragotri & Buist	27	5	−2.8	205	36	122	22
Potts Guy & TNO	64	11	−4.0	187	33	128	23
Mitragotri & TNO	75	13	−4.4	194	34	121	21
Cleek Bunge & TNO	275	49	−6.0	153	27	34	6
Total	564	100		564	100	564	100

Highlights in bold indicate best performing model(s).

Factors influencing models’ performance

The absorption process is known to be influenced by multiple factors. Some of them are associated with the physical-chemical properties of a substance, e.g. lipophilicity, with the substance’s direct or indirect environment, e.g. the formulation in which it is present, occlusion of the skin after contact to the substance, or with the biological condition of the skin, e.g. irritancy status.

Lipophilicity and molecular weight

Lipophilicity is described as the octanol/water partition coefficient. The relationships between the log prediction factor of each model (log PF) and the active substance’s lipophilicity (log k_ow) are shown in Figure 2. Cases with log PF<0 (indicated by the horizontal lines) tend to have lower log k_ow values. Relationships between the prediction factors and the active substance’s molecular weight (MW) are shown in Figure 3. No relationship between the MW and the prediction factor was identified. Pearson’s correlations were computed for both parameters and are presented in Table 3. Overall, there was a moderate, positive correlation between the log prediction factor and lipophilicity in the models Potts Guy & TNO, Mitragotri & TNO and Cleek Bunge & TNO. Poor to zero correlations existed between the examined variables in the remaining models. In addition, there was a poor correlation between the prediction factor and the molecular weight in all models except for COSMOS & TNO, COSMOS & Buist, and IH SkinPerm, where no correlation at all was apparent.

Figure 2. Relation of the prediction factor to the active substance’s lipophilicity (log k_ow). Dermal absorption was predicted using the models indicated in the figure and prediction factors (PF) were calculated by dividing by experimental values. The logarithm (log PF) is plotted against measured octanol-water partition coefficients (log k_ow). In figures (a,c) open circles represent the cases out of the applicability domain. The horizontal lines mark the limit, under which the predicted is lower than the experimental absorption and the vertical the limit under which a substance is considered hydrophilic.

Figure 3. Relation of the prediction factor to the active substance’s molecular weight. Dermal absorption was predicted using the models indicated in the figure and prediction factors (PF) were calculated by dividing by experimental values. The logarithm (log PF) is plotted against the molecular weight. In figures (a and c) open circles represent the cases out of the applicability domain. The horizontal lines mark the limit, under which the predicted is lower than the experimental absorption.

Table 3. Summary of the Pearson correlation coefficients calculated for the relationship between the log PF (log prediction factor) and lipophilicity (expressed as log k_ow) or molecular weight (MW). AD: values within applicability domain included only.

Model	Pearson correlation coefficient between log PF and log kow	p Value of the relationship between log PF and log kow	Pearson correlation coefficient between log PF and MW	p Value of the relationship between log PF and MW
Potts Guy & TNO	0.563	<0.01	0.157	<0.01
Potts Guy & TNO (AD)	0.533	<0.01	0.13	<0.01
Mitragotri & TNO	0.596	<0.01	0.131	<0.01
Cleek Bunge & TNO	0.746	<0.01	0.174	<0.01
Cleek Bunge & TNO (AD)	−0.046	0.652	−0.184	0.069
COSMOS & TNO	0.084	<0.05	0.04	0.375
Potts Guy & Buist	0.195	<0.01	0.132	<0.01
Mitragotri & Buist	0.189	<0.01	0.136	<0.01
Cleek Bunge & Buist	0.187	<0.01	0.129	<0.01
COSMOS & Buist	−0.014	0.739	0.064	0.129
IH SkinPerm	0.058	0.172	0.056	0.183

The dataset’s molecular weight correlated with the lipophilicity (r = 0.529, p < 0.01). Therefore, additional partial correlation testing was performed to statistically remove the influence of the molecular weight on the relationship between lipophilicity and log PF. The tests revealed that in most of the models the correlations between the log k_ow and the log PF are only slightly influenced by the molecular weight. The only model that presents a more considerable difference in the correlation between log k_ow and log PF when controlling for molecular weight was that of Cleek Bunge & TNO when used in its applicability domain. The correlation coefficient increased from −0.046 (p = 0.652) to 0.207 (p = 0.042), meaning that even though no statistically significant correlation existed initially, controlling for molecular weight yielded a poor correlation. This may be attributed to the fact that in the restricted sample that correctly applies to the model (n = 98) there was a stronger correlation between lipophilicity and molecular weight (r = 0.843).

Physical state

Figure 4 depicts the distribution of the prediction factors of all models based on the physical state of the tested product, liquid, paste, and solid. Significantly higher prediction factors result for products applied as solids than products applied as liquids or pastes in the models Potts Guy & TNO and Cleek Bunge & TNO. More specifically, for Mitragotri & TNO and IH SkinPerm, it was revealed that products applied as solids have significantly higher prediction factors than products applied as liquids. For COSMOS & TNO, Potts Guy & Buist, Mitragotri & Buist, Cleek Bunge & Buist, and COSMOS & Buist non-parametric testing showed that pastes and solids have significantly higher prediction factors than liquids.

Figure 4. Distribution of the prediction factor (log PF) in different physical states. The left whisker marks the minimum value and is defined as 1.5xIQR, the left part of the box represents the 1st quartile, the middle vertical line in the box represents the median, the right part of the box represents the 3rd quartile, and the right whisker the maximum value, defined as 1.5xIQR. The circles mark the outliers and the asterisks the extremes. The vertical line in the diagram marks the limit, under which the predicted is lower than the experimental absorption. Δ marks statistically significant differences to liquids; Π to pastes, as identified by non-parametric testing performed within each model, with alpha value set to 0.05.

Other factors: formulation type, skin thickness, irritancy, occlusion status after application

For Potts Guy & TNO, Mitragotri & TNO and Cleek Bunge & TNO, the prediction factor for suspensions or ‘other’ formulations was highest among all formulation types. Products formulated as emulsions or dispersions have similar log PF values to each other, yet lower than suspensions but higher than solutions. In general, solutions have the lowest log PF compared to the rest of the categories. For COSMOS & TNO, non-parametric testing revealed that products formulated as emulsions or solutions have significantly lower log PF values than products formulated as suspensions. Emulsions had lower prediction factors than products that belong to the category ‘other’. No other statistically significant effect was detected within the model COSMOS & TNO. For Potts Guy & Buist, Mitragotri & Buist, and Cleek Bunge & Buist the tests revealed that products formulated as emulsions have significantly lower log PF than products formulated as suspensions or fall into the category ‘other’. Solutions have significantly lower prediction factor than suspensions. The tests did not detect any other statistically significant differences. For COSMOS & Buist, the non-parametric tests revealed that products formulated as emulsions have significantly lower log PF than suspensions or products that belong to the category ‘other’. No other significant differences exist between the prediction factors of the examined formulation types. For the model IH SkinPerm, significance testing revealed that products formulated as emulsions have lower log PF than products formulated as suspensions, dispersions, or ‘other’.

Non-parametric testing between the skin preparations isolated epidermis and dermatomed skin showed no statistically significant difference between the log PFs in any model. No statistically significant differences existed among the log PFs of different skin irritancy types for the models Potts Guy & TNO, Mitragotri & TNO, Cleek Bunge & TNO and IH SkinPerm. For COSMOS & TNO, Potts Guy & Buist, Mitragotri & Buist, Cleek Bunge & Buist and COSMOS & Buist the log PF of slight irritant substances was statistically significantly higher than the log PF of non-irritants. No other significant difference was detected regarding the skin irritancy.

Comparison of the distribution of the prediction factors of the models of Potts Guy & TNO, Mitragotri & TNO and Cleek Bunge & TNO based on the occlusion state under which the experiments were conducted (non-occluded, semi-occluded, and occluded), revealed no significant differences among the occlusion states. For COSMOS & TNO, Potts Guy & Buist, Mitragotri & Buist, Cleek Bunge & Buist, and COSMOS & Buist the prediction factors tend to be higher under semi-occlusive or non-occlusive conditions, than the prediction factors calculated when the skin sample is occluded. Statistically significant differences exist in the log PF among all states of occlusion for IH SkinPerm. Under semi-occlusive or non-occlusive conditions, the prediction factors are higher compared to situations where the skin sample is occluded. The log PF of experiments conducted with non-occluded skin is lower than that of experiments conducted with semi-occluded skin.

All detailed results can be found in the supplementary material.

Discussion

Models’ performance

This study demonstrated that the examined models could not predict the experimental values adequately, even after statistical fitting. However, their inaccuracy does not necessarily mean that they could not offer useful support for risk assessment. In fact, IH SkinPerm overpredicted in many cases (38%) within a 1- to 10-fold, and underpredicted only 2% of the dataset, suggesting that it could provide a conservative, yet a reasonably accurate model for pesticides. At this point, it should be stressed that the model was not used as recommended by its developers. Almost no values for log k_ow and very few for water solubility were available at pH 5.5. All density values were taken as 1000 mg/cm³ after ensuring that even for higher deviations from 1000 mg/cm³, the results would still lie within the 90% confidence interval. Even though the developers recommend the use of 1 cm for the stagnant air layer, it was taken as 3 cm, to cover worst case scenarios for occlusion/semi-occlusion. It is therefore expected that this adjustment might lead to higher estimations than it would otherwise do in experiments where non-occluded skin was examined. COSMOS was the second-best performing model. The equation for the estimation of k_p is similar to that of Potts Guy or Mitragotri, and its derivation relies on the assumptions of Potts Guy. The advancement of this model is that it is based on a richer dataset where confidence scores were assigned to each datapoint.

The model of Cleek Bunge was applied as a correction of the Potts Guy model for highly lipophilic molecules to account for the impaired permeation arising from the resistance offered by the hydrophilic viable epidermis. According to its developers, correction with Cleek Bunge significantly affects the estimate of dermal absorption only when its applicability domain is considered. This means that the use of this model for substances not belonging to its applicability domain would alter the performance of Potts Guy only marginally. The Cleek Bunge model has been employed successfully in previous studies either as a stand-alone model or as a correction to the highly lipophilic substances [23,27,28]. However, the present work provides evidence that the model is only suitable as a correction within its applicability domain and cannot, in any case, be applied as a stand-alone model. Application of Cleek Bunge to the total sample resulted in an underprediction of the absorption of almost half of the cases. Restricting the dataset to the substances for which the model is recommended (17% of the whole dataset), resulted in more accurate predictions for 7.3% of the whole dataset.

Visual evaluation of the distribution of the predicted absorption values (Figure 2(e–h)) revealed that the results of any combinations of Buist were similar to one another, whereas when the k_p-estimating models are combined with the TNO formula, different distribution profiles are obtained. A possible explanation could be that the results of the Buist combinations are less dependent on the k_p than those of TNO. To prove this assumption, we assessed the sensitivity of the corresponding equations to changes of the k_p analytically, by comparing their partial derivatives with respect to k_p. After transformation and comparison of the corresponding equations we could show that the results from the combination of the TNO equation with the models of Potts Guy, Mitragotri, Cleek Bunge, and COSMOS represent differences in predicted k_p values, while combinations of the Buist equation are less sensitive to changes of the permeability coefficient (see section 4 of the supplementary material).

Permeation of hydrophilic molecules is usually described poorly

Analysis of the results in the models of Potts Guy & TNO, Mitragotri & TNO and Cleek Bunge & TNO showed that the less lipophilic an active substance is, the lower the prediction factor. These results indicate that the examined models are particularly dependent on lipophilicity and usually predict the absorption of hydrophilic solutes unsatisfactorily. Given these observations, the question arises as to whether there are more accurate ways to describe the transdermal transport of hydrophilic molecules. Wilschut and colleagues had undertaken the first step towards this direction. This group attempted to improve the predictivity of the state-of-the-art model of Potts Guy. However, even after adapting the underlying equation and achieving an overall better fit, the QSAR still failed to predict the permeability of highly hydrophilic molecules compared to the other models examined in the study [29]. Later on, Mitragotri proposed a mechanistic model based on analytical solutions of the diffusion equation [19]. In this model, four pathways were considered: free-volume diffusion through lipid bilayers, lateral diffusion along lipid bilayers, diffusion through pores, and diffusion through shunts. The results from the equation that describes the transport of lipophilic chemicals (log k_ow>1) as free volume diffusion through lipid layers (the model also examined in the present work) were in perfect agreement with the results from the Potts Guy QSAR. This outcome is consistent with the results of the current study. Mitragotri overcame the limitation of Potts Guy by showing that permeation of hydrophilic solutes can be explained mechanistically and mathematically by diffusion through aqueous pores (equation not evaluated in the present study).

Other possible routes through which hydrophilic solutes may cross the stratum corneum have also been suggested: the lateral and the transcellular pathways [30]. Kasting and colleagues were among the first to propose the transcellular route for the transport of hydrophilic chemicals [31]. Another study provided evidence for the existence of a transcellular route with preferential corneocyte partitioning as the possible pathway for hydrophilic molecules [32]. A drawback of this model is that it does not account for the anisotropic nature of the lipid bilayers. This limitation is absent from the microscopic model of Wang et al., which considers three transport routes in the stratum corneum: the lateral lipid diffusivity, D_lat, the transverse mass transfer for hopping across lamellar bilayers, k_trans, and the isotropic corneocyte-phase diffusivity, D_corn [33]. Their results support the theory that solutes, regardless of their lipophilicity, permeate the stratum corneum in a predominately transcellular manner along with transversal, transbilayer hopping. From a morphological point of view, the latter pathway makes sense since every reasonable representation of the lipid bilayer demonstrates that solutes must traverse the lipid bilayers, even though this diffusivity is lower than the lateral diffusivities. Also, corneocyte partitioning is of high importance when it comes to transporting material to the next set of lipid bilayers as the permeant crosses the stratum corneum. Models accounting for these processes may show improved predictivity when evaluated against our dataset.

Additional factors influencing absorption

Correlation analysis revealed that the predictivity was not influenced by the molecular weight of the active substances in the examined dataset. From a physiological point of view, dermal absorption is thought to decrease in particular when compounds larger than 500 g/mol are presented to the skin [34,35]. The examined dataset is restricted regarding the molecular weight range with only three compounds exceeding 500 g/mol. This could explain why the model predictions do not differ within the MW range covered here. The physical form of a chemical (liquid or solid) influences its effective concentration on the skin [36]. Taking into account the morphology of the stratum corneum and by hypothesising that the penetration of the stratum corneum occurs mainly through the intercellular lipids surrounding the keratinocytes, it is expected that only a fraction of particles on the skin surface is in contact with the intercellular lipid [37]. This may pose a rate-limiting factor to the penetration of solids since the distribution of particles on the skin surface is not as homogenous as of liquid formulations. In fact, previous statistical analysis of experimental data on pesticide absorption showed that overall solid formulations are absorbed less than liquids [38]. The results of the present study are consistent with that finding. Non-parametric testing revealed that, compared to liquids, significantly higher prediction factors result in all models for solids. This means that the models cannot account for the decreased absorption of solids, resulting in considerably higher predictions compared to the experimental values. This observation is expected for almost all models since they are trained with aqueous solutions. The only model that should account for the physical state of the applied chemical is IH SkinPerm. However, IH SkinPerm showed similar responsiveness as the other models, regarding the physical state. It is essential to bear in mind that due to the nearly 40-fold difference between the sample sizes (n_liquid = 533; n_solid = 13) this statistical result is to be interpreted with caution.

Non-parametric analysis revealed that the models respond differently to different formulation types. A vehicle that interacts with the stratum corneum lipids or proteins can either accelerate or slow down the compound’s permeation into the skin [2]. Many products contain additional chemicals like surfactants and solvents to stabilize the formulation. These can alter the barrier properties of the stratum corneum by increasing the lipid bilayer’s fluidity resulting in an enhancement of skin permeability [39]. In general, solutions showed the lowest prediction factors. Compared to solutions, higher log PF values resulted for emulsions, indicating that emulsions are generally absorbed less than expected. Emulsions contain active substances dissolved in a solvent (typically oil) and emulsifiers. Upon dilution, the concentrate will form an oil-in-water emulsion. The droplet size, the emulsifier, and the surfactant organisation (micelles, lyotropic liquid crystals) in the emulsion may affect percutaneous absorption [40]. However, due to the lack of information on the surfactant and the thermodynamic activity of the active substance in its solvent, it is difficult to identify which parameters may lead to decreased absorption compared to solutions. Suspensions exhibited the highest log PF values, meaning that the active substances included in the mixture typically permeated the skin less than calculated by the models. A possible explanation for that may be that in suspensions the active substances are present as solid particles and will not be able to permeate the skin as readily as from other formulations.

Occlusion of the skin may lead to higher permeability rates mainly due to enhanced hydration that can be caused by the reduced trans-epidermal water loss [36]. This is in agreement with the observation of lower prediction factors for experiments conducted with occluded skin compared to those conducted with semi-occluded or unoccluded skin. Interestingly, prediction of dermal absorption measured under occlusive conditions was slightly better for most models. It may be speculated that this is also due to the training dataset which was obtained mostly under infinite and thus occlusive conditions. Occlusion of the skin after application of volatile compounds has been shown to enhance the permeation of all the compounds compared to non-occluded skin [41]. Of the models examined, only IH SkinPerm allowed correcting for volatilisation. However, since we choose to assume occlusion status always as occluded/semi-occluded, it can be expected that IH SkinPerm overestimated absorption for experiments conducted with unoccluded skin to a greater extent. The results from the statistical tests are in agreement with this hypothesis.

The thickness of the skin may influence the extent of permeation and the distribution of solutes. Although data on the role of skin thickness in in vitro percutaneous permeation studies are inconsistent, it is believed that thinner skin may contribute to the enhanced flux of some compounds [42,43]. Therefore, a lesser degree of overestimation was anticipated for experiments conducted with epidermal membranes. Contrary to expectations, no statistically significant differences were identified between dermatomed skin and isolated epidermis, indicating that the models’ performances are not influenced by the thickness of the skin.

Compounds with the potential for skin irritation and/or corrosion may harm the skin barrier and thus contribute to the increase of dermal permeation. Skin irritants have been shown to damage the stratum corneum lipid bilayer and decrease the skin’s ceramide content [39]. After treatment of skin with irritants such as sodium lauryl sulphate (SLS), enhancement of dermal penetration could be observed both in vivo and in vitro [44,45]. The magnitude of the enhancement is dependent on the compound, the experimental procedure, and the severity of irritation. Contrary to these observations, the present study shows that the models did not underestimate the absorption of irritants or corrosives. Overall, no statistically significant differences were detected in the predictive quality of the models concerning skin irritancy. The most reasonable explanation for that is that the statistical power of the tests might have been compromised due to the underrepresentation of irritants and corrosives in the sample. Also, in this study, irritancy of the applied dose was inferred from properties of the active substance only, as respective information on identity and properties of the co-formulants was not available. A concise overview about the parameters affecting skin absorption, expected shortcomings of model predictivity related to these parameters as well as the observed relationship between the parameter and model output is provided in Table 4.

Table 4. Summary of the examined factors: expected vs. observed outcomes on the predictivity of the models.

Examined factor	Expected outcome	Observed outcome
Lipophilicity	The absorption of hydrophilic compounds may be underestimated by the models which account only for the lipophilic pathway.	For some models no effect. For some models the absorption of hydrophilic compounds is underestimated.
Molecular size	The absorption of large/small molecules may be over-/underestimated.	No effect.
Physical state	The absorption of solids may be overestimated, since most models are trained on liquids.	For all models the absorption from solids (and sometimes pastes) is overestimated to a greater extent than for liquids.
Formulation type	The absorption of dispersions, emulsions, and suspensions may be overestimated. The absorption of solutions may be overestimated only to a lesser extent.	For some models the absorption from emulsions or suspensions is overestimated to a greater extent. For some models the absorption from solutions is underestimated or overestimated to a lesser extent.
Skin thickness/type	The absorption on isolated epidermis may be underestimated.	No effect.
Irritation properties of active substance	The absorption of compounds with the potential to cause irritation may be underestimated.	For some models no effect. For some models the absorption of slight irritants was overestimated to a greater extent than non-irritants.
Occlusion after application	The absorption of compounds where the skin is occluded after application may be underestimated.	For some models no effect. For some models the absorption is overestimated to a greater extent for semi-occlusive or non-occlusive conditions compared to occlusive conditions.

Various efforts have been made to improve QSARs for dermal absorption by including more structural descriptors. These models could be included in a follow-up to our study, but the required input values are not always easily accessible. Abraham and co-workers applied the theory of linear free energy relationships (LEFR) relating solvatochromic properties to skin permeability through parameters, referred to as Abraham descriptors, such as hydrogen bond acidity $\sum {\alpha }_2^{\mathrm{H}}$ and basicity $\sum {\beta }_2^{\mathrm{H}}$, polarisability ${\pi }_2^{\mathrm{H}}$, excess molar refraction ${\mathrm{R}}_2$ and the McGowan volume ${\mathrm{V}}_{\mathrm{x}}$ [46]. In addition to the solvatochromic considerations of Abraham, Patel and co-workers [47] developed a multiple linear regression (MLR) QSAR which identified log k_ow, MW, sssCH (the sum of E-state indices for all methyl groups) and ABSQon (the sum of average charges on oxygen and nitrogen atoms) as key parameters. The importance of lipophilicity, size, and molecular stereochemical properties has also been underlined in Baert, where the authors statistically evaluated the contribution of 1630 parameters in dermal permeation [48]. Riviere and Brooks incorporated the mixture factor (MF) to account for the interactions within the vehicle relative to the compounds’ subsequent partitioning into the stratum corneum [49–51], which improved predictions in cases where chemicals are presented to the skin in non-aqueous vehicles/mixtures. Guth applied the MF approach and developed models explicitly targeted for the prediction of dermal absorption from agrochemical formulations using the Abraham descriptors to express the active substances’ chemical properties, log k_ow, the topological polar surface area (TPSA) and number of hydrogen bond acceptors/donors (HBA/HBD) to describe the properties for the formulation components, as well as an indicator variable for the species human and rat [52]. Nevertheless, this approach would still lack consideration of possible interactions of the mixture components with the skin, which could enhance absorption (e.g. through extraction or fluidization of the lipid bilayer) or could cause irritation and/or corrosion. A step forward in that direction could be the examination of models developed specifically for transdermal enhancers, as discussed in Tsakovska [53]. Unfortunately, these latter approaches could not be examined against our dataset due to a lack of detailed information on formulation composition, which is usually treated as confidential business information in the case of plant protection products. Finally, there have been attempts to improve model performance by including descriptors of exposure such as dose, state of occlusion, skin thickness, hydration and time [23,24,27,54].

Replacement/refinement of EFSA default values for regulatory risk assessment

In addition to the performance and factor analysis, we also wanted to examine whether the models can offer a useful replacement or refinement alternative to the current EFSA default values used in regulatory risk assessment of pesticides when experimental data is not available [5]. The first column of Table 5 shows the number of cases for which the application of the in silico model examined would refine the default value. That is, the predicted absorption would be lower than the default values as applicable according to the current EFSA guidance [5]. When predicted absorption was below the respective default value, we further examined whether estimates were above and below the experimental value. Predictions below the experimental value (refinement underestimates) may be problematic in the public health context. While the model combination Cleek Bunge & TNO would allow refinement of the default assumption in 77% (433 out of 564) of the cases in the examined dataset, application of this model provided a prediction below the experimental value in 60% of those cases. The only model with a low rate of refinement underestimates in a range that may be acceptable was IH Skin Perm with 9% of the refined cases being underestimates of the experimental value. However, when using this model, the refinement of the default assumption was possible in 65 out of 564 cases (11%) of the total dataset only.

Table 5. Applicability of the models as a refinement alternative to the current EFSA default values. The number of cases for which the predicted absorption would be lower than the default values as applicable according to the EFSA guidance [5] was counted. When, for those cases, predicted absorption was below experimental absorption, these were classified as refinement underestimates.

Model	Refinement of default values applicable [n out of a total of 564 cases]	Refinement underestimates [percentage]
IH SkinPerm	65	9%
COSMOS & TNO	32	22%
Mitragotri & Buist	106	23%
Cleek Bunge & Buist	86	24%
Potts Guy & Buist	88	25%
COSMOS & Buist	42	26%
Mitragotri & TNO	139	52%
Cleek Bunge & TNO	433	60%
Potts Guy & TNO	101	61%

Highlights in bold indicate models most suitable for refinement (Cleek Bunge & TNO) and with the lowest percentage of refinement underestimates (IH SkinPerm).

Ultimately, any model used in risk assessment should predict the situation in humans under real life conditions. As shown previously, in vivo data obtained in the rat is greatly overestimating absorption across the skin in humans, while in vitro human data is more appropriate [55]. Therefore, in vivo rat data was not used in the reference dataset for this study and preference was given to human in vitro (ex vivo) data that is now widely accepted [5]. Use of information calculated from rat in vivo, rat in vitro and human in vitro data by the so called triple pack approach [55] may provide a better description of the real-life situation in some cases, but would have significantly reduced the size of our database, while the impact on precision and accuracy of the reference dataset would have been unclear. The latter is generally an issue to be considered when interpreting analyses of any models’ predictivity. In regulatory practice, variability factors are thus included, with a factor of 3.2 assigned to intraspecies differences in toxicokinetics, and variability in dermal absorption testing is explicitly taken into account in evaluation of the data [5]. Further studies with a richer dataset may be required before recommendations on the use of computational models, including IH SkinPerm, for the replacement or refinement of established default value for dermal absorption of pesticides are issued.

Acknowledgements

The authors wish to thank Dr. Hans Mielke and Prof. Dr. Ursula Gundert-Remy for their support throughout the project. Data collection was funded by BfR SFP 1322-606.

Disclosure statement

No potential conflict of interest was reported by the authors.

Supplementary material

Supplemental data for this article can be accessed at: https://doi.org/10.1080/1062936X.2019.1644533.

Additional information

Funding

This work was supported by the Bundesinstitut für Risikobewertung [Data collection was funded by BfR SFP 1322-606.].