GARDskin dose–response: methodological description
The conventional GARDskin protocol assays chemicals at single concentrations and performs hazard classifications into the categories: skin sensitizer and non-sensitizer. Although results from the assay are binary, it has been noted that the observed effects seem to be dose-dependent, demonstrating an inverse correlation between concentrations at which chemicals have been assayed, and their relative sensitizing potency. However, by itself, the GARD input concentration may not be sufficiently predictive, as observed previously25. Therefore, to further explore the relationship between chemical exposure concentrations and GARDskin responses, and its subsequent association with sensitizing potency, dose–response measurements were incorporated into the GARDskin protocol. Specifically, it was hypothesized that the lowest concentration required to exceed the binary classification threshold (DV ≥ 0) and generate a positive classification, referred to as the cDV0 concentration, would be informative of sensitizing potency. As such, the proposed dose–response measurements would also be aligned with common toxicological principles for potency assessments, where the response value (DV), the binary threshold (DV ≥ 0), and the derived cDV0 concentrations may be viewed as an analogue to the response value (stimulation index (SI)), the binary threshold (SI ≥ 3), and the EC3 concentrations in the LLNA assay, respectively. To evaluate this premise, a total of 29 chemicals, comprised of 7 non-sensitizers and 22 skin sensitizers of varying potency, were tested in a titrated range of concentrations. Given that the experimental setup included several novel elements, prior knowledge of the dose–response relationships’ characteristics were lacking. Therefore, the exposure concentrations were selected to evenly sample several concentrations of each chemical instead of favoring a few doses with replicates. This approach was preferred due to the reduced risk of selecting poor concentrations which could have resulted in poor estimates of the dose–response relationships26,27. The standard GARDskin pipeline was used to generate response values, by applying the SVM model to assign DVs to each chemical and concentration, based on the expression levels of the genes in the GARDskin biomarker signature. For each tested chemical, the assayed concentrations were plotted against the responses, i.e., the DVs, as illustrated in individual dose–response plots in Fig. 1.
Dose–response relationships for the examined chemicals. Points represent measured response values (decision values) at respective concentrations. Points with white centers represent the GARD input concentrations. Dotted lines describe the log-logistic models, and the solid lines describe the linear interpolations. *The running median was used as positive interpolation point to reduce the potential impact of noise.
Dose–response relationships and cDV0 estimation
Following inspection of the acquired data, log-logistic models were fitted to each chemical to attempt to model the dose–response relationships. No-effect tests were applied to assess the presence of concentration-dependent increases in responses. The p-value for respective chemical and test is described in Table 1. As can be seen, the p-values for salicylic acid, xylene, 1-butanol, glycerol, octanoic acid, vanillin, trans-anethol, and kanamycin sulfate, agree with observations of the raw data (Fig. 1) and suggest that these chemicals did not induce increased responses in the assay over the considered concentrations. The results are further suggestive of weak responses for benzyl salicylate and benzocaine. Examination of the data showed that benzocaine did not exhibit any positive decision values, though an upward trend at the highest concentrations could potentially still be recognized, which would explain the relatively low p-value associated with the model’s fit. The two highest concentrations of benzyl salicylate generated positive responses, though with relatively low magnitudes. Remaining chemicals exhibited more distinct dose–response relationships with more significant p-values from the no-effect tests, suggesting quantifiable effects on the response values.
Next, cDV0 values were estimated for each chemical. Two different approaches were considered for the calculations. First, the concentrations were estimated using the already fitted log-logistic models by determining the curves’ x-axis intercepts. The log-logistic models also allowed for characterization of the uncertainty of respective cDV0 estimate, by calculation of approximate 95% confidence intervals, see Table 1. The log-logistic models do, however, make certain assumptions regarding the underlying data and the relationship between concentrations and response values, which might not necessarily be valid, while simultaneously requiring a relatively complex fitting procedure. Therefore, linear interpolation between two points on adjacent sides of the decision border was considered as an alternative, simpler and less constrained, method for estimation of cDV0. The derived concentration estimates from both approaches are compared in Fig. 2. As can be seen, both methods generally produced very similar concentrations. Three chemicals received inconsistent estimates when the two methodologies were compared. These were 2-hydroxyethyl acrylate, linalool, and isoeugenol. When examining the dose–response relationships for each of these chemicals, it can be seen that the observed discrepancy for 2-hydroxyethyl acrylate arose due to the two positive points present at the lowest exposure concentrations, giving an apparent U-shape to the overall response. Because of the unexpected appearance of these data, the chemical was excluded from further assessment. For linalool, a cDV0 concentration was only definable with linear interpolation. For isoeugenol, the conflicting estimates followed from the shape of the dose–response curve, which seemed to deviate from the sigmoidal shape expected by the log-logistic models. Instead, the data suggested that the chemical was able to render positive decision values at moderate magnitudes already at low concentrations, which was followed by a plateau until concentrations reached > 180 μM, at which a second increase in the decision values was observed. Given this appearance, it can be argued that the locality of the linear interpolation provides a better estimate of the true cDV0-value, since its fit is not constrained by a pre-determined shape of the curve. However, additional data would be useful for determining the induced dose–response curve’s shape with more certainty. Because of these observations, and the overall similarity between the results, the linear interpolation method was selected as the most appropriate technique for estimating cDV0, and the results described below are based on these interpolated values, which are described in Table 1.
Comparison between the cDV0 estimation procedures. The scatter plot describes cDV0 values that were determined using either log-logistic modelling (x-axis) or linear interpolation (y-axis). The dashed line represents the identity line.
Comparison of cDV0 with established potency measures
A visual representation of the experimentally derived dose–response curves for all chemicals are presented in Fig. 3, where the individual dose–response curves are colored according to their respective potency sub-category as defined by UN GHS (1A, 1B, no Cat) by the application of a 2% LLNA EC3 threshold to discriminate the stronger sensitizers from the weaker. As apparent by the figure, two of the LLNA non-sensitizers, benzyl alcohol and phenol, exhibited response-values sufficient to exceed the binary classification threshold (DV ≥ 0) at any of the assayed concentrations. The positive responses corresponded to several data points of benzyl alcohol, which is nevertheless generally considered a human skin sensitizer, and the highest evaluated concentration of phenol (500 μM). Further inspection of the figure and the generated dose–response curves clearly reveals a relationship between GHS potency and assayed concentrations, where stronger sensitizers (1A) required lower concentrations to exceed the binary classification threshold than did the weaker sensitizers (1B). Interestingly, with only a single exception for the chemical diethyl maleate, which is a CLP 1A/1B borderline chemical, a complete resolution between the GHS potency classes was observed.
Visualization of acquired dose–response data. Each line represents a generalized additive model fitted to the dose–response data of a single chemical. The colors of the lines describe the GHS potency sub-categories, into which chemicals were assigned based on their LLNA results (1A: EC3 ≤ 2%; 1B: EC3 > 2%; No Cat: NS).
To further evaluate the hypothesis that cDV0 values were informative of the chemicals’ relative potencies, the correlations between cDV0 and available potency measures were determined. First, the correlation with LLNA was considered. Figure 4 displays a scatter plot which suggests that chemicals ranked as strong sensitizers in the LLNA also seem to obtain lower cDV0 values. The linear association between the two metrics was quantified using Pearson’s correlation coefficient, which was determined to be 0.81 (p = 9.1 × 10–5; n = 17), which is indicative of a strong relationship, suggesting that GARDskin cDV0 was significantly associated with sensitizing potency, as defined by LLNA EC3 values.
Performance assessment of GARDskin Dose–Response. Scatter plots displaying the relationship between estimated cDV0 values and (a) LLNA EC3 values and (b) human NOEL values. The lines represent linear regression models fitted to the data, and the shaded areas describe the 95% confidence intervals of the fits. Encircled datapoints represent indirectly acting haptens. GARDskin cDV0 is given in weight-based concentrations to facilitate comparison between potency measures.
Second, the relationship between human NOEL values and cDV0 values was considered and Fig. 4 shows the relationship between these two metrics. Similar to the LLNA-figure, a clear correlation can be observed. Indeed, quantification of the strength of the linear association resulted in a Pearson correlation coefficient of 0.74 (p = 1.5 × 10–3; n = 15), further supporting the observation that cDV0 values were informative of sensitizing potency.
Finally, indirectly acting haptens, i.e., chemicals that require either abiotic or biotic activation prior to being able to induce sensitization reactions, are also displayed in Fig. 4. As can be seen, the data suggest that cDV0 values are also informative of potency for these types of chemicals, suggesting that the GARDskin dose–response methodology is also able to assess indirectly acting haptens.
