Transcriptional profiling of host gene expression in chicken embryo lung cells infected with laryngotracheitis virus

Background Infection by infectious laryngotracheitis virus (ILTV; gallid herpesvirus 1) causes acute respiratory diseases in chickens often with high mortality. To better understand host-ILTV interactions at the host transcriptional level, a microarray analysis was performed using 4 × 44 K Agilent chicken custom oligo microarrays. Results Microarrays were hybridized using the two color hybridization method with total RNA extracted from ILTV infected chicken embryo lung cells at 0, 1, 3, 5, and 7 days post infection (dpi). Results showed that 789 genes were differentially expressed in response to ILTV infection that include genes involved in the immune system (cytokines, chemokines, MHC, and NF-κB), cell cycle regulation (cyclin B2, CDK1, and CKI3), matrix metalloproteinases (MMPs) and cellular metabolism. Differential expression for 20 out of 789 genes were confirmed by quantitative reverse transcription-PCR (qRT-PCR). A bioinformatics tool (Ingenuity Pathway Analysis) used to analyze biological functions and pathways on the group of 789 differentially expressed genes revealed that 21 possible gene networks with intermolecular connections among 275 functionally identified genes. These 275 genes were classified into a number of functional groups that included cancer, genetic disorder, cellular growth and proliferation, and cell death. Conclusion The results of this study provide comprehensive knowledge on global gene expression, and biological functionalities of differentially expressed genes in chicken embryo lung cells in response to ILTV infections.


Background
Infectious laryngotracheitis virus (ILTV; gallid herpesvirus 1) is the only member of the Iltovirus genus of the Alphaherpesvirinae subfamily of the Herpesviridae family. ILTV includes~150 kb of linear dsDNA genome consisting of two unique regions (unique long; UL and unique short; US), inverted repeats (IR) and terminal repeats (TR) flanking the US region [1]. About 76 open reading frames (ORFs) have been shown to express viral proteins in ILTV [2]. The genome structure and gene contents of the ILTV genome clearly prove its classification as an alphaherpesvirus [3]. Infection of ILTV causes an upper respiratory disease in chickens during lytic infection, and ILTV can establish latency in the central nervous system.
Respiratory symptoms of ILTV infection include dramatically increased mucus formation in the trachea and tracheal hemorrhage that can cause up to 70% mortality. Currently, live attenuated vaccines developed from chicken embryo or cultured cells are commercially available to control ILTV disease [4]. However, vaccinal laryngtracheitis (VLT), resulting from reversion of vaccine virus to virulent form and spreading from vaccinated-to unvaccinated birds, is often associated with the use of live attenuated ILTV vaccines [5,6].
Microarray data sets can be interpreted further by clustering analysis. Many of the heuristic clustering methods have several shortcomings; these include the determination of the number of clusters which generally is unknown when there is no prior knowledge of the number or there no other information about the structure of the data to be clustered. A model-based clustering method can overcome these critical drawbacks by treating a clustering problem as a model selection problem over a variety of candidate models specified by different numbers of clusters and distribution, and by estimating the number of clusters in the clustering analysis. The best model is selected on the basis of a model selection criterion, simultaneously providing the optimal number of clusters and assigning cluster membership to observations.
To our knowledge, the effects of ILTV infection on changes in global gene expression in host cells have not been previously reported. ILTV is a special type of herpesvirus that causes acute respiratory disease in poultry. Thus, the objective of this study was to understand host responses to ILTV infection in cultured chicken embryo lung cells using microarray analysis. The microarray used in the current study contains 44 K chicken genes including functionally identified genes, predicted ORFs, ESTs, genomic contigs, chicken microRNAs and various control spots [27]. Importantly, functional analysis of differentially expressed genes should follow gene discovery research. To this end, a software program (Ingenuity Pathway Analysis), specifically developed to analyze large data sets such as microarray data for biological functionalities, gene networks, and physiological pathways [28], was used to assign biological functionalities and molecular interactions in chicken embryo lung cells after 1 to 7 days in response to ILTV infection.

Results and Discussion
Gene expression profile of lung cells infected by ILTV Primary chicken embryo lung cells at passage 1 were infected by the USDA reference strain of ILTV and cells were collected at 1, 3, 5, and 7 dpi. Cytopathic effects (CPE) were observed by 3 dpi, which became more severe by 5 dpi. Massive cell disruption was observed at 7 dpi ( Figure 1). Total RNA was isolated from both controls and infected lung cells at each dpi time point and subjected to microarray analysis.
To control for dye bias effects, spike-in control mixtures were utilized by mixing with RNA samples according to the manufacturer's recommendations (See Methods). The spike-in RNA controls consisted of two sets of synthetic RNA mixtures derived from the Adenovirus E1A genes with different concentrations in each set [29]. The Agilent chicken 4 × 44 K oligo gene expression array contains 320 spike-in indicator spots to be hybridized with spike-in controls of both A mix, which was hybridized with Cy 3, and the B mix hybridized with Cy 5 on each array. These spike-in sets were mixed with either uninfected control or infected samples and co-hybridized to arrays. The ratio of signal intensities for all spike-in spots were calculated, evaluated, and revealed no significant dye effects on all array slides (data not shown) as reported previously [29]. All raw and normalized data were deposited in the Gene Expression Omnibus (GEO; accession number: GSE20630).
Normalized signal intensities were subjected to statistical analysis to find differentially expressed genes during ILTV infection in cultured embryonic lung cells. The 44 K array revealed 11,491 genes with significant signal intensities that were sorted by signal to noise ratio (SNR) >3, meaning that real (forward) signals of the samples were three times greater than background signals. In order to discover time course change in gene expression patterns, a model-based method [30] was used for clustering the gene expression profiles. A key drawback in heuristic clustering techniques is that it is difficult to determine the number of clusters a priori. The method enables the number of clusters to be determined by estimating the number of components in a multivariate normal mixture model from which the data are generated. The clustering analysis resulted in three gene groups ( Figure 2). Group 1 included a total of 789 genes that showed significant differential expression in response to ILTV, Group 2 included 6,265 genes that displayed moderate alterations, and Group 3 included 4,437 genes that revealed no alterations during ILTV infection at four time points in chicken lung cells. Of the 789 genes in Group 1 exhibiting differential expression in response to ILTV (see Additional file 1), the top 10% (79 genes) were sorted by statistical analysis based on the highest value of standard deviations using the mean values of four different time points (Table 1). This approach highlights genes with more significant alterations in response to ILTV overtime. Out of the 789 genes, 390, 370, 320, and 422 genes were down-regulated, while 399, 419, 469, and 367 genes were up-regulated relative to uninfected cells at 1, 3, 5, and 7 dpi, respectively.

Quantitative reverse transcription-PCR (qRT-PCR)
To validate the microarray data, 20 genes were subjected to qRT-PCR using gene specific primer pairs and the same RNA samples as in the microarray analysis (Table 2). Results were analyzed by 2 -ΔΔCt method to determine relative levels of gene expression at each dpi time point compared to uninfected control [31]. There were no differences between microarray data and the qRT-PCR at any dpi time point (Table 3). However, it should be noted that fold change values for certain genes obtained by qRT-PCR analysis showed much greater expression levels than those observed in the microarray analysis. For example, the fold changes for the gene expression of matrix metalloproteinase (MMP) 27, interleukin (IL) 6, fatty acid binding protein (FABP) 4, IL8, and CXC chemokine K60 at 3-or 5 dpi showed much higher levels in qRT-PCR analysis compared to fold changes shown in microarray analysis (Table 3). Possibly, this qualitative difference between methodologies may be attributed to the upper detection limits of the fluorescent intensities for the array scanner. Based on quality control measures, such as the spike-in controls and the results of targeted qRT-PCR indicate that the microarray data sets for differential gene expression are valid to investigate genome-wide differential expression patterns for host responses during ILTV infection.

Expression clustering
The pattern of differential gene expression over time can provide insights into biologically functional relevance among genes. In the present study, a model-based clustering method [30] was used to cluster alteration patterns for the 789 differentially expressed genes in response to ILTV infection and revealed 7 gene clusters exhibiting distinct expression patterns ( Figure 3 and Additional file 2). The 287 genes placed in cluster (C) 1 showed only nominal increases at 3 and 5 dpi followed by decreased expression levels at 7 dpi that were similar to those at the onset of the experiment. The C2 representing 97 genes exhibited a dramatic increase in gene expression only at 7 dpi, whereas the expression levels of the 90 genes in C3 progressively declined at 5 and 7 dpi. Three genes in C4 showed higher expression during early infection (1 dpi), sharp increases at 3 and 5 dpi, followed by a slight decline at 7 dpi. Expression patterns of 9 genes in C5 showed slightly lower expression at 1 dpi relative to the other time points, a dramatically increase at 3 and 5 dpi, followed by decreased expression at 7 dpi. The 85 genes in C6 showed lower expression at 1 dpi followed by a progressive increase during the later time points, which was opposite to 218 genes in C7 that showed higher expression at 1 dpi followed by decreased expression at 3, 5, and 7 dpi. GenBank accession numbers for genes in each cluster are shown in the Additional file 2.
Interestingly, the genes in C4 that exhibited the highest expression during ILTV infection include cytokines (IL8 and IL1β) and a chemokine (CXC-K60), while in the C5, IL6 was most highly expressed. From these findings, it is reasonable to hypothesize that expression of functionally relevant genes such as cytokines and chemokines released in response to an immune challenge may be regulated similarly during a specific immunological challenge.

Functional Gene Ontology (GO)
Recently, new bioinformatics tools have been developed to facilitate efficient analysis of biological functionalities for large numbers of differentially expressed genes obtained from microarray analysis. By using the IPA program (http://www.ingenuity.com/), bioinformatics aspects of differentially expressed genes during ILTV infection were analyzed for the relevance of gene functionalities and gene networks. While 789 differentially expressed genes were used as the input number of genes, only 275 have been characterized with specific cellular functions according to the IPA program. Results obtained with the IPA program in terms of numbers of genes, biological functions of genes and categories were independently confirmed with a second bioinformatics tool, DAVID (Database for Annotation, Visualization and Integrated Discovery) version 6.7 (http://david.abcc. ncifcrf.gov/) (data not shown). The group of 275 differentially expressed genes was placed into 65 functional groups (see Additional file 3) and the top 25 functional groups of genes are displayed in Figure 4. The main categories for gene functionalities include diseases and disorders, molecular and cellular functions, and physiological system development. Gene information was repeatedly used in multiple groups of functionalities due to the multi-functional characteristics for designated genes. It can be seen that the cancer related function contains the highest number (140) of genes, while 125 genes were involved in genetic disorders, and 54 genes were grouped as inflammatory responses.

Gene network analysis
Gene network analysis, which represents the intermolecular connections among interacting genes based on functional knowledge inputs, was performed on the differentially expressed genes using IPA program (see Additional file 4). In this way, 21 possible gene networks for all days post ILTV infection were generated based on differential gene expression. Of these 21 gene networks, only 6 gene networks were identical at all time points following ILTV infection. Possibly, the reason why only 6 of 21 networks were identical might be because of differences in the sets of focus molecules which are generated from p-values and fold change values of differential gene expression that are used in IPA algorithms. The dynamics of alterations in gene expression for a subset of genes during the time course of ILTV infection can provide insights into biologically interacting genes within a network that display functional similarities. The most interactive network (network #1) is presented in Figure 5 whereas the other networks are shown in Additional file 5. Network #1 is closely associated with a signaling pathway of IL6, which is a cytokine known to be involved in cell proliferation and inflammatory responses [32]. The top functions related with genes in network #1 involve cancer, gastrointestinal disease, and the cell cycle. Interestingly, expression of certain genes in network #1 such as CDC20 (cell division cycle 20 homolog), PTTG1 (pituitary tumor transforming 1), CDC2, and Cyclin B, which are associated with cellular proliferation in cell cycle progression, appeared to be inversely related to IL6 expression. The dynamics of alterations in gene expression over time during ILTV infection suggest that ILTV infection elevates IL6 expression followed by the inhibition of cellular proliferation. In contrast, expression patterns of HPGD (hydroxyprostaglandin dehydrogenase 15-NAD), SOCS (suppressors of cytokine signaling), JAK (Janus kinase 1), and NASP (nuclear autoantigenic sperm protein) were independent of the IL6 expression pattern ( Figure 5A, B and Additional file 5A). JAK is known to enhance cellular proliferation through the signal transducer and activator of transcription (STAT) pathway that can be suppressed by IL6 signaling [33]. The consistent downregulation of JAK supports a role of JAK in the repression of cellular proliferation by ILTV infection. The top functions of genes in network #2 are involved with cellular compromise, connective tissue disorders, and post-translational modifications. Several heat shock proteins (HSP) were also focused in this network (see Additional file 5B). Heat shock proteins, especially the HSP70 family that serve as molecular chaperones, are known to interact with viral early immediate genes in HSV-1 genomic DNA replication [34]. Interestingly, since the expression of several HSPs in network #2 were downregulated, it is reasonable to hypothesize that the lower HSP through 5 dpi may lead to production of erroneous virion structures of ILTV that in turn results in low ILTV titers in tissue culture, which has been reported to barely exceed one infectious unit per cell [1,35]. Network #3 contains genes for growth factors and matrix metalloproteinases (MMPs), and genes have top functions associated with endocrine system function and development, carbohydrate metabolism, and digestive system function and development (see Additional file 5C). Expression for growth factors and MMPs increased on 3 dpi and remain elevated through 7 dpi. This observation is in agreement with other reports that infection with herpesviruses, such as HSV and HCMV, lead to an increase in growth factor expression and MMPs for extracellular remodeling, tissue invasion and angiogenesis [36][37][38]. Networks #4 and #5 contain cytokine genes (IFNβ and IL1β), chemokine genes (CCL20 and CCL4), and genes of the NF-kB families (NF-kB and NFIB) with top functions that are involved in organism injury and abnormalities, antigen presentation, cell mediated immune responses, lipid metabolism, small molecule biochemistry, and molecular transport (see Additional file 5D and 5E). Finally, genes in network #6 contain IL1, NF-kB, and ID1 that function in cardiac inflammation, cardiovascular disease and in the inflammatory response (see Additional file 5F). Interacting molecules found in networks #4, 5, and 6 are mostly focused on the host immune responses against pathogenic infections.
The network analysis suggests that a large number of biological pathways, regulated by various sets of genes, closely interact with each other in host responsiveness during ILTV infection. More detailed interactions among genes showing altered expression levels in each network are currently under investigation to identify host-response mechanisms that may occur in conjunction with general immunological reactions during ILTV infections.
The fold changes in gene expression of key molecules associated with cellular immune response, cell signaling, MMP molecules, cytokines, chemokines, and cell proliferation were plotted individually ( Figure 6). This was The gene expression levels of 20 genes at four different time points in microarray analysis were confirmed by qRT-PCR. The expression levels were presented by fold changes values in microarray analysis, while, for qRT-PCR, the values were calculated by 2 -ΔΔCT method [31], which would be comparable to fold changes in the microarray. All values are mean values determined by the calculation from three replicating assays. pathogenesis [39]. Similarly the viral oncoprotein meq in MDV is known to activate MMP3 transcription [40]. Furthermore, the balance between MMP9 and tissue inhibitor of metalloproteinases 1 (TIMP1) was altered in human macrophages in HCMV infection, such that MMP9 activity declined in response to HCMV infection. Moreover, it was shown that HCMV infection may affect atherogenesis in mice through the control of MMP9 expression [38]. Taken together, these results suggest that MMPs generally play a role in herpesvirus pathogenesis, but different isoforms of MMPs may be capable of responding to specific herpesviruses. Furthermore, in the present study with chicken lung cells, expression of surfactant protein A (SFTPA-1; GenBank accession -AF411083; Table 1 and 3) was downregulated by ILTV infection. Since the SFTPA-1 is a transcriptional indicator of EGFR (epidermal growth hormone receptor) signaling pathway [41], the reduction in SFTPA-1 expression suggests that the EGFR signaling pathway is suppressed by ILTV infection. This conclusion is consistent with the reported downregulation of EGFR functions in HCMV infected human lung [41] and foreskin fibroblastic cells [42]. Additionally, our findings are similar to a report in which the mRNA expression of SFTPA-1 declined after inoculation of influenza A virus H9N2 into chicken lung cells [43]. Another molecule, ovotransferrin (Otrf; GenBank accession-X02009; Table 1 and 3) has both iron transport-and antibacterial-activities. In the present study, the expression level of Otrf increased in ILTV infected cells ( Figure 6A). The finding is qualitatively similar to a possible role of Otrf in MDV infection that might protect the spread of MDV in chicken embryonic fibroblast cells [44]. Otrf also was found to accelerate the expression of immune response gene such as IFN-γ against MDV infection [45]. Furthermore, lactoferrin, which is a homologous form of transferrin in mammals, showed antiviral activity against canine herpesvirus [46]. It has been hypothesized that the anti-viral activity of Otrf may be similar to the anti-HSV capability of mammalian transferrins [47]. IL6 plays a role in both pro-inflammatory and anti-inflammatory responses [32]. The elevation of IL6 expression during virus infection has been wellcharacterized as one of the immune response activities in the pathogenesis of various viruses, such as Dengue virus in human [48], or HSV-1 in mice [49]. Also, the expression of IL6 is increased by KHSV-encoding miRNA [50]. Similarly, expression of IL8 and IL1β were upregulated in ILTV infected cells and also in MDV infected chicken lung [51]. These observations are consistent with previous findings of NF-kB activation via IL8 signaling pathway by human herpesvirus (HHV)-8 infection [52] and HHV-6 infection in peripheral blood mononuclear cell cultures [53], respectively.
In addition to the well-characterized host-virus interactions, a variety of unique responses detected with the microarray analysis of ILTV infection in the present study. For instance, the expression level of vasoactive intestinal peptide (VIP; GenBank accession number-U09350: Table 1 and 3) decreased dramatically in ILTV infection, but the functional consequences were not determined. Likewise, genes related to various metabolic enzyme functions such as acyl-CoA synthetase long chain-1 (ACSL1; GenBank Accession number-AJ851480; Table 1 and 3) were differentially expressed in ILTV infected cells, but the precise mechanisms in the host response have not been verified. Therefore, further investigations are being performed to identify unique and more deeply involved interactions between host lung cells and ILTV.

Conclusions
In this study, we have demonstrated changes in global gene expression in response to ILTV infection in chicken embryo lung cells using microarray analysis. A total of 789 differentially expressed genes were involved in a variety of molecular and cellular defense mechanisms of immune response, cell cycle regulation, cellular metabolism, and matrix metalloprotenases. Moreover, the bioinformatics studies, such as gene ontology and gene network analysis, using knowledge based bioinformatics tools (Ingenuity Pathway Analysis and DAVID) revealed biological functionalities and intermolecular connections among interacting genes associated with differentially expressed genes. Our study provides epigenetic insights into the pathogenesis of ILTV infection in chicken lung cells.

Cell culture and ILTV infection
Cell culture reagents were purchased from Invitrogen Life Technologies (Carlsbad, CA, USA). Chicken embryo lungs were isolated from 19 day old specific-pathogen free (SPF) chicken embryos (Charles River Laboratories, North Franklin, CT, USA). Lung tissues were homogenized and incubated in a 0.125% trypsin solution for 30 min at room temperature (25°C). Cells dissociated from lung tissues were suspended in a 1:1 ratio of mammary epithelial growth media (MEGM; Lonza, Rockland, ME, USA) and Dulbecco's Modified Eagle's Medium (DMEM, 0.45% glucose) plus 2% fetal bovine serum (FBS), 100 units/ml penicillin, 100 μg/ml streptomycin, and 2 mM L-glutamine in 10 cm tissue culture dishes (Sarstedt Inc., Newton, NC, USA) pretreated with 0.5% gelatin in PBS to improve cell adhesion. Cultured cells were grown at 39°C containing 5% CO 2 until cells reached confluent monolayers (2 to 4 days). The USDA reference strain of ILTV (National Veterinary Services Laboratories, Ames, IA, USA) was used to infect the chicken embryonic lung cells at a multiplicity of infection (m.o.i.) of 0.1. Infected cells were incubated at 37°C for 1 hr with rocking gently every 15 min. After the incubation, 10 ml of media, 1:1 MEGM/DMEM, were added to each culture dish, and the cells were incubated at 37°C in 5% CO 2 for up to 7 days. This research was performed under the permitted protocol approved by both the Institutional Biosafety Committee (IBC; permit number: 10007) of University of Arkansas and the Animal and Plant Health Inspection Service (APHIS; permit number: 102743) of United States Department of Agriculture (USDA).

Total RNA extraction
Total RNA was extracted from uninfected-or ILTV infected chicken embryonic lung cells at 1, 3, 5, and 7 dpi using TRIzol reagent (Invitrogen Life Technologies, Carlsbad, CA, USA) following the manufacturer's instructions. Total RNA was treated with DNase I (New England BioLabs Inc., Ipswich, MA, USA), and RNA was re-purified by the TRIzol reagent. The quality of RNA was checked by fractionation on an agarose gel (data not shown).

Probe labeling and microarray hybridization
A two color labeling microarray system was used to compare uninfected-and ILTV infected embryonic lung cells at 1, 3, 5, and 7 dpi. Fluorescently labeled complementary RNA (cRNA) probes were generated by using the Two Color Microarray Quick Labeling kit (Agilent Technologies, Palo Alto, CA, USA) and following the manufacturer's instructions. RNA spike-in controls were used to adjust possible dye effects following manufacturer's instructions. The Spike-in controls represent two sets of ten synthesized RNA mixtures derived from the Adenovirus E1A transcriptome with different concentrations in each set [29]. These spike-in sets were mixed with either uninfected control or infected samples and co-hybridized to arrays. Briefly, 2 μg of total RNA were mixed with Spike-ins and converted to cDNA using reverse transcriptase and oligo dT primers in which T7 promoter sequences were added. T7 RNA polymerase was used for the synthesis and labeling of cRNA with either Cy3 dye for the uninfected control or Cy5 dye for the ILTV infected samples. The fluorescently labeled cRNA probes were purified using the Qiagen RNeasy Mini Kit (Qiagen Inc., Valencia, CA, USA), and the concentration, fluorescent intensities, and quality of labeled cRNA probes were determined using a Nano-drop spectrophotometer (Thermo Scientific, Wilmington, DE, USA). An equal amount (825 ng) of Cy3 and Cy5 labeled cRNA probes were hybridized on a 4 × 44 K Agilent custom chicken oligo microarray (array ID: 017698). The hybridized slides were washed using a commercial kit package (Agilent Technologies, Palo Alto, CA, USA) and then scanned using a Genepix 4000B scanner (Molecular Devices Corporation, Sunnyvale, CA, USA) with the tolerance of saturation setting of 0.005%. Three biological replicates were conducted.

Microarray data collection and analysis
Background-corrected red and green intensities for each spot were used in the subsequent analysis. Global normalization based on local polynomial regression (loess) was applied to the intensities to remove effects that were due to undesirable systematic variations in microarray experiments rather than biological differences. The average values of the resulting normalized expression values in replicate hybridization sets were considered in the subsequent analysis. In order to identify a given set of genes that exhibited major alterations over time, a model-based clustering method [30] was employed, and the genes in the cluster were considered as differentially expressed over the time period. All analytic techniques were implemented in R (http://www.R-project.org).

Quantitative reverse transcription-polymerase chain reaction (qRT-PCR)
Reverse transcription was performed with 3 μg of total RNA using Superscript II reverse transcriptase (Invitrogen Life Technologies, Carlsbad, CA, USA) with oligo dT [12][13][14][15][16][17][18] primers (Invitrogen Life Technologies, Carlsbad, CA, USA) following the manufacturer's instructions. The reverse-transcribed cDNA were diluted by 1:10 ratio and a portion (1 μl) of each product was subjected to qPCR under the following conditions: 40 cycles of 95°C for 30 s, gene-specific annealing temperature for 62°C for 1 min, extension for 30 s at 72°C, and a final extension at 72°C for 10 min. A non-template control and endogenous control (chicken GAPDH) were used for the relative quantification. The differential expression levels for the ILTV infected group were compared by the 2 -ΔΔCT method against the uninfected controls [31]. Primers for qRT-PCR were designed using Primer3 software (http://frodo.wi.mit.edu/cgi-bin/primer3/pri-mer3.cgi) with these parameters: amplicon length, 95-100 bp; primer length, 18-27 nucleotides; primer melting temperature, 60-64°C; primer and amplicon GC content, 20-80%; difference in melting temperature between forward and reverse primers, 1-2°C. Primers were synthesized by Integrated DNA Technologies (Coralville, IA, USA). Primer information is listed in Table 2.

Bioinformatics
Functional interpretation of differentially expressed genes was analyzed in the context of gene ontology and molecular networks using the Ingenuity Pathways Analysis (IPA) 6.5 software (Ingenuity Systems®; http://www.ingenuity. com). The differentially expressed genes were categorized, compared to genetic categories in the IPA database, and ranked according to p-values [54]. The IPA analysis determined the subcategories within each category which is supplied with an appropriate p-value and the number of genes identified. Since the size of the created network could potentially be enormous, the IPA software limited the number of molecules in the network to 35, leaving only the most important ones based on the number of connections for each focus gene (focus genes = a subset of uploaded significant genes having direct interactions with other genes in the database) to other significant genes [55].