Mechanistic Interrogation of Cell Transformation In Vitro: The Transformics Assay as an Exemplar of Oncotransformation

The Transformics Assay is an in vitro test which combines the BALB/c 3T3 Cell Transformation Assay (CTA) with microarray transcriptomics. It has been shown to improve upon the mechanistic understanding of the CTA, helping to identify mechanisms of action leading to chemical-induced transformation thanks to RNA extractions in specific time points along the process of in vitro transformation. In this study, the lowest transforming concentration of the carcinogenic benzo(a)pyrene (B(a)P) has been tested in order to find molecular signatures of initial events relevant for oncotransformation. Application of Enrichment Analysis (Metacore) to the analyses of the results facilitated key biological interpretations. After 72 h of exposure, as a consequence of the molecular initiating event of aryl hydrocarbon receptor (AhR) activation, there is a cascade of cellular events and microenvironment modification, and the immune and inflammatory responses are the main processes involved in cell response. Furthermore, pathways and processes related to cell cycle regulation, cytoskeletal adhesion and remodeling processes, cell differentiation and transformation were observed.


Introduction
Cancer is a consequence of a complex process requiring the alteration and modification of several biological traits that includes changes at the molecular, cellular, tissue and organ level. In humans, this process takes years to become phenotypically evident and the diagnosis is often made when the organ function is already impaired.
Addressing this complexity by using appropriate experimental in vitro models that can resemble the biological process is still a challenge.
Cancer development has been associated with chemical exposure since Sir Percival Pott described, in 1775, the occurrence of scrotum cancer in chimney sweepers [1,2]. This first observation paved the way to the discovery of several classes of chemicals acting as carcinogens in human, most (but not all) of which have been demonstrated to act through genotoxic mechanisms, involving a direct damage to DNA [2]. At the beginning of the 21st century, the identification of cancer hallmarks added new perspectives in the understanding of cancer onset and development, following chemical exposure.

Statistical Analysis of Microarray Experiment
The t-test unpaired (p ≤ 0.05, Benjamini-Hochberg) 0.02 µg/mL vs. 0.1% DMSO was performed using GeneSpring software, this indicated an extensive gene perturbation with B(a)P treatment on BALB/c 3T3 cells. Approximately 800 differently expressed genes were identified, among which 459 genes were over-expressed and 347 under-expressed. No fold change cut-off was applied.

Biological Interpretation of Microarray Experiments
The GeneSpring output dataset of 800 significantly modulated genes were analyzed further through the MetaCore V 6.34 and Single Experiment Workflow, which resulted in a total of 577 modulated genes included in the Metacore gene set analysis, corresponding to 667 network objects. Application of MetaCore Enrichment Analysis resulted in several Pathway Maps, Process Networks and Gene Ontology (GO) cellular processes with a statistically significant false discovery rate (FDR ≤ 0.05) (Tables 1-3).

Statistical Analysis of Microarray Experiment
The t-test unpaired (p ≤ 0.05, Benjamini-Hochberg) 0.02 μg/mL vs. 0.1% DMSO was performed using GeneSpring software, this indicated an extensive gene perturbation with B(a)P treatment on BALB/c 3T3 cells. Approximately 800 differently expressed genes were identified, among which 459 genes were over-expressed and 347 under-expressed. No fold change cut-off was applied.

Biological Interpretation of Microarray Experiments
The GeneSpring output dataset of 800 significantly modulated genes were analyzed further through the MetaCore V 6.34 and Single Experiment Workflow, which resulted in a total of 577 modulated genes included in the Metacore gene set analysis, corresponding to 667 network objects. Application of MetaCore Enrichment Analysis resulted in several Pathway Maps, Process Networks and Gene Ontology (GO) cellular processes with a statistically significant false discovery rate (FDR ≤ 0.05) (Tables 1-3).
The statistically significant data were quantitatively dominated by the involvement of immune response and inflammatory processes. Additional key pathways and processes that were statistically significant were related to cell cycle regulation (proliferation/apoptosis), cytoskeletal adhesion and remodeling processes, cell differentiation and transformation ( Figure S1). The statistically significant data were quantitatively dominated by the involvement of immune response and inflammatory processes. Additional key pathways and processes that were statistically significant were related to cell cycle regulation (proliferation/apoptosis), cytoskeletal adhesion and remodeling processes, cell differentiation and transformation ( Figure S1).

Pathway and Processes Related to Immune Response
The immune response was characterized by the involvement of several cytokinedependent pathways, interconnected by common modulated genes (Tables 1 and 2, Figure 3). Specific pathways of high relevance include: "Immune response_IFN-alpha/beta signaling via JAK/STAT", "Immune response_IFN-alpha/beta signaling via MAPKs" and "Immune response_IFN-alpha/beta signaling via PI3K and NF-kB pathways". These pathways are generally activated as a consequence of alpha and beta interferons (IFNs) binding to the respective transmembrane receptors, which, through the activation of intermediate transcription factors such as STAT and IRF, determine inflammatory, immune and proliferative processes (Figure 3). To better address the characterization of single genes, the gene set "Process Network: Inflammation_Interferon signaling" was scanned (Table S1). Table 1. The t-test unpaired (p ≤ 0.05, Benjamini-Hochberg) 0.02 µg/mL treatment vs. vehicle 0.1% DMSO was performed using GeneSpring software. No Fold Change cutoff was applied. In total, 800 significantly modulated genes were analyzed through the MetaCore Single Experiment Workflow. The first 50 statistically significant Pathway Maps in the order of significance are listed (FDR ≤ 0.05). "Gene Ratio" indicates the number of significantly altered genes matching the objects in those specific pathways, out of the total number of genes involved in different specific pathways of the MetaCore database. Finally, up-regulated and down-regulated endpoints indicate individual and/or family genes positively or negatively altered in those specific pathways.      Table 3. Significantly modulated genes resulted from GeneSpring analysis were analyzed further through the MetaCore Single Experiment Workflow. The first 10 statistically significant Gene Ontology (GO) cellular processes are listed (FDR ≤ 0.05). "Gene Ratio" indicates the number of significantly altered genes matching the objects in those specific GO Process out of the total number of genes involved. The treatment with B(a)P appeared to induce the activation of the complement system, as highlighted by the significant modulation of the related Pathway Maps, in which a considerable over-expression of the C3 gene is shown (Figure 3).

GO
Furthermore, other immune response modulated pathways included interleukindependent pathways, associated with pro-inflammatory functions. Of particular interest is the over-expression of the interleukins Il-22 and Il-18 (Table 1; Pathway n. 15; 47) ( Figure 3).
Up-regulation of some negative regulators of cytokine signaling and immune responses, which prevent genes from being expressed, were also detected. In particular, the over-expression of Nfkbia, which encodes a member of the NF-kappa-B inhibitor family (I-kb), and Socs3, which encodes the suppressor of cytokine signaling member 3, were highlighted ( Figure 3). Both these genes may act in a classic negative feedback loop to regulate cytokine and JAK-STAT signal transduction [28,29].
Process Networks and GO processes analysis also confirmed the importance of the cytokine-mediated inflammation signals (Tables 2 and 3).
Among the different significantly modulated Process Networks (p-values and FDR provided in Table 2), the "Inflammation inflammasome" was of high interest, with the conspicuous presence of INF-related genes in the differential expressed gene list (Table S2). It is also worth noting the Chil1 (Chitinase-like 1, also known as Chi3l1) differential gene expression (FC: +20.65; Table 3: GO processes i.e., n.1-3-5). This gene modulation has been already proposed in other studies as a potential link between the molecular microenvironment related to inflammation and cells transformation [30,31]. Moreover, this gene has been shown to be secreted by Cancer Associated Fibroblasts (CAFs) in breast carcinomas during both early and advanced stages, and to play a central role in tumor promotion, producing a pro-inflammatory and immunosuppressed microenvironment [30,31].

Pathways and Processes Related to EMC Modification and EMT/MET
Several statistically significant pathways and processes were linked to the "Epithelial to Mesenchymal Transition" (EMT), to processes of differentiation and modification of extracellular matrix (ECM) and cells adhesion (Tables 1-3).
In addition, the modulation of "Signal Transduction_TGF-beta, GDF and Activin signaling" and "Signal transduction_NOTCH signaling" are presented in the Process Networks list (   The treatment with B(a)P appeared to induce the activation of the complement s tem, as highlighted by the significant modulation of the related Pathway Maps, in wh a considerable over-expression of the C3 gene is shown (Figure 3).
Furthermore, other immune response modulated pathways included interleukin-d pendent pathways, associated with pro-inflammatory functions. Of particular interes . Canonical and non-canonical interferon-mediated signaling pathways. Type 1 interferons (IFN) can trigger signals through the JAK-STAT-mediated canonical pathway and the formation of transcriptional-activator complexes, which translocate into the nucleus and activate the interferonstimulated response element (ISRE) or IFN-γ-activated site (GAS). The activation of ISRE leads to the transcription of IL-1β, which, in turn, activates the complement cascade. Type 1 IFNs, specifically IFNα, stimulate the activation of IL-22, through IL-22R, which can efficiently stimulate STAT1 and its downstream pathway. IFNs' non-canonical pathways include IFNα/B signaling via MAPKs pathway, involved in chromatin rearrangement, and IFNα/B signaling via PI3K and NF-kB pathway, involved in survival signals, angiogenesis and cell cycle. The interplay with Toll like receptor (TLR)-mediated pathway is already shown. The activation of T:R by pathogen-associated or damage-associated molecular patterns (PAMPs, DAMPs) triggers the pro-inflammatory NF-kB-mediated signaling pathway, which is inhibited by SOCS3, through the myddosome degradation. Toll-IL1-R receptor activates IL18, which support the activation of NF-kB pathway and interplays with IL-22 to activate MAPK-signaling pathway.
Furthermore, the up-regulation of Dcn (Decorin) ( Table 2 and Table S4.) is proposed to inhibit the TGF-β.
The down-regulation of Wisp1 and Jag1 could be also evidence of NOTCH1 signaling pathway depression, an event which may occur in sustaining fibroblast proliferation [33,34].
Moreover, the up-regulation of Lif (leukemia inhibitory factor, FC +2.55) ( Table 1; pathway n. 40) is a cytokine which belongs to the IL-6 superfamily and has been noted to be an important regulatory factor [35][36][37].

Genes Involved in the Metabolism of Benzo (a) Pyrene
Extensive literature evidence defines B(a)P (also 3-MCA and 2,3,7,8-tetrachlorodibenzop-dioxin) as a bifunctional AhR inducer [38,39], which, upon receptor ligand binding and heterodimerization with the heterodimeric transporter protein ARNT, is able to induce both Phase 1 oxidative enzymes (such as the cytochrome P450 enzymes CYP1A1 and CYP1B1) and Phase 2 conjugating enzymes (such as glutathione transferase and UDPglucuronosyltransferase [39][40][41] In this study, a remarkable up-regulation of the Cyp1a1 (FC = +8.17) and Cyp1b1 (FC = +2.76) genes were observed. Positive regulation of the Cyp1a1 and Cyp1b1 genes was also identified in the significantly modulated Process Network "Signal transduction ESR1 nuclear pathway" ( Table 2; Process Network n. 25), as was Ugt1a6b over-expression; this gene encodes for UDP-glucuronosyl-transferase (UGT).
With regard to inducible enzymes, the up-regulation of Ncf4 (Neutrophil cytosolic factor 4, also known as P40phox, FC: +5.070; Table 2; Process Network n. 16) is notable. As a member of the NADPH oxidase subunit, this enzyme plays an important role in NADPH enzyme regulation of oxidase activity and the subsequent production of reactive oxygen species (ROS). A direct transcription regulation of P40phox by AhR has been proposed [42].

Discussion
Here, we have utilized the "Transformics Assay" to develop an evaluation and analysis of the effect of the lowest transforming concentration of B(a)P using the BALB/c 3T3 A31-1-1 cell line. In this CTA, the 0.02 µg/mL B(a)P concentration resulted in the lowest concentration related to the foci formation.
The study results reported here consolidate upon those previously reported in the BALB/c 3T3 model validation [43], where the validation report showed an increment in malignant foci formation at the concentration range 0.01-0.02 µg/mL.
The study aimed to identify the molecular and cellular events at the tested threshold concentration following 72 h of exposure.
A deep perturbation of a large set of immune and inflammatory pathways and processes was observed, showing cross-talk interactions and reciprocal control, essentially 'immune-editing'. Overall, the data obtained are consistent with principal key events related to the biological process "inflammation", regardless of the tissue type, specifically (i) cell activation and phenotypic modification, which includes alterations in their secretory activity, activation of biosynthesis pathway, production of pro-inflammation molecules and changes in metabolism and morphology; (ii) increase of molecules with pro-inflammatory action; (iii) recruitment/activation of leukocytes [44] (see Figure S1 for further information about some involved Differential Expressed Genes (DEGs)).
It is postulated that the first steps in the inflammation signaling could be: mitogen activated kinases (MAPKs) activity and Pattern Recognition Receptors (PPRs) signaling activation [45]. PPRs have been reported to recognize both pathogens (PAMPs, Pathogen-Associated Molecular Patterns) and endogenous ligands released by distressed or damaged cells, generally called DAMPs (Damage-Associated Molecular Patterns) [46,47].
This result is consistent with other studies which investigated some B(a)P (and other AhR ligands) mechanisms of toxicity mediated by cross-talk between AHR and MAPK [48][49][50].
Many of the inflammation related pathways and processes observed are associated with a "cellular antiviral state", highlighting overlaps in some sterile and non-sterile stress responses. Indeed, it has been proposed that DNA damage response sensors can mediates the type I interferon (IFN) response and the inflammasome activation [51].
Consistent with this assumption and the known B(a)P genotoxicity, the "Transcrip-tion_P53 signaling pathway" (Table 1 pathway n. 39) can be highlighted.
Moreover, in this study, the modulation of the "Immune response_TLR5, TLR7, TLR8 and TLR9 signaling pathways (Toll-like receptors) ( Table 1, Pathway n. 42) and the upregulation of the Rig-I (Retinoic acid-inducible gene 1) together with the modulation of the Process Networks "inflammation inflammasome" and "inflammation interferon signaling" have been observed ( Table 2, Tables S1 and S2 and Figure 4).
The activation of such Pattern Recognition Receptors (RIG-I and TLRs) could represent the initiation step for IRF-and NF-kB-dependent signals activation, leading to the release of type 1 interferons and other cytokines and, consequently, to the JAK/STAT signaling cascade, leading to transactivation of ISG (interferons-stimulated genes) promoters, generating the INF response ( Figure 4).
It is also worth noting the concomitant up-regulation of RIG-I and some ISGs (such as Oas1a), which have been linked to a positive-feedback/feed-forward mechanisms and the formation of an auto-amplification loop RIG-I/IRF [29,49]. Moreover, the over-expression of RIG-I has been connected to the relative inflammasome activation and mature proinflammatory interleukins release (such as IL-18) [52,53].
As the results show, several significantly modulated pathways are either directly or indirectly related to the cell transformation process, including EMT-related pathway modulation, ECM molecular signaling remodeling and TGF-b and WNT signaling modulations.
Results from the Enrichment Analysis, as well as the directly observable morphological change in spindle-shaped cells, suggest that fibroblasts are undergoing several modifications, and there are some overlaps with the "CAFs phenotype" [36,37]. Taken in combination, these background signaling events can be considered a composite molecular signal of ECM remodeling, as indicated by the modulation of the relative Pathway maps (Table 1; pathways n 5;10) and Process Networks (Tables 2 and S4.) Moreover, the scored Lif gene up-regulation has been noted to be associated with a constitutive activation of the JAK1/STAT3 signaling pathway, this sustains the contractile and pro-invasive abilities of fibroblasts [35][36][37].
Other gene modulations suggest the inhibition of NOTCH 1 and TGF-β pathways and the deregulation of the WNT signaling pathway (Figure 4). The WNT and TGF-β pathways are both directly and indirectly involved in tissue homeostasis processes of self-renewal and proliferation, differentiation, cell movement and adhesion [55][56][57][58].
The product of Wnt-5b, which is up-regulated in our cell model, is a ligand regulating fibroblasts activation via non-canonical signaling. It occurs in inflammation responses mediating IL-6 and CXCL8 release and JNK, p38 and p65 NF-κB activation [59].
While the putative down-regulation of TGF-β signaling has been previously considered controversial, with respect to the classical "molecular phenotype" of myofibroblasts related to CAFs, more recently, Biffi et al. [36] have highlighted differences in CAFs subpopulation in pancreatic ductal adenocarcinoma, discerning between two subtypes characterized by either myofibroblastic or inflammatory phenotypes, with the latter expressing Lif and JAK/STAT activation [36]. Interestingly, the LIF/STAT3 and WNT/β-catenin, as well as the inhibition of TGF-β pathways, are also all involved in the maintenance of pluripotency [54]. Considering the final end-point of oncotransformation and not forgetting that EMT can be triggered in a tissue-specific manner, we consider that these combined molecular data illustrate the overall molecular dynamics of the "evolving tumor-microenvironment".
From this work, we can observe that the AhR-mediated signaling pathway has a major role in triggering the molecular and cellular response to B(a)P in this cell model. The AhR transactivation is supported by the noticeable over-expression of phase I Cyp1a1 and Cyp1b1 genes and of the p40phox enzyme [23,26,42].
The AhR mediates the xenobiotic genotoxicity, but it is also functionally involved in cell adhesion, spreading, migration and modulation of ECM composition, either in physiological cell conditions or following the exposure to xenobiotics, and in particular to B(a)P While this is a very specific example that so far has not been reported for other cancer types, the fact that this has been observed contributes to an improved weight of evidence for the down-regulation of TGF-β. It is also noteworthy that TGF-β interplay with AhR signaling has been reported [60][61][62].
Interestingly, the LIF/STAT3 and WNT/β-catenin, as well as the inhibition of TGF-β pathways, are also all involved in the maintenance of pluripotency [54]. Considering the final end-point of oncotransformation and not forgetting that EMT can be triggered in a tissue-specific manner, we consider that these combined molecular data illustrate the overall molecular dynamics of the "evolving tumor-microenvironment".
From this work, we can observe that the AhR-mediated signaling pathway has a major role in triggering the molecular and cellular response to B(a)P in this cell model. The AhR transactivation is supported by the noticeable over-expression of phase I Cyp1a1 and Cyp1b1 genes and of the p40phox enzyme [23,26,42].
The AhR mediates the xenobiotic genotoxicity, but it is also functionally involved in cell adhesion, spreading, migration and modulation of ECM composition, either in physiological cell conditions or following the exposure to xenobiotics, and in particular to B(a)P [41]. Moreover, the AhR is involved in autophagy and EMT processes; these functions can be altered by AhR transmigration and depletion in cytoplasm [63,64].
Our data show the inflammatory/immune response after 72 h of exposure as the main cellular process (Tables 1-3 and Figure 4). Results reported in the paper from Mascolo et al. [23] showed a prolonged inflammation state in time. Particularly in this previous study, transforming concentration of 3-MCA treatment had caused the IL-6 production after 28 days from the cellular exposure (end time point of the CTA). Worth noting, the IL-6 cytokines could represent the downstream product of a number of modulated pathways in this B(a)P study. Otherwise, the IL-6 receptor up-regulation can be highlighted ( Table 1, Pathway n. 40). It is also worth noting that the inflammation initiation may influence the B(a)P metabolic fate supporting the activation of steps leading to genotoxic damage/injury. This hypothesis is supported by a number of reports showing that exposure to inflammation fostering molecules can slow down the B(a)P metabolism and generate higher levels of B(a)P-DNA adducts [24,[65][66][67].
Particularly the phase II enzymes, such as UGTs, are necessary for the transformation of non-polar, exogenous and endogenous compounds into their hydrophilic derivatives. With respect to B(a)P, the activity of UGTs is fundamental for preventing the formation of epoxide BPDE (B(a)P-7, 8-dihydrodiol-9, 10-oxide), a potent mutagenic metabolite [68,69]. Our group has previously noted that the loss in expression of UGTs associated with cell transformation is a crucial event to discriminate the adaptive from the adverse response that leads to transformation [23].
Unlike UGT1A7, UGT1A9 and UGT1A10, UGT1a6 (up-regulated in this study) is unable to detoxify the intermediate metabolite BPD (B(a)P-trans-7,8-dihydrodiol), and thus, is unable to prevent its oxidation and consequent epoxide formation [68,69]. Moreover, down-regulation of some UGTs (but not UGT1A6) have been observed in human hepatocarcinogenesis, while the up-regulation of Ugt1a6 has been observed in the rat hepatocarcinogenesis process [68,69]. Altogether, these data show adverse effects of carcinogens on wide-ranging aspects of the cellular system. Particularly the inflammation signals stimulation, together with oxidative stress effects and inhibition of healing, can build a self-sustaining inflammation system and a transforming microenvironment. At the same time, these immune responses may alter the phase I/phase II metabolism of the carcinogen. Metabolic derived debris (products of oxidative stress and/or unstable metabolites activity) may accumulate in cells and induce the innate immune response to DAMPs, forcing the response to become chronic and maladaptive [21].
According to the multiple roles of the AhR receptor, its transmigration and gene expression activation could be involved in the aberrant immune response [21,23,70], triggering a reactive metabolite cascade that leads to genotoxicity and also promotes fibroblast motility [41] (Figure 4).
However, further, more specific studies are needed to clarify the complex role of AhR in such processes and eventual ligand-specific effects.
Furthermore, with the identification of a LTC of 0.02 µg/mL B(a)P, we have a preliminary threshold for the CTA when identifying whether a result for an AhR bin of test chemicals is positive or negative, within an IATA.
Interestingly, in vivo, endpoints predictive of carcinogenicity induced by B(a)P are not the most sensitive toxicity outcomes, and the LTC reported here improves sensitivity for the prediction of oncotransformation, as compared to the RCB. For example, with B(a)P in vivo rodent studies, the critical point of departure of 0.025 mg B(a)P/kg-bw-day for neurotoxicity is used, not the lowest observed adverse level of 0.54 mg B(a)P/kg-bw-day for rodent forestomach carcinogenicity [71].

Experimental Design
set of plates for each treatment (NT, 0.1% DMSO, 0.02 µ g/mL B(a)P, 0.2 µ g/mL B(a)P, 4 µ g/mL 3-MCA) were maintained in culture for the CTA and for the RNA isolation. Cells were treated 24 h after seeding. Total RNA was isolated after 72 h of exposure. RNA was extracted from a set of plates treated with the Lowest Transforming Concentration (LTC) 0.02 µ g/mL B(a)P (16 total plates: 4 technical replicate and 4 biological replicates) and a set of plates treated with the vehicle 0.1% Dimethyl Sulfoxide (DMSO) (8 total plates: 2 technical replicates and 4 biological replicates), as control. The analysis of the entire transcriptome by microarrays has been performed. In parallel the CTA experiment was conducted and finalized 31-32 days after seeding, when cells were fixed and stained. Particularly, M10F medium was used for routine culture: MEM supplemented with 10% Fetal bovine serum FBS and 1% penicillin 10,000 U/mL/streptomycin 10 mg/mL; DF212F medium was used for the late stage of transformation assay: DMEM/F12 with 2 µg/mL insulin, 2% Fetal bovine serum FBS and 1% penicillin 10,000 U/mL/streptomycin 10 mg/mL solution. Ten replicates for each test sample were performed. (C) Integrated experimental protocol: Cells were treated 24 h after seeding. Total RNA was isolated after 72 h of exposure. RNA was extracted from cells treated with the Lowest Transforming Concentration (LTC) 0.02 µg/mL B(a)P (16 total plates: 4 technical replicate sand 4 biological replicates) and cells treated with the vehicle 0.1% Dimethyl Sulfoxide (DMSO) (8 total plates: 2 technical replicates and 4 biological replicates) as control. The analysis of the entire transcriptome by microarrays was performed. The microarray experiment was conducted with a glass slide containing 8 × 60 K formatted arrays (Agilent's Whole Mouse Genome Oligo microarray), four of which were hybridized with the treated cells lysate (four biological replicates), and four with the control lysate (three biological replicates and one technical replicate). The quality control analysis and data normalization were performed by Agilent Feature Extraction. The GeneSpring GX software (Agilent Technologies, Santa Clara, CA, USA) was employed for the statistical analysis. A preliminary concentration range finding CTA was performed to select the treatment concentrations to be tested in the Transformics Assay. The Transformics Assay experimental design is shown in Figure 5. Cell culture experiments were run in parallel. A set of plates for each treatment (NT, 0.1% DMSO, 0.02 µg/mL B(a)P, 0.2 µg/mL B(a)P, 4 µg/mL 3-MCA) were maintained in culture for the CTA and for the RNA isolation. Cells were treated 24 h after seeding. Total RNA was isolated after 72 h of exposure. RNA was extracted from a set of plates treated with the Lowest Transforming Concentration (LTC) 0.02 µg/mL B(a)P (16 total plates: 4 technical replicate and 4 biological replicates) and a set of plates treated with the vehicle 0.1% Dimethyl Sulfoxide (DMSO) (8 total plates: 2 technical replicates and 4 biological replicates), as control. The analysis of the entire transcriptome by microarrays has been performed. In parallel the CTA experiment was conducted and finalized 31-32 days after seeding, when cells were fixed and stained.

Cells
The original stock of BALB/c 3T3 cells, clone A31-1-1, was obtained from the Health Science Research Resource Bank (Osaka, Japan). Cells were tested, characterized and authenticated by the cell bank for the species origin.
Working cultures were expanded from the original cryopreserved stock. Two vials were thawed at the same time and mixed to create one cell suspension. Cells were grown in M10F, Minimum Essential Medium supplemented with 10% fetal bovine serum (FBS, SigmaAldrich®, St. Louis, MO, USA). Only subconfluent cells were used, and the target cells were not maintained beyond the third passage after thawing. Cultures were maintained in a humidified incubator with an atmosphere of 5% CO 2 in air at 37 • C. DMSO was used at a 0.1% concentration as the vehicle for the two chemicals and as the negative control.

Chemicals and Solutions
The stock solution of 3-MCA (4 mg/mL) was delivered to the culture medium to reach the final working concentration of 4 µg/mL.
The B(a)P working solutions were prepared by diluting sequential stocks solutions in the culture media to obtain concentrations ranging from 0.0025 µg/mL to 0.04 µg/mL, within the preliminary CTA, and the concentration of 0.02 µg/mL in the subsequent Transformics Assay.

CTA and Cytotoxicity Assessment
The test was performed according to the EURL-ECVAM (EU Reference Laboratory European Centre for the Validation of Alternative Methods, European Commission, Joint Research Centre, I-21027 Ispra (VA), Italy) validated protocol [27,72,73]. This protocol includes cytotoxicity testing within the CTA. Details are reported in Supplementary Information, available Online.
A preliminary CTA was performed in order to identify the LTC, by testing B(a)P at concentrations of 0.0025 µg/mL, 0.005 µg/mL, 0.01 µg/mL, 0.02 µg/mL and 0.04 µg/mL (Figure 2). While all the concentrations tested were cytotoxic, they did not all incur transformation. The LTC observed for B(a)P was 0.02 µg/mL, and therefore, this concentration was selected for testing in the CTA for the Transformics Assay. A concentration of 0.2 µg/mL B(a)P was used in the subsequent Transformics CTA, as a positive control.
Cells were treated 24 h after seeding. Following 72 h exposure, the medium was removed, and the cell cultures were incubated in a humidified incubator at 37 • C with 5% CO 2 . The time of exposure was chosen according to EURL-ECVAM validated protocol [27,72,73].

Total RNA Extraction
Total RNA was isolated after 72 h of exposure. RNA was extracted from cells treated with the Lowest Transforming Concentration (LTC) 0.02 µg/mL B(a)P (16 total plates: 4 technical replicate sand 4 biological replicates) and cells treated with the vehicle 0.1% Dimethyl Sulfoxide (DMSO) (8 total plates: 2 technical replicates and 4 biological replicates), as control. Total RNA was isolated by using TRIzol Reagent (Invitrogen, San Diego, CA, USA) followed by purification with Rneasy affinity column (Qiagen, Valencia, CA, USA), according to the manufacturer's instructions. RNA quality was assessed by Agilent Bioanalyzer 2100 (RNA 6000 Nano LabChip, Agilent Technologies, Palo Alto, CA, USA). Four type 1 biological replicates were obtained for each treatment (0.02 µg/mL B(a)P and 0.1% DMSO).

Total RNA Labeling and Hybridization
The cRNA was labeled, purified and hybridized on oligonucleotide slides (SurePrint G3 Mouse Gene Expression 8x60K Microarray Kit) using the Low Input Quick Amp Labeling Kit, version 6.9.1, December 2015 (Agilent Technologies, Santa Clara, CA, USA) (www.genomics.agilent.com HYPERLINK http://www.chem.agilent.com/, accessed on 13 May 2018). Four arrays were hybridized with the treated cells lysate (four biological replicates), and four with the control lysate (three biological replicates and one technical replicate).

Statistical Analysis of Microarray Data
Raw data were filtered for intensity and quality signal by using GeneSpring GX (Agilent Technologies). A t-test unpaired p ≤ 0.05 corrected by Benjamini-Hochberg, was used to select differentially expressed genes after B(a)P treatment in comparison with 0.1% DMSO. No fold change cut-off was applied.

Tools of Biological Interpretation
The list of differentially expressed genes was evaluated by using Metacore™ software V6.34 (Calrivate Analytics), an integrated 'knowledge-based' platform with a manually annotated database of protein interactions verified by 'small experiment' data and gene-disease associations (https://portal.genego.com/, accessed on 13 May 2018 ). An enrichment analysis was performed using the Analyze Single Experiment workflow, to identify Pathways Maps Process Networks and GO cellular processes altered after 72 h of 0.02 µg/mL B(a)P exposure.

Conclusions
We have demonstrated how the application of the Transformics method to the CTA can support the elucidation as to how the molecular environment of exposed cells can be resembled, and how it can enable the analysis of several molecular events involved in toxicity and cell transformation. By taking this approach, we have shown how the test method has the ability to highlight diverse mechanisms and their plausible interactions leading to cell malignant transformation. This type of approach could be suitable to enable an interlinked understanding of chemical specific modes of action along the key events of the NGTxC IATA process.
Some perceived limitations regarding this method can be ascribed to (i) species extrapolation and comparison with the rodent cancer bioassay and (ii) the variability of transcriptomics data and ways in which such data are reported.
For the former consideration, it is acknowledged that the traditional BALB/c 3T3 CTA showed a concordance of 68%, sensitivity of 75% and specificity of 53%, positive predictivity of 77% and false positive and negative rates of 47 and 25%, respectively, compared to the rodent bioassay [74]. However, it is important to emphasize that the toxicogenomic approach described herein can support the human relevance for data interpretation.
With regard to limitations that may derive from the transcriptomics approach, it is recognized that the transcriptome is of a highly dynamic nature with consequent difficulties in the biological interpretation and that cellular post-transcriptional events are lacking in representation. To overcome this, the Transformics Assay can be performed at any time during the in vitro transformation process, giving the possibility to compare different timepoint profiles [23], and furthermore, internationally agreed OECD transcriptomics reporting formats guidance is underway [75]. The latter underlines the wish of the international regulatory community to see greater applications of 'omics tools for regulatory purposes', such as that described here.
It is, therefore, proposed that the Transformics Assay can support the use of the CTA to model the carcinogenesis process in vitro within the IATA. Further studies will be conducted to confirm its potential for NGTxC hazard assessment [10].