{'fraction_non_alphanumeric': 0.028760172331258976, 'fraction_numerical': 0.018305409286740067, 'mean_word_length': 5.141345249294449, 'pattern_counts': {'":': 0, '<': 13, '<?xml version=': 0, '>': 6, 'https://': 0, 'lorem ipsum': 0, 'www.': 2, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
System-based Identification of Toxicity Pathways Associated With Multi-Walled Carbon Nanotube-Induced Pathological Responses

The fibrous shape and biopersistence of multi-walled carbon nanotubes (MWCNT) have raised concern over their potential toxicity after pulmonary exposure. As in vivo exposure to MWCNT produced a transient inflammatory and progressive fibrotic response, this study sought to identify significant biological processes associated with lung inflammation and fibrosis pathology data, based upon whole genome mRNA expression, bronchoaveolar lavage scores, and morphometric analysis from C57BL/6J mice exposed by pharyngeal aspiration to 0, 10, 20, 40, or 80 µg MWCNT at 1, 7, 28, or 56 days post-exposure. Using a novel computational model employing non-negative matrix factorization and Monte Carlo Markov Chain simulation, significant biological processes with expression similar to MWCNT-induced lung inflammation and fibrosis pathology data in mice were identified. A subset of genes in these processes was determined to be functionally related to either fibrosis or inflammation by Ingenuity Pathway Analysis and were used to determine potential significant signaling cascades. Two genes determined to be functionally related to inflammation and fibrosis, vascular endothelial growth factor A ( vegfa ) and C-C motif chemokine 2 ( ccl2 ), were confirmed by in vitro studies of mRNA and protein expression in small airway epithelial cells exposed to MWCNT as concordant with in vivo expression. This study identified that the novel computational model was sufficient to determine biological processes strongly associated with the pathology of lung inflammation and fibrosis and could identify potential toxicity signaling pathways and mechanisms of MWCNT exposure which could be used for future animal studies to support human risk assessment and intervention efforts.


Introduction
Nanotechnology is an emerging discipline in both industrial and medical fields, which necessitates the development of nanotoxicology to determine the biological effects of occupational and commercial nanoparticle exposure (Oberdorster et al., 2005). Multi-walled carbon nanotubes (MWCNT) are fibrous nanoparticles consisting of multiple concentric cylindrical carbon tubes that are appealing for both industrial and medical purposes due to their efficient electronic conductivity, great strength, and strong capillary forces, while maintaining a small size, light weight, high surface area to mass ratio, and low density (Iijima, 1991;Ajayan, 1999;Castranova, 2011). The physical attributes of MWCNT, while useful from an engineering standpoint, make them easily aerosolized and a potential inhalation hazard during synthesis, product use, and disposal. In vitro studies of MWCNT exposure determined toxicity to both lung epithelial and microvascular endothelial cells with increases in reactive oxygen species (ROS) production, NF-ƙB signaling, cytokine release, cytoskeletal reorganization, and endothelial cell permeability (Walker et al., 2009;Ye et al., 2009;He et al., 2011;Srivastava et al., 2011;Pacurari et al., 2012). Mouse and rat in vivo studies determined that MWCNT can reach the alveolar region of the lung after pharyngeal aspiration and inhalation, respectively, and induce a transient inflammatory reaction followed by a progressive fibrotic response (Muller et al., 2005;Mercer et al., 2010;Mercer et al., 2011).
Although chronic inflammation has been suggested as the underlying mechanism governing the progression to fibrosis, this does not appear to hold true for MWCNT as the initial inflammatory response to MWCNT diminishes before the progressive fibrotic response begins Porter et al., 2010;Mercer et al., 2011). Therefore, it is essential to uncover the significant biological processes directing MWCNT-induced inflammation and fibrosis so as to determine potential outcomes and hallmarks of exposure. We hypothesize that the identification of transcription-related biological processes and pathways, which match the patterns of BAL quantification  for inflammatory pathology and morphometric scoring of collagen (Mercer et al., 2011) for fibrosis in MWCNT-treated mice, could identify critical toxicity pathways and potential mechanisms of MWCNT-induced lung inflammation and fibrosis for early identification and intervention.
Recently, our group conducted an in vivo dose-response time-course study of MWCNT exposure in C57BL/6J mice to determine the ability of MWCNT to induce pulmonary inflammation, damage, and fibrosis . Mice were exposed to 0, 10, 20, 40, or 80 µg of MWCNT by pharyngeal aspiration with endpoints monitored at 1, 7, 28, and 56 days post-exposure . The results indicated that a transient inflammatory response occurred 1 day post-exposure with peak activity 7 days post-exposure. A fibrotic response was noted 28 days post-exposure, which progressed through 56 days post-exposure . Nevertheless, MWCNT-induced toxicity pathways and mechanisms underlying these observed in vivo pathological responses remain unknown. We hypothesize that systematic analyses of gene expression profiles and pathological data could identify transcription-related biological processes correlated with the observed pathological patterns of lung inflammation and fibrosis, which could reveal MWCNT-induced toxicity pathways and pathogenesis. The current study sought to use a novel computational system to identify transcription-related biological processes and pathways associated with these MWCNTinduced pathology responses in a comprehensive systematic evaluation. A novel computational model, previously reported by Dymacek (2011) was applied to genome-wide mRNA expression profiles and pathological analysis of mouse lungs taken at these respective time points so as to determine biological processes significantly correlated with inflammation (bronchoalveolar lavage fluid [BAL] score ) or fibrosis (morphometric analysis of alveolar interstitial fibrosis (Mercer et al., 2011)). These biological processes were then analyzed through Ingenuity Pathway Analysis (IPA) to determine gene subsets functionally related to inflammation or fibrosis. In vitro gene and protein expression data of two genes functionally related to inflammation and fibrosis, vascular endothelial growth factor A (vegfa) and C-C motif chemokine 2 (ccl2), were validated through cell culture studies.
This study determined that a novel computational model was sufficient to identify transcription-related biological processes strongly associated with lung inflammation BAL scores and fibrosis morphometric analysis. Potential toxicity signaling pathways of MWCNT exposure were determined and validated in vitro. The use of these toxicogenomics data and in vivo animal-model based gene expression profiling integrated with in vitro verification may allow for successful toxicity profiling of MWCNT as well as the identification of potential signaling pathways involved in the etiology of MWCNT-induced injury.

MWCNT
MWCNT used in both mouse and cell studies were obtained from Mitsui & Company (MWCNT-7, lot #05072001K28) and have been previously characterized . Briefly, the bulk MWCNT exhibit a distinctive crystalline structure with the number of walls ranging from 20 to 50 walls. Overall, MWCNT trace metal contamination was 0.78%, including sodium (0.41%) and iron (0.32%) with no other trace metal contamination over 0.02%. Transmission electron microscopy (TEM) micrographs of MWCNT dispersed in dispersion medium (DM) demonstrated that DM promotes significant dispersion of MWCNT. The quantitative analysis of TEM micrographs revealed that the median length of the MWCNT sample was 3.86 µm (GSD 1.94) and the count mean width was 49 ± 13.4 (SD) nm. The zeta potential of the MWCNT in the DM was determined to be −11mV.

Animals
Animal studies were performed as previously described . Briefly, male C57BL/6J mice (7 weeks old) were obtained from Jackson Laboratories (Bar Harbor, ME). Individual mice were housed one per cage in polycarbonate isolator ventilated cages and provided HEPA-filtered air with fluorescent lighting from 0700 to 1900 hours. Autoclaved Alpha-Dri virgin cellulose chips and hardwood Beta-chips were used as bedding. Mice were monitored to be free of endogenous viral pathogens, parasites, mycoplasms, Helicobacter, and CAR Bacillus. Mice were maintained on Harlan Teklad Rodent Diet 7913 (Indianapolis, IN) and tap water was provided ad libitum. Animals were allowed to acclimate for at least 5 days before use. All animals in this study were housed in an AAALAC-accredited, specific pathogen-free, and environmentally controlled facility. All animal studies and procedures were approved by the National Institute for Occupational Safety and Health ACUC.

MWCNT Pharyngeal Aspiration Exposure
Suspensions of MWCNT were prepared in DM and administered as previously described (Porter et al., 2008;Porter et al., 2010). In brief, each treatment group consisted of 8 mice, which were anesthetized with isoflurane (Abbott Laboratories, North Chicago, IL). When fully anesthetized, the mouse was positioned with its back against a slant board and suspended by the incisor teeth using a rubber band. The mouth was opened and the tongue gently pulled aside from the oral cavity. A 50 µl aliquot of sample was pipetted at the base of the tongue, and the tongue was restrained until at least 2 deep breaths were completed (but not for longer than 15 seconds). Following release of the tongue, the mouse was gently lifted off the board, placed on its left side, and monitored for recovery from anesthesia. Mice received either DM (vehicle control), or 10, 20, 40 or 80 µg MWCNT suspended in DM .

Tissue RNA Extraction
Total RNA was extracted from frozen mouse lung tissue samples (−80° C) in RNAlater using a RNeasy Fibrous Tissue Mini Kit according to manufacturer's protocol (Qiagen, USA) as previously described (Pacurari et al., 2011). Total RNA was eluted in RNase-free water and stored at −80° C until further analysis. The quality and concentration of each RNA sample were determined using a NanoDrop-1000 Spectrophotometer (NanoDrop Tech, Germany).

Microarray Expression Profiling
Extracted RNA was analyzed for expression profiling using Agilent Mouse Whole Genome Arrays (Agilent, Santa Clara, CA). A universal reference design was employed using Stratagene Universal Mouse Reference RNA -Cat. No. 740100 (Agilent) as the reference RNA. Total RNA quality was determined on an Agilent 2100 Bioanalyzer with all samples having RNA integrity numbers (RIN) greater than 8. Total RNA (250ng) was used for labeling, using the QuickAmp labeling kit (Agilent). RNA extracted from each mouse was labeled with cyanine (Cy)-3-CTP (PerkinElmer, Waltham, MA) and reference RNA with (Cy)-5-CTP. Following purification of labeled cRNAs, 825ng of Cy3-and Cy5-labeled cRNAs were combined and hybridized for 17 h at 65° C in an Agilent hybridization oven. Microarrays were then washed and scanned, using an Agilent DNA Microarray Scanner.

Pathological Datasets
Inflammatory datasets were obtained by analysis of BAL fluid taken from MWCNTexposed mice at 1, 7, 28 and 56 days post-exposure as previously described . Mice were euthanized with an i.p. injection of sodium pentobarbital (>100 mg/kg body weight) followed by exsanguination. A tracheal cannula was inserted and BAL was performed through the cannula using ice cold Ca2 + andMg2 + -free phosphate buffered saline, pH 7.4, supplemented with 5.5mM d-glucose (PBS). The first lavage (0.6 ml) was kept separate from the rest of the lavage fluid. Subsequent lavages, each with 1ml of PBS, were performed until a total of 4ml of lavage fluid was collected. BAL cells were isolated by centrifugation (650×g, 5 min, 4 °C). An aliquot of the acellular supernatant from the first BAL (BAL fluid) was decanted and transferred to tubes for analysis of lactate dehydrogenase (LDH) and albumin. The acellular supernatants from the remaining lavage samples were decanted and discarded. BAL cells isolated from the first and subsequent lavages for the same mouse were pooled after resuspension in PBS, centrifuged a second time (650×g, 5 min, 4 °C), and the supernatant decanted and discarded. The BAL cell pellet was then resuspended in PBS and placed on ice. Total BAL cell counts were obtained using a Coulter Multisizer 3 (Coulter Electronics, Hialeah, FL) and cytospin preparations of the BAL cells were made using a cytocentrifuge (Shandon Elliot Cytocentrifuge, London). The cytospin preparations were stained with modified Wright-Giemsa stain and cell differentials were determined by light microscopy. Fibrosis datasets were obtained by morphometric analysis of Sirius Red staining for connective tissue in MWCNT-exposed mice at 1, 7, 28 and 56 days post-exposure as previously described (Mercer et al., 2011). Briefly, mice were euthanized by an overdose of pentobarbital (> 100 mg/kg body weight, i.p.) followed by transection of the abdominal aorta to provide exsanguination. To accomplish lung fixation, the trachea was cannulated and the lungs removed from the chest cavity. The lungs were then inflated with 1 ml of 10% neutral buffered formalin over a 1 minute period and the trachea tied off. After 4 to 5 hours, the lungs were trimmed and processed overnight in a tissue processor. For each animal, the left lung lobe was placed in the embedding carrier with a consistent apex to base orientation and embedded in paraffin. For morphometric studies, paraffin sections of the left lung (5 µm thick) were cut. A new region of the disposable knife blade was used to section each block and the water bath was changed frequently in order to prevent potential cross-contamination that might result from MWCNT passage on the knife between sections. The sections were then deparaffinized and rehydrated with a xylene-alcohol series to distilled water. To enhance the contrast between tissue and MWCNT, lung sections were stained with Sirius Red (Junqueira et al., 1979). Sirius Red staining consisted of immersion of the slides in 0.1% Picrosirius solution (100 mg of Sirius Red F3BA in 100 ml of saturated aqueous picric acid, pH 2) for 1 -2 hours followed by washing for 1 minute in 0.01 N HCl. Sections were then briefly counterstained in freshly filtered Mayer's hematoxylin for 2 minutes, dehydrated, and mounted on a slide with a coverslip. Quantitative morphometric methods were used to measure the average thickness of Sirius Red positive connective tissue fibers in the alveolar regions. Volume and surface density were measured using standard morphometric analyses (Weibel, 1980a;Weibel, 1980b). This consisted of basic point and intercept counting. Volume density was determined from counting the number of points over all tissues in the alveolar regions and points over Sirius Red positive connective tissue. Surface density of the alveolar wall was determined from intercepts between a line overlay and the alveolar wall. These point and intercept counts were made using a 121-point/11-line overlay graticule (12.5 mm square with 100 divisions), at 100× magnification, taken at six locations equally spaced across each section (one section per animal). This process was repeated twice for each animal. In order to limit the measurements to alveolar parenchyma, areas containing airways or blood vessels greater than 25 mm in diameter were excluded from the analysis. Average thickness of the Sirius Red positive connective tissue fibers of the alveolar wall was computed from two times the ratio of volume density of point to the surface density of the alveolar wall. The collagen fiber content of granulomatous lesions in the airspaces was assessed by a separate tabulation of points over Sirius Red positive connective tissues in granulomas and expressed as a percentage of total alveolar collagen. Mean linear intercept, a measure of the average size of the alveolar/alveolar duct airspaces in the alveolar region, was computed from the ratio of volume density to surface density (Weibel, 1980a).

Microarray Data Preprocessing and Filtering
Data were exported from the Agilent DNA Microarray Scanner using Feature Extraction v10 as tab-delimited text files after background subtraction, log transformation, and lowess normalization and reported as log or relative expression of the sample compared to the universal reference. Data were read from each file into R using a custom script. For each array, values for control spots, spots which were saturated on either channel, and spots which were not well above background on at least one channel were considered unreliable and/or uninformative and were replaced by "NA". Values were collated into a single table, and probes for which fewer than 10 present values were available were removed. For probes spotted multiple times on the array, values were averaged across replicate probes. The resulting table is available as a series matrix in the NCBI Gene Expression Omnibus repository with accession number GSE29042 . A web-interface (http:// www.mwcnttranscriptome.org) was developed to visualize the expression pattern of every gene in the whole genome in each MWCNT treatment condition.

Computational System
The computational system ( Figure 1B) was divided into four main components: a preprocessing component, the Pattern Finding component, the Coefficient Expander (CE) component, and the Functional Process Evaluation (FPE) component. Source code for the computational system can be found at http://sourceforge.net/projects/megpath. First, the preprocessing step was used to identify probes with significant changes in expression. Missing data were imputed using the K-means nearest neighbor algorithm as implemented by the impute.knn function in the impute R package from Bioconductor (Seattle, WA). Using the Bioconductor package, a set of differentially expressed genes for each dose and time point were identified by performing a two-class unpaired Significance Analysis of Microarrays (SAM) between the treated samples and the dose zero samples from the corresponding time point. A threshold delta value was chosen to produce a false discovery rate of 1% using the find Delta function from the same package. The list of probes called as significant was subsequently filtered by restricting those probes which were at least 1.5 fold up-or down-regulated. Fold changes were computed from the data before imputation of missing values.
Additionally, a linear model was fit to the data, modeling the log expression of each gene as a function of time, dose, and the interaction of time with dose. The t-statistic associated with the dose and interaction parameters following the SAM algorithm was moderated and a threshold set to control for a false discovery rate of 0.1%, thus generating a list of genes whose expression values were significantly dependent on dose and a list of gene whose expression values were significantly dependent on dose in a time-dependent fashion. The combined list of probes was described by Guo et al. (Guo et al., 2012) and was used by the Pattern Finding component in the current study.
Second, the Pattern Finding component was based on a Non-negative Matrix Factorization algorithm. This algorithm attempted to find a set of non-orthogonal basis vectors (patterns), which could be linearly combined to reconstruct the original probe expression data. In addition to finding the patterns, the algorithm also found coefficients relating each probe to each pattern. These coefficients could be used to describe how closely a probe matches a pattern. Most importantly, the Pattern Finding algorithm allowed for a probe to be associated with multiple patterns. In this way, patterns could be thought of as functions, such as fibrosis or inflammation, hence probes may be involved in multiple functions. The Pattern Finding algorithm worked as a Monte Carlo Markov chain with each location in the coefficient and pattern matrices having an associated probability density function. Quantitative histopathology data were used as one pattern.
The third step was to apply the Coefficient Expander component. This step attempted, through the use of simulated annealing, to find optimal coefficients for each probe in the genome-wide microarray data after preprocessing from the patterns found in the Pattern Finding step, and, therefore, enabled the probe's expression to be reconstructed from the patterns with minimal error.
The final step was to calculate the Functional Process Evaluation (FPE) score for a given pathway of genes. The FPE score was based on the enrichment score used in Gene Set Enrichment Analysis (Subramanian et al., 2005). Each gene's coefficients were normalized to obtain the relative importance of each pattern on the gene. Genes which were not common to both the biological process and microarray data after preprocessing were ignored and not included in the computation. Biological processes with fewer than 15 genes were excluded from further analysis. If a gene had multiple probes, the probe which could be reconstructed with least error was chosen. A biological processes' p-value was found by comparing its FPE score to the score of 1,000 randomly generated gene sets with the same number of genes. After p-values had been calculated for all processes, they were adjusted for multiple hypotheses testing by using the Benjamini and Hochberg method. A process with a p value less than 0.05 was considered significant. The leading set of a process is defined as the subset of genes, which was used to compute the Pathway Evaluation score. Genes are not restricted to being in only one leading set allowing for genes to be influenced by multiple patterns and used in multiple functions. The average expression of the genes in the leading set will strongly resemble the original pattern. However, a gene in the leading set is not required to look exactly like the pattern letting allowing for both known biological information and expression patterns to be incorporated.
Leading sets were found from the gene sets of the C2 Canonical Pathways and C5 Gene Ontology databases in MSigDB (Subramanian et al., 2005). The C2 Canonical Pathways database consists of 880 curated sets of genes corresponding to metabolic and signaling pathways. The C5 database consists of 1454 gene sets derived from the Gene Ontology project (http://geneontology.org).

Ingenuity Pathway Analysis
Data were analyzed through the use of Ingenuity Pathway Analysis (IPA) (Ingenuity® Systems, www.ingenuity.com). A network/My Pathway is a graphical representation of the molecular relationships between molecules. Molecules are represented as nodes, and the biological relationship between two nodes is represented as an edge (line). All edges are supported by at least one reference from the literature, from a textbook, or from canonical information stored in the Ingenuity Knowledge Base. Human, mouse, and rat orthologs of a gene are stored as separate objects in the Ingenuity Knowledge Base but are represented as a single node in the network. Nodes are displayed using various shapes that represent the functional class of the gene product. A total of 773 significant inflammation genes identified in the computational system were subjected to an Inflammatory Response -Inflammation overlay to determine which genes in the significant inflammation leading set were directly involved in inflammation according to IPA (Table 1). A total of 890 significant fibrosis genes were subjected to an Organismal Injury and Abnormalities -Fibrosis overlay to determine which genes in the significant fibrosis leading set were directly involved in fibrosis according to IPA (Table 2). To determine the interactions between genes which have only been experimentally observed in the lung, the Build-Trim tool of IPA was used. Direct and indirect interactions were trimmed to a Confidence Level of Experimentally Observed, and Tissue & Cell Lines included both Organ Systems of Lung and LungCell Lines.

Cell Culture
Small airway epithelial cells (SAEC) were cultured in SABM media (Lonza) supplemented with a SingleQuot Kit (Lonza). Cells were maintained at 37°C with 5% CO 2 .

Enzyme-Linked Immunosorbent Assay (ELISA)
SAEC were plated at 60,000 cells per well in a 24-well dish and grown at 37° C for 48 hours. Cells were serum starved overnight followed by exposure to 1 µg/ml or 2.5 µg/ml MWCNT for 24 hours. Conditioned media were collected and assayed for vascular endothelial growth factor A (VEGFA) and C-C motif chemokine 2 (CCL2) protein expression levels using DuoSet ELISA Development Systems from R&D Systems (Minneapolis, MN) according to manufacturer's protocol. Statistical analysis was done using a two-sample t-test assuming unequal variances.

Cellular RNA isolation
RNA was isolated from SAEC using RNAprotect Cell Reagent and an RNeasy Mini Kit from Qiagen according to the manufacturer's protocol (Qiagen, Valencia, CA). RNA concentrations were determined using a NanoDrop 1000 Spectrophotometer (NanoDrop Technologies, Wilmington, DE) and RNA quality was assessed using an Agilent 2100 Bioanalyzer (Agilent Technologies, Santa Clara, CA).

Real-Time Polymerase Chain Reaction
Total RNA (1µg) was converted into complementary DNA (cDNA) using a High Capacity cDNA Reverse Transcription Kit from Applied Biosystems (Life Technologies, Carlsbad, CA). All quantitative real-time PCR (qRT-PCR) reactions were performed on a 7500 Real-Time PCR system from Applied Biosystems. Each treatment group consisted of three biological replicates. qRT-PCR analysis for each biological replicate was performed in triplicate, and the Ct values obtained were normalized to the 18S housekeeping gene. Validated gene expression assays from Applied Biosystems were employed to carry out the mRNA expression profiling. The following gene expression assays were used: VEGFA (Hs00900055_m1); CCL2 (Hs00234140_m1); and 18S (Hs99999901_s1). Thermal cycling conditions were as follows: 50° C for 2 minutes, 95° C for 10 minutes, followed by 40 cycles of 95° C for 15 seconds and 60° C for 10 minutes.

Overview of an in vivo MWCNT exposure study and in vitro validation
A schematic of the overall method of determining and validating relevant processes related to lung inflammation and the progression to fibrosis after MWCNT exposure is depicted in Figure 1A. A total of 480 mice were randomized into three groups 1) gene expression profiling from snap frozen lung tissues, 2) BAL collection for inflammation assessment, and 3) lung tissue fixation for pathological analysis of fibrosis as determined by morphometric analysis of Sirius Red staining of lung tissue for collagen at 1, 7, 28, and 56 days postexposure to MWCNT by pharyngeal aspiration Mercer et al., 2011). Each time point consisted of 8 mice exposed to 0 (dispersion media [DM] control), 10, 20, 40, or 80 µg of MWCNT dispersed in DM for each animal group (Porter et al., 2008;Porter et al., 2010). Genome-wide mRNA expression profiles were analyzed by microarray through mRNA samples purified from the collected tissue and run on an Agilent Mouse Whole Genome Array. BAL was evaluated for the presence of polymorphonuclear leukocytes to assess inflammation, and morphometric analysis of Sirius Red staining for collagen in the alveolar walls was used to evaluate fibrosis Mercer et al., 2011). To determine significant transcription-related biological processes and genes with expression corresponding to MWCNT-induced lung inflammation or fibrosis patterns, the computational system ( Figure 1B) evaluated 41,059 probes on the microarray and established biological processes by incorporating the pathological data as input patterns in the simulation, as described in the Materials and Methods. For each significant biological process (BH adjusted p < 0.05), the corresponding leading set consists of the genes from the process which are most strongly related to the input pattern. Not all genes in the leading set will exactly resemble the pattern, but the average expression of the leading set will. The identified significant genes (SAM analysis; p < 0.05; FDR < 1%; fold change > 1.5) in the leading sets were then entered into IPA to identify genes functionally associated with inflammation and fibrosis and to depict molecular interactions in the lung. Based on the comprehensive evaluation, vegfa and ccl2 were selected for in vitro validation.

Identification of biological processes with expression patterns resembling MWCNTinduced inflammation or fibrosis pathology
A computational system was used to identify genes and biological processes with transcriptional activities, which matched the observed pathological patterns of lung inflammation or fibrotic collagen in the alveolar wall in the MWCNT-exposed mice ( Figure  1B). The preprocessing step found 2,996 unique probes which were significantly (p < 0.05; FDR < 1%; fold change > 1.5) up-regulated ( Figure 2A) or down-regulated ( Figure 2B) using Significance Analysis of Microarrays (SAM) or a linear model showing significant (p < 0.05; FDR < 0.1%) dose-response or dose and time interactions ( Figure 2C). Using this set of 2,996 genes, quantitative BAL and pathological data of MWCNT-induced inflammation or quantitative morphometric analysis of fibrosis were used as input patterns to find gene coefficients for reconstruction of the gene expression. Specifically, results for 3 sets of data were found, 2 sets relating to fibrosis (morphometrically determined changes in collagen within the alveolar wall) (Mercer et al., 2011) and 1 relating to inflammation (BAL) . Pathology data for fibrosis at dose 80 µg across the 4 time points was fit as an input pattern. The computational system found 64 total significant (BH adjusted p < 0.05) leading sets, the subset of genes which was used to compute the Functional Process Evaluation (FPE) score, representing the level of correlation with the fibrosis morphometric data for each biological process in the databases. Morphometric data for fibrosis occurring on day 56 across 4 doses was fit in the computational system with 84 significant (BH adjusted p < 0.05) leading sets found. Lastly, inflammation BAL scores at dose 40 µg across 4 time points was used, and 110 leading sets were found to be significantly (BH adjusted p < 0.05) correlated with the inflammation pattern.
Example results for each of the pathology data are shown in Figure 3. The average of the mRNA expression of genes in the leading set closely resembled the pathology data, indicating that in general, the transcriptional activities of the leading set genes correlated with changes in the pathology. The leading sets REACTOME GPCR LIGAND BINDING ( Figure 3A) and REACTOME HEMOSTASIS ( Figure 3B) were found in the C2 Canonical Pathways database and consisted of 156 genes (Supplemental Table 1) and 147 genes (Supplemental Table 2), respectively. The leading set of Immune System Process ( Figure  3C) was found in the C5 database and consisted of 163 genes (Supplemental Table 3). CCL2 ( Figure 3D) was contained in the leading set of REACTOME GPCR LIGAND BINDING. Although the CCL2 expression does not exactly follow the pattern, the average of all gene expression in the leading set does. The same can be seen for VEGFA ( Figure 3E). Importantly, our computational system does not constrain genes to being in only one leading set, allowing for genes to be involved in multiple processes. For instance, CCL2 was found to be involved in both MWNCT-induced fibrosis ( Figure 3A) and inflammation ( Figure 3F).

Determination of genes functionally involved in inflammation and fibrosis
To determine which genes were significantly altered in response to MWCNT exposure, leading set genes which attained a fold change of 1.5-fold or greater were input into Ingenuity Pathway Analysis (IPA) to determine if they were functionally involved in inflammation or fibrosis according to currently accepted literature.
The inflammation and fibrosis biological processes consisted of 773 and 890 unique genes, respectively, identified to be significantly altered (p < 0.05; fold change >1.5) after MWCNT exposure (significant inflammation) with a false discovery rate (FDR) of 1% in SAM analysis. Of the 773 significant inflammation genes, 67 were determined to be directly involved in inflammation by IPA (Supplemental Table 4). Of the 890 significant fibrosis genes, 69 were determined to be directly involved in fibrosis by IPA (Supplemental Table  5).
Using IPA and these 67 inflammation genes and 69 fibrosis genes, we determined those genes which were significantly involved in IPA Function and Disease Annotations associated with MWCNT-induced fibrosis. A recent report by Mishra et al. (Mishra et al., 2012) determined that low, physiologically relevant doses of MWCNT equivalent to those in our mouse study could significantly elevate the levels of transforming growth factor β (TGF-β) and matrix metalloproteinase-9 (MMP-9) in lung epithelial cells, as well as increase mechanisms of collagen production and cellular activation. Therefore, we used IPA to determine which genes in our significant (SAM analysis; p<0.05; FDR<1%; fold change >1.5) inflammation and fibrosis gene sets were involved in these processes (Supplemental Tables 4 and 5) Many inflammation genes were involved in general cell activation by functional association with the IPA function and disease annotations, including Cell Movement, Proliferation of Cells, and Morphology of Cells (Supplemental Table 4). Genes found in the significant inflammation set were also involved in the function and disease annotations, including Injury of Lung (ccl2,cd14,il6,il1r1,olr1,ptgs1,ptgs2,selp,sphk1,and tnf),Degradation of Connective Tissue (fn1,il6,il1b,il1rn,osm,ptgs1,ptgs2,and tnf), as well as the signaling pathway VEGF Signaling (pgf) (Supplemental Table 4). No significant inflammatory genes were found in the TGF-β signaling pathway according to IPA. Many fibrosis genes were also involved in the general cell activation function and disease annotations, such as Cell Movement, Proliferation of cells, and Morphology of Cells (Supplemental Table 5). Several genes in the significant fibrosis set were involved in the function and disease annotations, including Injury of Lung (adra2a,c3,ccl2,hif1a,il5,il6,il1r1,mmp12,ptgs2,selp,tnf,and vegfa),Degradation of Connective Tissue (fcgr2b,fn1,il6,il1b,il1rn,mmp13,osm,ptgir,tnf,and tnfrsf1b), as well as the signaling pathway VEGF Signaling (actc1, hif1a, ptk2, and vegfa) (Supplemental Table 5). Interestingly, 3 genes in the significant fibrosis set, bmpr2, smad4, and smurf2, were involved in the IPA TGF-β Signaling pathway, again suggesting that TGF-β signaling may play an important role in the progression of fibrosis and that the computational system was efficient in determining those biological processes which were functionally related to MWCNT-induced inflammation and fibrosis. An additional analysis of the significant inflammation ( Figure 5A) and fibrosis ( Figure 5B) genes by IPA determined those genes that have been experimentally shown to have an interaction specifically in the lung (Tables 1 and 2).

VEGFA and CCL2 in vivo and in vitro RNA expression
The inflammation leading set genes ( Figure 5C) and fibrosis leading set genes ( Figure 5D) were ranked by their frequency of inclusion in the biological processes significantly correlated with the pathological data. Two genes, ccl2 and vegfa, were selected for in vitro validation. Ccl2 was the top ranked gene that was involved in the most biological processes correlated with the inflammation ( Figure 5C) and among the top 20 genes involved in the most biological processes correlated with the fibrosis ( Figure 5D). Consistently, in the IPA lung interaction networks ( Figure 5A and 5B), ccl2 is in a hub that interacts with both TNF and IL1β hubs in the inflammation and fibrosis networks. Vegfa was found to be functionally associated with the fibrosis leading set and is integral to angiogenesis, or the formation of new blood vessels (Ferrara and Davis-Smyth, 1997). Angiogenesis is necessary for the formation of fibrotic tissue, and VEGF has been suggested as a serum biomarker for ranking the severity of idiopathic pulmonary fibrosis (Thannickal et al., 2004;Strieter and Mehrad, 2009;Ando et al., 2010). In a separate study, angiogenesis was observed after MWCNT exposure in human endothelial cells and in a coculture of both human epithelial and endothelial following epithelial exposure (Snyder-Talkington, In Review). Both genes were functionally validated in IPA analysis as involved in inflammation and/or fibrosis. Based on these results, ccl2 and vegfa were analyzed for their in vitro mRNA and protein expression levels following MWCNT exposure to validate the in vivo results. The top ranked gene for fibrosis [chemokine (C motif) receptor 1 (XCR1)] ( Figure 5D) was not shown to be functional involved in fibrosis in IPA analysis. Therefore, it was not selected for in vitro validation.
In vivo mRNA levels of vegfa showed stable expression levels across all days and doses with a significant decrease in expression on day 56 at dose 40 µg ( Figure 5B, Figure 6) and closely resembled the time-course of the morphometric collagen score data and leading set average of the biological process Reactome Hemostasis ( Figure 3B). Ccl2 showed a consistent dose-dependent increase in mRNA expression on all days with significant increases at all doses on day 1, doses 20, 40, and 80 µg on day 7 and doses 40 and 80 µg on day 56 ( Figure 5A, Figure 6). Ccl2 in vivo mRNA expression data closely resembled the fibrosis day 56 dose-response morphometric analysis and leading set average of biological process Reactome GPCR Ligand Binding ( Figure 3A) and was similar to the inflammation BAL pattern and leading set average for Immune System Process ( Figure 3C).
To assess the ability of MWCNT to induce similar RNA expression changes in vitro, SAEC were exposed to MWCNT at either 1 µg/ml (approximately equivalent to the in vivo dose of 20-40 µg ) or 2.5 µg/ml (approximately equivalent to the in vivo dose of 80 µg ) for 24 hours, and their mRNA expression levels analyzed. MWCNT exposure at both 1 and 2.5 µg/ml exposure levels induced modest but significant increases in vegfa mRNA expression in vitro in a dose-dependent manner ( Figure 7A). MWCNT exposure at both 1 µg/ml and 2.5 µg/ml levels induced an increase in ccl2 mRNA expression with a significant increase at 1 µg/ml ( Figure 7B).

VEGFA and CCL2 in vitro protein expression
To determine if the change in in vitro mRNA expression levels after exposure to MWCNT resulted in an increase in protein expression, conditioned media from cells exposed to either 1 or 2.5 µg/ml MWCNT for 24 h was collected and analyzed by ELISA for VEGFA and CCL2 protein expression. VEGFA showed significant increases in protein expression levels over control after 24 h of MWCNT exposure ( Figure 7C). CCL2 also showed significant increases in protein expression levels after 24 h of exposure ( Figure 7D). This demonstrated that the increase in mRNA expression levels of VEGFA and CCL2 after MWCNT exposure in vitro resulted in a concordant increase in protein expression and indicated that a similar increase may occur after in vivo exposure.

Discussion
Integrated in vivo and in vitro studies and in silico analysis is a recent endeavor in toxicological sciences. Novel methods for the analysis of current in vivo data are needed to develop predictive in vitro models so as to determine the toxicity profile of multiple material variants, such as various types of CNT. Our computational system was sufficient to identify potentially activated functions and pathways, which match inflammatory BAL scores and morphometric alveolar interstitial fibrosis data. By identifying the leading gene sets of the significant functions and pathways, our system can extract genes which are strongly associated with BAL markers and morphometry and that have potential involvement in inflammation and collagen production. The employment of IPA allowed for global analysis of our leading sets throughout the body of accepted scientific literature so as to target our results to those genes known to be involved in inflammation and fibrosis. The comparable results between gene expression profiles of our targeted genes in vivo with those found after MWCNT exposure in vitro suggest that our computational system is sufficient to determine potential outcomes of MWCNT exposure. This analysis may therefore identify potential signaling pathways and mechanisms that may be studied in vitro to determine potential in vivo outcomes and prognostic indicators of MWCNT exposure.
The epithelial lining of the lung is the first physical barrier to inhaled particles, and inflammation is a necessary process for the response to and recovery from lung injury. Upon injury, inflammatory mediators are released to trigger an immune response so as to remove the invading pathogen and allow for wound healing, involving processes such as reepithelialization and angiogenesis (Thannickal et al., 2004;Wynn, 2008). To replace the injured alveolar epithelial cells, alveolar type II cells dedifferentiate and move into the wounded area where they differentiate into new alveolar type I cells. Fibroblast-like cells in the lung secrete extracellular matrix (ECM) along which the alveolar type II cells move (Shi et al., 2009). Upon resolution of the injury, this ECM is typically reabsorbed and normal lung function and architecture is maintained; however, when an inflammatory response becomes chronic due to persistent injury or uncontrolled signaling, the inflammatory process can become pathogenic. Fibrosis, excessive collagen production, in the lung occurs when the deposition of ECM is poorly controlled and there is a loss of normal lung function and architecture (Shi et al., 2009). Because many of the regulators of normal development and inflammation also govern the process of fibrosis, there are multiple hypotheses on what actions must occur for the response to switch from normal inflammation to pathogenic fibrosis (Strieter and Mehrad, 2009). The exposure of mice by pharyngeal aspiration to 0, 10, 20, 40, and 80 ug of MWCNT was suggested to represent an exposure that could be compared to human occupational exposures . Sampling of lungs at 1, 7, 28, and 56 days post-exposure allowed the determination of gene expression changes that occur in the initial inflammatory stage as well as in later fibrotic stages of disease. These studies suggested that a single exposure to MWCNT can not only induce inflammation, but that MWCNT are biopersistent and induce a delayed fibrotic response as determined by increased collagen in the alveolar walls. Using a novel computational system, the correlation of global mRNA expression profiles to the changes in BAL score and morphometric analysis was analyzed. This identified transcription-related biological processes with expression patterns resembling the pathological patterns of inflammation and fibrosis in MWCNT-exposed mice, allowing for the identification of critical toxicity pathways and potential mechanisms for intervention. The results showed that this systematic analysis could identify relevant genes and pathways in MWCNT-induced lung injury from in vivo studies, which were further validated in in vitro experiments.
Previous studies used in vivo or in vitro genome-wide mRNA expression data to infer toxicity in carbon nanotube-exposed rats (Ellinger-Ziegelbauer and Pauluhn, 2009;Alazzam et al., 2010;Peng et al., 2010). In addition, a combination of microarray data, benchmark dose methods, and Gene Ontology annotations were used to identify potentially adverse biological processes in toxicity (Thomas et al., 2007;Burgoon and Zacharewski, 2008). Our novel computational system allowed for the discovery of non-parametric patterns, which could be used to reconstruct microarray data, incorporated quantitative pathological data, and was capable of working on both time-series and dose-dependent data. Unlike traditional clustering techniques (Tamayo et al., 1999;Waring et al., 2001;Yeung and Ruzzo, 2001), our system enabled genes to be included in multiple coexpression groups and be involved with multiple patterns. As a result, 23 genes were found to be involved in both MWCNTinduced inflammation and fibrosis in this study (Tables 1 and 2). The leading sets could be thought of as the genes in a function or pathway, which were most strongly influenced by a particular pattern. Although the expression of an individual gene in the leading set may not exactly match a pattern, the average of all the gene expression in the leading set will, indicating that the process as a whole responds similarly to the pattern. In this study, only genes in leading sets related to inflammation or fibrosis and from gene sets which were found to be significantly represented were studied. It should be noted that similar to the Bayesian Decomposition method, our computational system uses data collected from multiple conditions, in this study either 4 time points or 4 dose conditions, for each pattern matching and gene expression reconstruction. This computational model would not be applicable to modeling experimental data collected with less than 3 treatment conditions. In addition, since our computational model is non-parametric, it does not make any inference of unobserved experimental conditions. The use of IPA to determine if genes significantly altered in the leading sets were involved in inflammation or fibrosis allowed for an in depth analysis based upon data derived from relationships between genes and disease states taken from the currently accepted literature knowledge base. These analyses were rooted in and verified by experimental results collated from numerous sources. A total of 67 significantly altered genes were determined by IPA to be directly involved in the inflammatory process while 69 significantly altered genes were determine by IPA to be directly involved in fibrosis. Of the significantly altered genes, two genes, ccl2 and vegfa, were chosen to determine their in vivo and in vitro expression levels due to their roles in the cell during the development of inflammation and fibrosis as well as their rankings during gene profiling.
The formation of new blood vessels is an early response to tissue injury and a continuous process in the formation of fibrosis. Angiogenesis is necessary to sustain the tissue with oxygen, and increased angiogenic potential has been seen in patients with fibrosis (Thannickal et al., 2004;Strieter, 2008). VEGFA is a predominant angiogenic factor that acts upon endothelial cells for the proliferation of new blood vessels (Ferrara and Davis-Smyth, 1997). Because angiogenesis is integral to the formation of excessive ECM, we chose to determine if MWCNT had the ability to increase the expression of vegfa for both inflammatory and fibrotic processes. Additionally, we chose to determine the expression levels of ccl2 as an indicator of the inflammatory process. CCL2 is a known stimulator of the immune response, initiating chemotaxis in a variety of cell types, such as monocytes, lymphocytes, and basophils, as well as inducing the production of collagen from fibroblast cells. Due to these processes, CCL2 is suggested to play a role in inflammatory diseases (Rose et al., 2003). Interestingly, CCL2 has also been suggested to play a role in angiogenesis and upregulates the expression of VEGFA while, in turn, VEGFA has also been suggested to increase the expression of CCL2 (Hong et al., 2005;Yadav et al., 2010).
The dose-dependent increase in ccl2 mRNA expression at all days and doses in vivo suggests its role in the initial inflammatory process. Although the in vivo mRNA levels of vegfa remained relatively constant across all days and doses, the in vivo protein levels are unknown and may enhance collagen production. In vitro levels of ccl2 and vegfa mRNA also increased with increasing dose, reflecting what is seen in the in vivo analysis. In vitro analysis of the protein levels of CCL2 and VEGFA suggest that even modest increases in mRNA levels were able to significantly upregulate protein expression, and a similar increase in protein expression may occur in vivo. The analogous changes to vegfa mRNA levels in vitro, with subsequent increases in protein levels, suggests that MWCNT may have a similar effect in vitro to that seen in vivo. This may allow for potentially significant cellular processes to be identified by computational means and for the analysis of the mechanisms and signaling cascades behind MWCNT-induced effects to be validated in an in vitro manner.

Conclusions
A novel computational model was presented which was sufficient to determine transcription-related biological processes strongly associated with BAL and morphometric markers of lung inflammation and fibrosis, respectively, following exposure to MWCNT in mice. The biological processes were analyzed through IPA to determine genes and signaling pathways functionally involved in lung inflammation and fibrosis. Concordance of expression in two representative functionally involved genes in the in vivo analysis was confirmed in vitro, and the novel computational model was validated as a useful method to identify potential toxicity pathways. The use of these toxicogenomics data and in vivo animal model-based gene expression profiling integrated with in vitro verification may allow for successful in vitro toxicity profiling of MWCNT as well as the identification of potential signaling pathways involved in the etiology of MWCNT-induced injury.

Supplementary Material
Refer to Web version on PubMed Central for supplementary material.  Step 1: Preprocessing to identify significantly changing genes and identify potential interesting genes.
Step 2: Find patterns and coefficients from the interesting genes to reconstruct the gene expression data.
Step 3: Find coefficients for the entire genome.
Step 4: Using the patterns, coefficients, and pathways/functions, identify significant pathways.  Three leading sets found to be significant in a search of the C5 and C2 Canonical Pathways databases using pathological data. Computations were based on the observed experimental data points only; lines have been added to emphasize the patterns used in the computational system. For each pathway, (A) REACTOME GPCR LIGAND BINDING, (B) REACTOME HEMOSTASIS, and (C) IMMUNE SYSTEM PROCESS, the average of all the genes in the leading set shows strong similarity to the pathology data. The genes for these leading sets are listed in Supplementary Tables 1-3. Expression fold change values are shown for CCL2,  (A) IPA analysis of the 67 significant inflammation genes to determine those interactions, which specifically occur in the lung. (B) IPA analysis of the 69 significant fibrosis genes to determine those interactions, which specifically occur in the lung. (C) Ranking of significant inflammation genes by their frequency of appearance in biological processes significantly correlated with histopathological data. (D) Ranking of significant fibrosis genes by their frequency of appearance in biological processes significantly correlated with histopathological data. In vitro gene and protein expression data of VEGFA and CCL2. (A) In vitro vegfa gene expression in SAEC after exposure to 1 µg/ml and 2.5 µg/ml MWCNT for 24 hr. (B). In vitro ccl2 gene expression in SAEC exposure to 1 µg/ml and 2.5 µg/ml MWCNT for 24 hr. (C) ELISA results of VEGFA protein expression after DM (82 ± 11 pg/ml), 1 µg/ml (239 ± 23 pg/ml), and 2.5 µg/ml (408 ± 18 pg/ml) MWCNT exposure for 24 hr. (D) ELISA results Table 1 Gene information for 34 inflammation genes significantly altered in MWCNT-treated mouse lungs above 1.5 fold change with an FDR of 1% in SAM analysis and with functional molecular interactions in inflammation in the lung in IPA analysis. These genes were also strongly correlated with lung BAL scores in MWCNTtreated mice.