A genomic biomarker signature can predict skin sensitizers using a cell-based in vitro alternative to animal tests

Background Allergic contact dermatitis is an inflammatory skin disease that affects a significant proportion of the population. This disease is caused by an adverse immune response towards chemical haptens, and leads to a substantial economic burden for society. Current test of sensitizing chemicals rely on animal experimentation. New legislations on the registration and use of chemicals within pharmaceutical and cosmetic industries have stimulated significant research efforts to develop alternative, human cell-based assays for the prediction of sensitization. The aim is to replace animal experiments with in vitro tests displaying a higher predictive power. Results We have developed a novel cell-based assay for the prediction of sensitizing chemicals. By analyzing the transcriptome of the human cell line MUTZ-3 after 24 h stimulation, using 20 different sensitizing chemicals, 20 non-sensitizing chemicals and vehicle controls, we have identified a biomarker signature of 200 genes with potent discriminatory ability. Using a Support Vector Machine for supervised classification, the prediction performance of the assay revealed an area under the ROC curve of 0.98. In addition, categorizing the chemicals according to the LLNA assay, this gene signature could also predict sensitizing potency. The identified markers are involved in biological pathways with immunological relevant functions, which can shed light on the process of human sensitization. Conclusions A gene signature predicting sensitization, using a human cell line in vitro, has been identified. This simple and robust cell-based assay has the potential to completely replace or drastically reduce the utilization of test systems based on experimental animals. Being based on human biology, the assay is proposed to be more accurate for predicting sensitization in humans, than the traditional animal-based tests.


Background
Allergic contact dermatitis (ACD) is a common inflammatory skin disease characterized by eczema and recurrent episodes of itching [1]. The disease affects a significant proportion of the population, with prevalence rates of 7.2% to 18.6% in Europe [2,3], and the incidence is increasing due to repeated exposure to sensitizing chemicals. ACD is a type IV delayed-type hypersensitivity response caused mainly by reactive T helper 1 (Th1) and interferon (IFN)γ producing CD8 + T cells, at site of contact with small chemical haptens in previously exposed, and immunologically sensitized, individuals [4].
Dendritic cells (DC) in the epidermis initiate the immune reactions by responding to haptens bound to self-molecules subsequently activating T cell-mediated immunity.
The REACH (Registration, Evaluation, and Authorization of Chemicals) regulation requires that all new and existing chemicals within the European Union, involving approximately 30.000 chemicals, should be tested for hazardous effects [5]. As the identification of potential sensitizers currently requires animal testing, the REACH legislation will have a huge impact on the number of animals needed for testing. Further, the 7th Amendment to the Cosmetics Directive posed a ban on animal tests for the majority of cosmetic ingredients for human use, to be in effect by 2009, with the exceptions of some tests by 2013. Thus, development of reliable in vitro alternatives to experimental animals for the assessment of sensitizing capacity of chemicals is urgent. To date, no non-animal replacements are available for identification of skin sensitizing chemicals, instead the preferred assay is the mouse Local Lymph Node Assay (LLNA) [6], followed by the Guinea pig maximization test (GPMT) [7]. An in vitro alternative to these animal models should exhibit improved reliability, accuracy and importantly correlate to human reactivity.
DCs play key roles in the immune response by bridging the essential connections between innate and adaptive immunity. Upon stimulation, they can rapidly produce large amounts of mediators that affect chemotaxis and activation of other cells at the site of inflammation, and can selectively respond to various pathogens and environmental factors, by fine-tuning the cellular response through antigen-presentation. Thus, exploring and utilizing the immunological decision-making by DCs during stimulation with sensitizers, could serve as a potent test strategy for the prediction of sensitization.
Factors that complicate and impede the use of primary DCs as a test platform include adaptable phenotypes and specialized functions of different DC subpopulations, in addition to their wide and sparse distribution. Thus, the development of assays based on the predictability of DC function must rely on alternative cell types or mimics of in vivo DCs. For this purpose, a cell line with DC characteristics would be advantageous, as it constitutes a stable, reproducible and unlimited supply of cells. MUTZ-3 is an unlimited source of CD34 + DC progenitors. Upon differentiation, MUTZ-3 can acquire phenotypes comparable to immature DCs or Langerhans-like DCs [8], present antigens through CD1d, MHC class I and II and induce specific T-cell proliferation [9]. Differentiated MUTZ-3 can also display a mature transcriptional and phenotypic profile upon stimulation with inflammatory cytokines [10].
In this report, we present a novel test principle for the prediction of skin sensitizers. To simplify the assay procedures and increase reproducibility, we employed progenitor MUTZ-3 cells, without further differentiation, and subjected the cells to stimulation with a large panel of sensitizing chemicals, non-sensitizing chemicals, and controls. The transcriptional response to chemical stimulation was assessed by genome-wide profiling. From data analysis, a biomarker signature of 200 transcripts was identified, which completely separated the response induced by sensitizing chemicals vs. non-sensitizing chemicals and the predictive power of the signature was illustrated, using ROC curves. The biomarker signature includes transcripts involved in relevant biological pathways, such as oxidative stress, DC maturation and cytokine responses, which further could shed light on molecular interactions involved in the process of sensitization. In conclusion, we have identified a biomarker signature with potent predictive power, which we propose as an in vitro assay for the identification of human sensitizing chemicals.

Results
The cellular rationale of the in vitro cell culture system DCs are essential immunoregulatory cells of the immune system demonstrated by their unique property to recognize antigen for the initiating of T cell responses, and their potent regulatory function in skewing immune responses. This makes them obvious targets for assay development. However, primary DCs constitute a heterogeneous and minor population of cells not suited for screening and the choice would be a human DClike cell line, with characteristics compared to primary DCs. Since no leukemic cell line with DC-like properties has been reported [11], the generation of human DClike cell lines relies on available myeloid leukemia cell lines. MUTZ-3 is a human acute myelomonocytic leukemia cell line with a potent ability to mimic primary human DCs [11]. Similar to immature primary DCs, MUTZ-3 progenitors express CD1a, HLA-DR and CD54, as well as low levels of CD80 and CD86 ( Figure 1). The MUTZ-3 population also contains three subpopulations of CD14 + , CD34 + and double negative cells, previously reported to be transitional differentiation steps from a proliferative CD34 + progenitor into a Figure 1 Phenotype of MUTZ-3 cells prior to stimulation with sensitizing and non-sensitizing chemicals. Cell surface expression levels of CD14, CD1a, CD34, CD54, CD80, CD86 and HLA-DR were assessed with flow cytometry. Gates were set to exclude debris and dead cells, and quadrants were established by comparing with relevant isotype controls. Results are shown from one representative experiment out of six. non-proliferative CD14 + DC precursor [8]. Consequently, constitutively differentiating progenitor MUTZ-3 cells were used as the basis for a test system. CD86 surface expression in response to sensitizer stimulation CD86 is the most extensively studied biomarker for sensitization to date, using e.g. monocyte derived dendritic cells (MoDCs) or human cell lines and their progenitors, such as THP-1, U-937 and KG-1. Thus, as a reference, cell surface expression of CD86 was measured with flow cytometry after 24 h stimulation, using 20 sensitizers and 20 non-sensitizers, as well as vehicle controls ( Figure 2). CD86 was significantly up-regulated on cells stimulated with 2-aminophenol, kathon CG, 2-nitro-1,4-phenylendiamine, 2,4-dinitrochlorobenzene, 2-hydroxyethyl acrylate, cinnamic aldehyde, p-phenylendiamine, resorcinol, potassium dichromate, and 2-mercaptobenzothiazole.
Hence, an assay based on measurement of a single biomarker, such as CD86, would give a sensitivity of 50% and a specificity of 100%. Consequently, CD86 cannot classify skin sensitizers, using a system based on MUTZ-3 cells.
Analysis of the transcriptional profiles in chemically stimulated MUTZ-3 cells The genomic expression arrays were then used to test the same 20 sensitizers and 20 non-sensitizers, in triplicates. The vehicle controls, such as DMSO and distilled water, were included in twelve replicates. In total, a data set was generated based on 144 samples. RMA normalization and quality controls of the samples revealed that the oxazolone and cinnamic aldehyde samples were significant outliers and had to be removed, or they would have dominated the data set prohibiting biomarker identification (data not shown). In addition, one of the Figure 2 Changes in CD86 expression following stimulation with sensitizing and non-sensitizing chemicals. Cell surface expression levels of CD86 were monitored after stimulation with chemicals for 24 h. A). Chemical-induced up regulation of CD86, in terms of changes in frequency of positive cells, were determined by flow cytometry, as exemplified by the comparison of 2-aminophenol-stimulated cells (right dotplot) and unstimulated controls (left dot plot). Results are shown from one representative experiment out of three. Gates were set to exclude debris and dead cells, and quadrants were established by comparing with relevant isotype controls. B) Compilation of frequencies of CD86positive cells after 24 h of stimulation. Statistical analysis was performed using Student's t test. *p < 0.05, # p < 0.01.
replicates of potassium permanganate had to be removed due to a faulty array. This left a data set consisting of 137 samples, each with data from measurements of 29,141 transcripts. In order to mine the data set for information specific for sensitizers vs. non-sensitizers, the software Qlucore Omics Explorer 2.1 was used, which enable real time principal component analysis (PCA) analysis. The input genes were at the same time sorted after desired criteria, i.e. sensitizers and non-sensitizers, based on ANOVA p-value selection. Two different ANOVA analyses were performed ( Figure  3). First, Figure 3A and 3B show PCA plots based on 1010 transcripts with a p-value of ≤ 2.0 × 10 -6 , from a one-way ANOVA analysis, comparing sensitizing vs. 20 sensitizing and 20 non-sensitizing chemicals were assessed with transcriptomics, using Affymetrix Human Gene 1.0 ST arrays. Structures and similarities in the gene expression data set were investigated, using principal component analysis (PCA) in the software Qlucore. A) PCA of genes differentially expressed in cells stimulated with sensitizing (red) versus non-sensitizing (green) chemicals (1010 genes identified with one-way ANOVA). B) PCA of genes differentially expressed in cells stimulated with sensitizing versus non-sensitizing chemicals (1010 genes), but now samples are colored by the compound used for stimulation. C) PCA of genes differentially expressed when comparing the different stimulations with 2-way ANOVA (1137 genes). Samples are colored according to sensitizing (red) and non-sensitizing (green) chemicals. D) PCA of genes differentially expressed when comparing the different stimulations with 2-way ANOVA (1137 genes), but now samples are colored by the compound used for stimulation. P, p-value from ANOVA. Q, p-value corrected for multiple hypothesis testing. non-sensitizing chemicals. As can be seen in Figure 3A, a clear discrimination can be made between the two groups, with non-sensitizers forming a condensed cloud in the lower part of the figure (green), while sensitizers stretch upwards in various directions (red). However, a complete separation is not achieved between the two groups at this level of significance. From Figure 3B, now colored according to stimulating agent, it is evident that one or more replicate of glyoxal, eugenol, hexylcinnamic aldehyde, isoeugenol, resorcinol, penicillin G and ethylendiamine grouped together with the control group. In addition, one replicate or more of the non-sensitizers tween 80, octanoic acid and phenol grouped closely with the sensitizers. Secondly, Figure 3C and 3D show PCA plots based on 1137 genes, with p-values ≤ 7.0 × 10 -21 , from a multi-group ANOVA analysis, comparing each individual stimulation. Identifying this large number of genes at this level of significance provided strong indications of the power in the data set. In Figure 3D, it is clear that the replicates group together, indicating high quality data. The triplicate samples of potassium dichromate have a discrete profile, which demonstrate a substantial impact of the cells compared to non-sensitizers. Furthermore, 2-hydroxyethyl acrylate, 2-aminophenol, kathon CG, formaldehyde, 2-nitro-1,4-phenylendiamine, 2,4dinitrochlorobenzoic acid, p-phenylendiamine, 2-mercaptobenzothiazole, cinnamic alcohol and resorcinol have replicates that group together, separate from the negative group. Still, as can be seen in Figure 3C as well as in 3A, complete separation is not achieved with neither of the gene signatures of 1010 and 1137 genes both selected on p-values.
Backward elimination identifies genes with the most discriminatory power Even though the data set contains genes with p-values down to 1 × 10 -17 , lowering the p-value cutoff did not achieve complete separation between sensitizers and non-sensitizers. Gene signatures entirely selected on pvalues does not provide the best possible predictive power, since the information is per se not orthogonal. To further reduce the number of transcripts for a predictive biomarker signature, we employed an algorithm for backward elimination ( Figure 4A). The algorithm removes genes one by one while taking into account not only the impact of genes individually, but how they perform collectively with the entire selected gene signature. For each gene eliminated, the Kullback-Leibler divergence (KLD) value is lowered, until a breakpoint is reached, at which point 200 genes remained. Continuing eliminating genes at this point causes the KLD to rise again, indicating that information is being lost ( Figure  4A). Therefore, the 200 genes with lowest KLD value were selected for further analysis. PCA of the 200 analytes now revealed that they have the ability to completely separate sensitizers from non-sensitizers, indicating that these transcripts can be used as predictors for sensitizing properties of unknown samples ( Figure 4B). Importantly, by coloring the samples in the PCA by their potency, according to LLNA, it is clear that potency can also be predicted ( Figure 4C), as extreme and strong sensitizers tend to group further from the non-sensitizers, while moderate and extreme sensitizers group closer to non-sensitizers. The 200 genes are termed the "Prediction Signature" and their identities are listed in Table 1. In addition, the transcriptional profiles of the differentially expressed genes are presented in a heatmap ( Figure 5).

Interrogation of the analysis used to identify the Prediction Signature
To validate the predictive power of our signature, we used a machine learning method called the Support Vector Machine (SVM) [12], which maps the data from a training set in space in order to maximize the separation of gene expression induced by sensitizing and nonsensitizing chemicals. As training set, 70% of the data set was selected randomly and the entire selection process was repeated. Starting with 29,141 transcripts, the signature was reduced to 200 transcripts, termed "Test Gene Signature", using ANOVA filtering and backward elimination, as described above. The remaining 30% of the data set was used to test each signature. The partitioning of the data set into subsets of 70% training data set and 30% test data set was done in a stratified random manner, while maintaining the relation of sensitizers and non-sensitizers. Thereafter, the Test Gene Signature was used to train an SVM model with the training set, and the predictive power of the model was assessed with the test set. This entire process was iterated 20 times. The frequency by which each gene in the Prediction Signature was included in the Test Gene Signatures is reported in table 1. Figure 6A shows a PCA plot based on the Test Gene Signature from one representative iteration. Clearly, the separation between sensitizers and non-sensitizers resembles the one observed for the Prediction Signature in Figure 4B. In Figure 6A, the samples of the sensitizing and non-sensitizing chemicals in the test set have been colored dark red and dark green respectively, indicating that they are not contributing to the principal components of the plot, but are merely plotted based on their expression values of the selected Test Gene Signature. As can be seen, sensitizers from the test set group with sensitizers from the training set, while non-sensitizers from the test set group with non-sensitizers from the training set. The final outcome of the SVM training and validation can be seen in Figure 6B, where the areas under the ROC curve are plotted for each iteration. The average area under the ROC curve of 0.98 confirms the ability to discriminate sensitizers from control samples. Based on this average, the estimated prediction performance of the assay reveals an accuracy of 99%, sensitivity of 99% and specificity of 99%. While this experiment does not validate the prediction power of the Prediction Signature per se, it does indeed validate the method by which it has been selected, supporting the claim that the Prediction Signature is capable of accurately predicting sensitizing properties of unknown samples.

Interactome, molecular functions and canonical pathways involving the Prediction Signature
Using Ingenuity Pathways Analysis (IPA, Ingenuity Systems Inc.), 184 of the 200 molecules in the signature were characterized with regard to the interactome, known functions and (canonical) pathways. The remaining 16 molecules could not be mapped to any unique IPA entries. The dominating functions identified were small molecule biochemistry (39 molecules), cell death (33), lipid metabolism (25), hematological system development (18), cell cycle (18), molecular transport (17), cellular growth and proliferation (16), and carbohydrate metabolism (15) ( Table 2).
Pathways possibly invoked by the molecules in the signature were also investigated using IPA. Those most highly populated involved NRF2-mediated oxidative response (10), xenobiotic metabolism signaling (8), protein ubiquitination pathway (7), LPS/IL-1 mediated inhibition of RXR function (6), aryl hydrocarbon receptor signaling (6) and protein kinase A signaling (6). These pathways are known to take part in reactions provoked by foreign substances, xenobiotics, which supports a relevant biology behind the genomic signature.

Discussion
Allergic contact dermatitis (ACD) is an inflammatory skin disease caused by an adaptive immune response to normally innocuous chemicals [13]. Small molecular weight chemicals, so-called haptens, can bind self-proteins in the skin, which enables internalization of the protein-bound allergenic chemical by skin dendritic cell (DC). DCs, under the influence of the local    microenvironment, process the protein-hapten complex, migrate to the local lymph nodes and activate naïve T cells. The initiation and development of allergen-specific responses, mainly effector CD8+ T cells and Th1 cells, and production of immunoregulatory proteins, are hallmarks of the immune activation observed in ACD. ACD is also the most common manifestation of immunotoxicity observed in humans [13] and hundreds of chemicals have been shown to cause sensitization in skin [14]. The driving factors and molecular mechanisms involved in sensitization are still unknown even though intense research efforts have been carried out to characterize the immunological responses towards allergenic chemicals. The REACH legislation requires that all chemicals produced over 1 ton/year are tested for hazardous properties such as toxicity and allergenicity [5], which increase the demand for accurate assays with predictive power for hazard identification. Additionally, the 7th Amendment to the Cosmetics Directive (76/768/EEC) poses a complete ban on using animal experimentation for testing cosmetic ingredients by 2013 if a scientifically reliable method is available. Thus, there is a significant need for predictive test methods that are based on human cells. Today, the identification of potential human sensitizers relies on animal experimentation, in particular the murine local lymph node assay (LLNA) [6]. The LLNA is based upon measurements of proliferation induced in draining lymph nodes of mice after chemical exposure [15]. Chemicals are defined as sensitizers if they provoke a three-fold increase in proliferation compared to control, and the amount of chemical required for the increase is the EC3 value. Thus, the LLNA can also be used to categorize the chemicals based on sensitization potency. However, LLNA is, besides the obvious ethical implications, also time consuming and expensive. Human sensitization data often stem from human maximization tests (HMT) [16] and human patch tests (HPT). In an extensive report from the Interagency Coordinating Committee on the Validation of Alternative Methods (ICCVAM), the performance characteristics of LLNA were compared to other available animal-based methods and human sensitization data (HMT and HPT) [17]. The LLNA performance in comparison to human data (74 assessments) revealed an accuracy of 72%, a sensitivity of 72% and a specificity of 67%. Various human cell lines and primary cells involved in sensitization have been evaluated as predictive test system, such as epithelial cells, dendritic cells and T cells, however, no validated test assay is currently available. THP-1, U937, KG-1 and MUTZ-3, naive or differentiated, are among the human myeloid cell lines most extensively evaluated as platforms for DC-based in vitro assays, as reviewed in [18]. These cells are easy to grow  and enable standardization of protocols. U937 and THP-1 are currently being evaluated in pre-validation stage for prediction of skin sensitization. The Human Cell Line Activation Test (hCLAT) is based upon analysis of CD86 and/or CD54 expression on THP-1 cells after chemical stimulation [19,20]. The Myeloid U937 Skin Sensitization Test (MUSST) also involves analysis of CD86 [21]. These assays are thus very limited in Figure 6 Validation of selection procedure of Prediction Signature. The method by which the Prediction Signature was constructed was validated by repeating the process on 70% randomly selected data (training set   readout. As CD86 is among the markers most extensively studied, we evaluated the expression level of this marker in our assay. We demonstrated its relevance but also its insufficient predictive power (Figure 2), since only 10 out of 20 sensitizing chemicals induced a significant up regulation of CD86. Various other single biomarkers have been suggested to be up regulated upon stimulation with sensitizing chemicals, such as CD40, CD80, CD54, CXCL8, IL-1β, MIP-1β, p38 MAPK, as reviewed in [18], yet single-handedly, none of them have enough predictive power to discriminate between sensitizing and non-sensitizing chemicals. The analysis of biomarker signatures, i.e. combination of biomarkers, has been shown to be superior in molecular diagnostic of cancer and superior to any single biomarker. Consequently, we therefore utilized the power of global transcriptomics and screened the gene regulation induced by a large set of well-defined chemicals and controls in search of predictive biomarker combinations. The large number of differentially expressed genes in MUTZ-3 cells stimulated with sensitizing chemicals vs. non-sensitizing controls revealed that MUTZ-3 indeed had a capacity to differentiate between these two groups. Efforts have previously been done to create assays based on genome analysis in various cell systems, such as e.g. CD34 + -progenitor cells-derived DCs [22][23][24]. While such assays might provide in vivo like environments, primary cells are not well suited for a high-throughput format considering both donor-dependent variations as well as ethical aspect of such cell sources. Furthermore, previous efforts within in vitro assay development for sensitization that rely on full genome analysis have used a limited set of testing compounds.
The present study utilized in all 40 compounds and efforts were made to divide these compounds into two subsets, for training and testing respectively. While these experiments have resulted in successful predictions (data not shown), it is our experience that sensitizing compounds differ greatly in their induced gene expression profile, as can be seen in Figure 3D. In this perspective, we strived to include as many training compounds as possible when identifying our Prediction Signature, and did not exclude any compounds for validation. Instead, we validated the method by which the Prediction Signature was identified, by subdividing the samples into training and test sets at random, using unseen data for validation, to avoid overfitting. At present, the Prediction Signature consists of 200 transcripts, based on Figure 4A. Continuing the elimination process beyond 200 transcripts causes loss of information, as seen by the rise of KLD. Experiments have shown that correct classifications are possible even with further reduced signatures, down to 11 genes (data not shown). A reduction of signature size could be assessed in conjunction with validation of the assay, using untested positive and negative compounds in a new test set. By reducing the signature size at this point, the risk of biasing the signature towards this data set increases, making it harder to correctly classify unknown samples. Additional test compounds will also serve to assess the frequency of extreme transcriptional profile outliers, such as Oxazolone and Cinnamic aldehyde, which had to be removed from the analysis performed in this study. A number of reasons may be attributed to the fact that these compounds were not compatible with the assay, such as solubility in the cell media or extreme toxic effects. In those cases, other in vitro alternatives may complement this assay, so that the safety assessment of chemicals for sensitization includes a battery of in vitro assays. Naturally, an additional data set with blinded compounds is essential to validate whether the assay truly performs as estimated by the random subdivisions into training and test sets.
Of note, our Prediction Signature is able to predict the potency of sensitizing compounds, as defined by the LLNA (Figure 4C). However, the potency predicted by LLNA and that of our classifier do not match for all samples. Notably, the moderate sensitizer 2-hydroxyethyl acrylate showed resemblance to strong and extreme sensitizers with respect to gene expression profile. Similarly, the moderate sensitizers ethylendiamine, hexylcinnamic aldehyde, and glyoxal grouped together with weak sensitizers. These findings support the fact that sensitizing potency, as defined, may need revising.
By studying the identity of the transcripts and their involvement in intracellular signaling pathways, we were also able to confirm the biological relevance of the Prediction Signature. Using IPA, we found that the most highly populated pathways were nuclear factor-erythroid 2-related factor 2 (NRF2) mediated oxidative response, xenobiotic metabolism signaling, protein ubiquitination pathway, LPS/IL-1 mediated inhibition of Retinoic X receptor (RXR) function, aryl hydrocarbon receptor (AHR) signaling, and protein kinase A (PKA) signaling. These pathways are all known to take part in reactions provoked by xenobiotics, and several were associated with oxidative stress. Furthermore, Toll-like receptor (TLR) signaling is among the top pathways found in IPA. Recent studies on assay development for prediction of sensitization in vitro have to a large extent focused on how danger signals are provided to antigen-presenting cells, inducing pro-inflammatory cytokines and chemokines, as well as co-stimulatory molecules needed for a specific T-cell response. We hypothesize that these signals are provided through the innate immune responses, in analogy with infections, as reviewed in [25].
The primary pathways found in this study involved NRF2 signaling. This is a pathway activated by Reactive Oxygen Species (ROS), and is a defense mechanism to xenobiotics and response to cellular stress. In the resting cell, NRF2 is bound by kelch-like ECH-associated protein 1 (KEAP1) and located in the cytosol. In the response to ROS activity, KEAP1 is targeted for ubiquitination and protesomal degradation, resulting in the translocation of NRF2 to the nucleus, where it activates transcription of genes containing anti-oxidant response elements (ARE) in their promoter region [26]. The functions of genes transcribed by NRF2 association to ARE include regulation of inflammation, migration of DC and anti-oxidant defense enzymes, such as NADPH quinone oxidoreductase 1 (NQO1) and glutathione S-transferases (GST) [27,28], genes found in the Prediction Signature. Furthermore, the NRF2/KEAP1/ARE pathway has previously been described as activated in response to skin sensitizers, inducing maturation of dendritic cells [29].
Similarly, AHR is a transcription factor in the cytosol that is activated by binding to ligands, which includes a wide range of xenobiotic chemicals, such as halogenated aromatic hydrocarbons, polyphenols and a number of pharmaceuticals [30]. In the absence of a ligand, AHR is bound by a complex of chaperon proteins, keeping it in the cytosol. Upon ligand binding, AHR is translocated to the nucleus, where it dimerizes with aryl hydrocarbon receptor nuclear translocator (ARNT) [30]. The ARNT/ AHR heterodimer then binds to xenobiotic response elements (XRE) in promoter regions of target genes. The typical target genes for XRE include enzymes for drug metabolism, such as the cytochrome P450 (CYP) superfamily, as well as cytoprotective enzymes mediating defense against oxidative stress, such as NQO1 [31]. Interestingly, while NQO1 is under control of both NRF2 and AHR, with both ARE and XRE in the promoter region, it has also been shown that AHR is among the target genes for the activated NRF2 pathway and vice versa [32]. Thus, a battery of protective enzymes are induced in response to a variety of xenobiotics, possibly through a number of signaling pathways, ultimately leading to the maturation of dendritic cells, as also indicated by the present data. The protein ubiquitination pathway is involved in degradation of short-lived or regulatory proteins involved in many cellular processes, such as the cell cycle, cell proliferation, apoptosis, DNA repair, transcription regulation, cell surface receptors and ion channels regulation, and antigen presentation. Of note, both NRF2 and AHR are in the resting cell bound by proteins that are targeted for ubiquitination upon ligand binding.
RXR is a nuclear receptor, with retinoic acid as the most prominent natural ligand [33]. It has previously been described as important for xenobiotics recognition and glutathione homeostasis, with cytoprotective enzymes as target genes [34,35].
TLR signaling is known to play a major role in dendritic cell maturation, as they activate transcription of a number of pro-inflammatory cytokines, chemokinereceptors for homing to lymph nodes and co-stimulatory molecules [36][37][38]. While TLR6 and TLR9 are present in our Prediction Signature, others have reported TLR4 as a crucial mediator of contact allergy to nickel [39]. As these receptors all signal through nuclear factor kappa-light-chain-enhancer of activated B cells (NF-B), it is not surprising that different compounds activate different receptors, considering the chemical diversity of skin sensitizers, as discussed above.
Lastly, PKA signaling is a vastly versatile pathway activated by numerous stimuli, and, to the best of knowledge, this pathway has not previously been reported in association with skin sensitization. However, individual species of CYPs are known to be phosphorylated by PKA, in response to elevated levels of cyclic adenosine monophosphate (cAMP), triggered by xenobiotics. In addition, cAMP levels influence the nuclear translocation of AHR, connecting these two pathways and their impact on CYP activity [40].

Conclusion
In this paper, we have demonstrated the predictive power of a genomic biomarker signature, which correctly classifies sensitizers and non-sensitizers. The biomarker signature was derived from the human DC-like cell line MUTZ-3, which was challenged with a panel of 40 reference chemical compounds. The biomarker genes were shown to be biologically relevant, as demonstrated by their involvement in cytoprotective mechanisms and pathways triggered by xenobiotic substances, supporting their relevance as predictor genes for skin sensitization. The findings reported in this paper might impact the development of in vitro assays for assessment of skin sensitization, which is crucial in order to replace the animal models currently in use.

Relative viability = fraction of viable stimulated cells fraction of viable unstimulated cells · 100
For toxic compounds, the concentration yielding 90% relative viability (Rv90) was used. For non-toxic compounds, a concentration of 500 μM was used. For nontoxic compounds that were insoluble at 500 μM in medium, the highest soluble concentration was used. For compounds dissolved in DMSO, the final concentration of DMSO in each well was 0.1%. The vehicle and concentrations used for each compound are listed in Table 4.

Chemical exposure of the cells
The human myeloid leukemia-derived cell line MUTZ-3 (DSMZ, Braunschweig, Germany) was maintained in α-MEM (Thermo Scientific Hyclone, Logan, UT) supplemented with 20% (volume/volume) fetal calf serum (Invitrogen, Carlsbad, CA) and 40 ng/ml rhGM-CSF (Bayer HealthCare Pharmaceuticals, Seattle, WA), as described [10]. Cultures were maintained at 200.000 cells/ml during expansion, with a media change every 3-4 days. No differentiating steps were performed. Instead, the proliferating progenitor MUTZ-3 was used for stimulations, as delivered by the supplier. Prior to each experiment, the cells were immunophenotyped using flow cytometry as a quality control. Cells were seeded in 6-well plates at 200.000 cells/ml. Stock solutions of each compound were prepared in either DMSO or distilled water, and were subsequently diluted so the in-well concentrations corresponded to the Rv90 value, and in-well concentrations of DMSO were 0.1%. Cells were incubated for 24 h at 37°C and 5% CO 2 . Thereafter, cells were harvested and analyzed by flow cytometry. In parallel, harvested cells were lysed in TRIzol reagent (Invitrogen) and stored at -20°C until RNA extraction. Stimulations with chemicals were performed in three individual experiments, so that triplicates samples were obtained.

Phenotypic analysis with flow cytometry
All cell surface staining and washing steps were performed in PBS containing 1% BSA (w/v). Cells were incubated with specific mouse mAbs for 15 min at 4°C. The following mAbs were used for flow cytometry: FITC-conjugated CD1a (DakoCytomation, Glostrup, Denmark), CD34, CD86, and HLA-DR (BD Biosciences), PE-conjugated CD14 (DakoCytomation), CD54 and CD80 (BD Biosciences). Mouse IgG1, conjugated to FITC or PE were used as isotype controls (BD Biosciences) and PI was used to assess cell viability. FACS-Diva software was used for data acquisition with FACSCanto II instrument (BD Bioscience). 10,000 events were acquired and gates were set based on light scatter properties to exclude debris and nonviable cells. Further data analysis was performed using FCS Express V3 (De Novo Software, Los Angeles, CA).

Preparation of cRNA and gene chip hybridization
RNA isolation and gene chip hybridization was performed as described [41]. Briefly, RNA from unstimulated and chemical-stimulated MUTZ-3 cells, from triplicate experiments, were extracted and analyzed. The preparation of labeled sense DNA was performed according to Affymetrix GeneChip™ Whole Transcript (WT) Sense Target Labeling Assay (100 ng Total RNA Labeling Protocol) using the recommended kits and controls (Affymetrix, Santa Clara, CA). Hybridization, washing and scanning of the Human Gene 1.0 ST Arrays were performed according to the manufacturer's protocol (Affymetrix). The microarray data have been deposited in the Array Express database http://www.ebi. ac.uk/arrayexpress/ with accession number E-MTAB-670.

Microarray data analysis and statistical methods
The microarray data were normalized and quality checked with the RMA algorithm, using Affymetrix Expression Console (Affymetrix). Genes that were significantly regulated when comparing sensitizers with non-sensitizers were identified using one-way ANOVA, with false discovery rate (FDR) as a correction for multiple hypothesis testing. In order to reduce the large number of identified significant genes, we applied an algorithm developed in-house for Backward Elimination of analytes [42]. With this method, we train and test a Support Vector Machine (SVM) model [12] with leave-one out cross-validation, with one analyte left out. This process is iterated until each analyte has been left out once. For each iterative step, a Kullback-Leibler divergence (KLD) is recorded, yielding N KLDs, where N is the number of analytes.
The analyte that was left out when the smallest KLD was observed is considered to provide the least information in the data set. Thus, this analyte is eliminated and the iterations proceed, this time with N-1 analytes. In this manner, the analytes are eliminated one by one until a panel of markers remain that have Chlorobenzene CB - [14] Diethyl phthalate DP - [48] Dimethyl formamide DF - [46] Ethyl vanillin EV - [52] Glycerol GLY - [48] Isopropanol IP - [48] Lactic acid LA - [ [55] List of sensitizers and non-sensitizers used in assay development. 1) HMT, Human Maximization Test; HPTA, Human Patch Test Allergen. Information is derived from [17]. 2) False positives in LLNA.
been selected based on the ability of each analyte to contribute with orthogonal information for the discrimination of skin sensitizers vs. non-sensitizers. The selected biomarker profile of 200 transcripts were designated the "Prediction Signature". The scripts for Backwards Elimination and Support Vector Machines were programmed for R [43], with the additional package e1071 [44]. ANOVA analyses and visualization of results with Principal Component Analysis were performed in Qlucore Omics Explorer 2.1 (Qlucore AB, Lund, Sweden). Hierarchical clustering for the heatmap was performed in R. List of concentrations and vehicles used for each testing compound. 1) Kathon CG is a mixture of the compounds MC and MCI. The concentration of this mixture is given in %.

Interrogation of the method for identification of the Prediction Signature
The data set was divided into a training set and a test set, consisting of 70% and 30%, of the chemical compounds, respectively. The division was performed randomly, while maintaining the proportions of sensitizers and non-sensitizers in each subset at the same ratio as in the complete data set. A biomarker signature was identified in the training set, using ANOVA filtering and Backward Elimination, as described above. This test signature was used to train an SVM, using the training set, which was thereafter applied to predict the samples of the test set. The process was repeated 20 times and the distribution of the area under the Receiver Operating Characteristic (ROC) curve [45] was used as a measurement of the performance of the model.

Assessment of biological functions of Prediction Signature using pathway analysis
In order to investigate the biological functions the gene profile of the 200 genes derived from the Backward Elimination was analyzed, using the Ingenuity Pathway Analysis software, IPA, (Ingenuity Systems, Inc. Mountain View, USA). The gene profile was analyzed using the 'Build' and 'Path Explorer' functions to build an interactome of the core genes from the Prediction Signature together with connecting molecules, as suggested by IPA. The molecules of the signature were connected using the shortest known paths. In this process only human data from primary cells, cell lines and epidermal tissue was used. Public identifiers were used to map genes in IPA. All molecules except for endogenous and chemical drugs were allowed in the network and all kinds of connections were allowed. Known 'Functions' and 'Canonical Pathways' from IPA were mapped to the signature using the 'Overlay' function. The most densely populated pathways and functions were reported. All were significant, using the built in IPA statistical measures (p-values for functions and -log(p-values) for pathways).