Gene Expression in the Spinal Cord in Female Lewis Rats with Experimental Autoimmune Encephalomyelitis Induced with Myelin Basic Protein

Background Experimental autoimmune encephalomyelitis (EAE), the best available model of multiple sclerosis, can be induced in different animal strains using immunization with central nervous system antigens. EAE is associated with inflammation and demyelination of the nervous system. Micro-array can be used to investigate gene expression and biological pathways that are altered during disease. There are few studies of the changes in gene expression in EAE, and these have mostly been done in a chronic mouse EAE model. EAE induced in the Lewis with myelin basic protein (MBP-EAE) is well characterised, making it an ideal candidate for the analysis of gene expression in this disease model. Methodology/Principal Findings MBP-EAE was induced in female Lewis rats by inoculation with MBP and adjuvants. Total RNA was extracted from the spinal cords and used for micro-array analysis using AffimetrixGeneChip Rat Exon 1.0 ST Arrays. Gene expression in the spinal cords was compared between healthy female rats and female rats with MBP-EAE. Gene expression in the spinal cord of rats with MBP-EAE differed from that in the spinal cord of normal rats, and there was regulation of pathways involved with immune function and nervous system function. For selected genes the change in expression was confirmed with real-time PCR. Conclusions/Significance EAE leads to modulation of gene expression in the spinal cord. We have identified the genes that are most significantly regulated in MBP-EAE in the Lewis rat and produced a profile of gene expression in the spinal cord at the peak of disease.


Introduction
On the basis of pathological findings and genetic studies that show that the majority of genes that predispose to multiple sclerosis (MS) are associated with T cell activation [1], MS is thought to be primarily a T-cell mediated autoimmune disease. B cells and antibodies are also thought to play a role in pathogenesis [2]. The symptoms of MS are due to demyelination and axonal loss, and possibly to temporary loss of nerve conduction due to circulating factors [3]. Experimental autoimmune encephalomyelitis (EAE) is the best available animal model of MS, and can be induced by inoculation of susceptible animals with a range of central nervous system (CNS) antigens. EAE was first induced by inoculation of whole brain tissue [4;5]. Myelin basic protein (MBP) was the first encephalitogenic protein to be studied in detail. Later studies have used myelin proteolipid protein (PLP) [6;7] and myelin oligodendrocyte glycoprotein (MOG) [8] to induce EAE. The susceptibility to EAE varies among animal strains and with the antigen used to induce the disease. There have been studies of the genetic loci that contribute to this susceptibility. The most important genes are the MHC class II, but other genes have been implicated, including other immune related genes [9][10][11]. Studies have also identified genetic loci that contribute to different EAE traits such as the severity of disease, maximum clinical score and the presence of demyelination [9;10].
In this study we have investigated MBP-EAE, which is characterized by weight loss and ascending paralysis [12] followed by spontaneous recovery. It has been clearly shown that recovery from weakness is associated with restoration of electrical conduction [13]. The pathology of MBP-EAE involves infiltration of the spinal cord and nerve roots with inflammatory cells [14], followed by demyelination of axons. In MBP-EAE there is apoptosis within the spinal cord [15;16] and there is also axonal damage [17]. The metabolic processes that cause damage in EAE include glutamate toxicity [18] and dysregulation of ion channels including calcium channels [19]. Because of the detailed knowledge of the pathogenesis of this EAE model, it is of interest to determine how gene expression correlates with the known features of this disease. The advantages of using EAE for studies of gene expression in the CNS are that tissue can be easily obtained at the peak of disease, without a significant post-mortem delay, allowing good preservation of RNA for expression studies.
Because of the extensive characterization of MPB-EAE in the Lewis rat, it is an excellent candidate for studies aimed at defining the gene expression profile of the disease. There have been no previous gene array studies of MBP-EAE in the Lewis rat; however, the expression profile of selected cytokines, chemokines, adhesion molecules, cell markers, matrix metalloproteinases and other genes relevant to disease etiology have been examined in this model by real-time PCR (RT-PCR) [20]. There have been microarray studies of gene expression in MBP-EAE in SJL/J mice [21] and of the effects of estrogen treatment on gene expression in MBP-EAE in BV8S2/AV4 double transgenic mice [22]. There have also been microarray studies of EAE induced with myelinoligodendrocyte glycoprotein (MOG) in the spinal cords of mice [23] and rats [24]. The present study was carried out to investigate changes in gene expression in MBP-EAE compared to control Lewis rats, and found that there were significant differences from controls.

a) Differential Gene Expression in MBP Induced EAE in the Lewis Rat
Rats were inoculated with MBP to induce EAE. The typical clinical course of this disease is shown in Figure 1A, which demonstrates that disease is generally monophasic with spontaneous recovery. Rats were sacrificed on day 13 post injection which corresponds approximately with the peak of clinical disease in our model. The course of disease of the rats that were used in the micro-array study is shown in Figure 1B.
RNA quality analysis was carried out on the BioRad Experion automated electrophoresis system. All preparations used in both assays had values of values of .9.5 for RNA quality indicator (RQI) [25]. This is illustrated in Figures S1 and S2. After array hybridization with the Affimetrix Rat Exon 1.0 ST array, which contains probe sets for 27,324 genes, data was uploaded for analysis with the Partek genomics suite. Transcripts from 8,793 genes were detected. Of these, we identified 2,350 genes that were significantly regulated between the healthy rats and rats with MBP-EAE (t test, (p-value ,0.05). The overall gene expression changes in MBP induced EAE are illustrated in the volcano plot in Figure 2, which shows a large number of genes that were highly regulated compared to healthy controls.
Of the 2,350 differentially regulated genes, 2,265 mapped to known locations in the Partek data-base and the remainder were unmapped. Of these mapped genes, 1,665 were mapped to named genes. Of the 1,665 differentially regulated genes, 998 were found to be down-regulated and 667 were up-regulated. These are listed in Tables S1 and S2.
The 2,530 genes that were differentially regulated were uploaded into the Ingenuity Pathway analysis software package (IPA) (Ingenuity H Systems, www.ingenuity.com ) for functional and pathway analysis. After filtering to exclude genes that were not of relevance to the nervous system, we produced an annotated set of 1,190 differentially regulated genes. These genes were used for our further analysis of functions, pathways and upstream regulators.
To produce a list of the genes that are most significantly and highly regulated in MBP-EAE, we used a False Discovery Rate (FDR) [26] of 5% to further refine the data-set. This produced a list of the most significantly regulated genes (p values #0.0045). These are shown in Table 1, which lists the most significantly upregulated genes (fold-change limit + or -4.0) and Table 2 which contains the most significantly down-regulated (fold change limit + or -2.0) transcripts in the data-set. A number of these genes mapped to loci that predispose to EAE in the rat (http://rgd.mcw. edu), or have been associated with one or more of the commonly used EAE models, as shown in Tables 1 and 2.

b) Quantitative Real-Time PCR (rt PCR) Validation of Micro-array Data
To confirm our micro-array findings, expression levels of a the genes caspase 1 (Casp1), oligodendrocyte myelin glycoprotein (OMG), alpha subunit voltage-gated sodium channel type 1 (Scn1a), Fas (TNF receptor superfamily member 6) and superoxide dismutase 2, mitochondrial (Sod2) were validated using RT-PCR in 8 rats in each group. These genes were selected as being regulated in EAE and of biological interest. Gene transcripts were amplified from cDNA preparations. These were obtained from the original total RNA samples used in the micro-array analysis and an additional 4 normal rats with MBP-EAE and an additional 4 healthy rats. The results are shown in Figure 3. All transcripts Figure 1. Clinical course of EAE. A. The figure shows the results of a clinical course experiment (n = 10) in which female Lewis rats (10-12 wks old) were injected with 50 mg Guinea Pig MBP emulsified in incomplete Freund's Adjuvant containing 4 mg/ml heat inactivated M.Butyricium. The peak of disease generally occurs between 12 and 14 days post injection and resolves, usually within 7-10 days later. Mild relapses occur in a few animals between 28 and 40 dpi. B. Spinal cords used to prepare RNA for microarray analysis were collected on day 13 when the disease course was either at its peak or in the phase of very early resolution. Clinical course scores were recorded using a 12 point scale [43]. doi:10.1371/journal.pone.0048555.g001 amplified showed similar expression patterns in the RT-PCR assay compared to the micro-array analysis.

c) Biological Functions Most Significantly Regulated in MBP-EAE
The functional predictions that are generated from IPA are based on the direction of expression of a number of downstream genes which have been previously shown to be associated with that function. The functions that are most significantly altered in MBP-EAE are shown in Tables 3 and 4. The functions that are predicted to be up-regulated include activation of central nervous system cells, cell death or apoptosis (of leukocytes, microglia and dopaminergic neurons), formation of amyloid fibrils, and proliferation of oligodendrocyte precursor cells. The functions that are predicted to be down-regulated include long term-potentiation, neurotransmission, synaptic transmission, quantity of vesicles and delay in the death of neurons.

d) Canonical Pathways Regulated in MBP-EAE
A number of Canonical pathways were identified in IPA as being significantly regulated in MBP-EAE, based on the expression profile of genes in MBP treated animals compared to healthy controls. These are shown in Figure 4. These include immune related pathways and neural pathways. The details of some of the pathways that were regulated are shown in Figures S3, S4, S5, S6, S7, S8, S9, S10, S11, S12, S13, S14, S15, S16, S17, S18, S19, S20.

e) Upstream Factors Predicted to be Regulated in MBP-EAE
From IPA we have identified the upstream regulators that are predicted to be activated or inactivated in MBP-EAE. These are shown in Table 5. These include Sirtuin 2, stat1 and presenilin 2.

Discussion
In this study, using female Lewis rats, we have investigated the expression of genes in MBP-EAE compared to controls. We have studied rats at day 13, which is the peak of clinical disease, when processes that lead to recovery are under way. The pathology of MBP-EAE involves infiltration of the CNS with macrophages and activated T lymphocytes [27]. There are morphological changes in microglia and astrocytes and there is demyelination of axons and axonal transaction. Consistent with this, we have found that MBP-EAE is associated with changes in expression of genes and biological functions of the immune system and also of the nervous system. We found markedly increased expression of MHC class II molecules and other immune related genes. This is consistent with the pathological findings of upregulation of MHC class II molecules and infiltration with immune cells [14]. Micro-array analysis of C57/BL6 mice with EAE induced with myelinoligodendrocyte protein also found upregulation of numerous genes [23] including MHC class II and immune related genes. We note increased expression of lipocalin 2, which is a protein that can be involved in the immune response and that has been found to be upregulated in MOG-EAE in DA rats [28] and in C57/BL mice [29].
Genes involved in the inflammatory response and in antigen processing and presentation were some of the most highly regulated genes identified. Cell-cell signalling and interaction and cell movement were also highly represented amongst the most differentially regulated gene transcripts identified. A number of chemokines and their receptors, annexin 3A and fibronectin were all up-regulated more than four-fold in the diseased spinal cords. Cellular growth factors, notably IL2-a and the common gamma chain were highly up-regulated as was the inflammatory cytokine IL1-b which is expressed in activated microglial cells [30]. IL18, the IL1 receptor agonist IL1rn, the interleukin receptors IL1r1, Figure 2. Volcano plot. This shows the significance (y axis) versus the fold change (x axis) for genes in rats with EAE compared to healthy controls. doi:10.1371/journal.pone.0048555.g002 Ccl7 chemokine (C-C motif) ligand 7 6.7

2.65E-04
EAE18b Cellular movement Ccl7 maps to a gene region associated with MOG EAE in the rat [53] Ccr5 chemokine (C-C motif) receptor 5 6.5 2.25E-05 -Cellular movement Ccr5 was upregulated in the spinal cords of rats with MOG induced EAE [54] Cx3cr1 chemokine (C-X3-C motif) receptor 1 5.6 1.19E-03 -Cellular movement Cx3cr1 is involved in neuroinflammation and upregulated in phagocytic cells in the rat brain in MOG induced EAE.
[ 55;56] Cxcl9 chemokine (C-X-C motif) ligand 9 113.50 2.51E-06 -Cellular movement Cxc19 is involved in neuroinflammation. [55] Cxcl16 chemokine (C-X-C motif) ligand 16 20 IL4r, IL6r, and the interferon gamma inducible cytokine IL18 were all moderately up-regulated (Table S1). The regulation of key components of the classical (C1q, C4, C3b, Serping1) and alternative (C3, Cfb) complement pathways and the complement component 3a receptor 1 (3ar1), but not molecules in the terminal pathway, is in accordance with previous studies which indicate that the common terminal pathway (culminating in the membrane attack complex) is not important for the development of fulminant demyelinating EAE [31].
Genes associated with apoptosis which were highly regulated in our data set include the apoptosis related serine peptidases (Caspases) which have all been associated with disease in various EAE models [32;33] and the toll like receptors 2 and 4 (Tlr2, Tlr4) which have previously been found up-regulated in active lesions of MS and MOG induced EAE in the DA rat [34;35]. As mentioned above, the apoptotic inducer lipocalin 2 (Lcn2) [36] which has previously been identified as up-regulated in MOG induced EAE in DA rats [28] and C57BL mice [29] was also highly expressed in MBP-EAE in the present study.
Down-regulation of gene transcripts in MBP induced EAE in the Lewis rat was much more subtle than up-regulation however almost a thousand genes were found to be significantly downregulated in the model (Table S2). Eighteen of these were downregulated more than 2-fold, notably including key enzymes in the cholesterol biosynthesis pathway and molecules involved in cellular growth and proliferation, signalling and interaction.
Of the most highly regulated genes in the spinal cord (Tables 1  & 2), many are associated with processes associated with induction of disease, such as cell-cell signalling, chemotaxis antigen presentation, activation of both infiltrating cells such as Tlymphocytes and macrophages and resident astrocytes and microglia. There is also regulation of genes associated with resolution of disease, for example apoptosis of immune infiltrating cells and activated glial cells, and growth and proliferation of oligodendrocyte precursor cells which migrate to the sites of demyelination and differentiate into mature oligodendrocytes able to carry out the remyelination of axonal processes damaged during the inflammatory phase of the disease [37].
EAE can only be induced in susceptible animals. In studies of different mouse and rat strains, some genetic loci have been identified as being important in the predisposition to EAE. In our study, many of the genes identified as being differentially regulated found to map to those regions (http://rgd.mcw.edu/). From IPA, the cellular functions that are regulated can be predicted. The upregulated functions in MBP-EAE included activation of CNS cells, cell death or apoptosis (of immune cells, microglia and dopaminergic neurons), formation of amyloid fibrils, and proliferation of oligodendrocyte precursor cells. Apoptosis is known to be an important function in EAE. Apoptosis of T cells is important in the regulation of MBP-EAE [15;16]. Apoptosis of oligodendrocytes could also contribute to pathology in EAE. Consistent with this, we have found upregulation of genes associated with apoptosis. Our finding of increased proliferation of oligodendrocyte precursor cells is consistent with the observation that recovery in MBP-EAE is due to restoration of conduction in demyelinated fibres.
The functions that are predicted to be downregulated include long term-potentiation, neurotransmission, synaptic transmission, quantity of vesicles and delay in the death of neurons. These downregulated functions presumably reflect damage to neural tissue.
From IPA we also obtained a list of the canonical pathways that were regulated in MBP-EAE. There were changes in immune pathways and neural pathways. The immune related pathways Grb14 growth factor receptor-bound protein 14

5.55E-04
Neuinf1 Cell-cell signalling and Interaction Grb14 maps to a gene region regulating the inflammatory microglial response in the rat. [75] Hnmt histamine N-methyltransferase

1.74E-03
Neuinf 4 Cell-cell signalling and Interaction Ppp1r14c maps to a gene region regulating the inflammatory microglial response in the rat. [75] Trhr thyrotropin-releasing hormone receptor include antigen presentation, T cell differentiation and complement pathways, all of which would be expected to be activated in EAE. The neural related pathways include dopamine feedback of cAMP signalling, Creb signalling in neurones, biosynthesis of steroids and ALS signalling. From IPA we have identified the upstream regulators that are predicted to be activated or inactivated in MBP-EAE. The greatest inactivation was predicted in sirtuin 2 (sirt2). Sirt 2 is is a member of the sirtuin family of histone deacetylases and has been shown to have a role in promoting longevity [38]. With proteomic studies, sirtuin 2 has also been found to be downregulated in EAE [39]. Resveratrol which is an agonist of sirt2 has been shown to reduce the severity of EAE [40]. Stat1 is a member of the JAK/STAT group of proteins [41]. Stat 1 regulates genes involved in interferon signalling which is shown in this study to be upregulated in MBP-EAE. Presenilin 2 is a member of the gamma-secretase family and is predominantly expressed in neurones where it plays a role in the processing of amyloid precursor protein [42] and was predicted in this study to be downregulated.
In conclusion, we have shown that MBP induced EAE alters gene expression in the spinal cord. We have produced a gene expression profile for this model in the Lewis rat which correlates well with the documented characteristics of the disease. Further examination of the pathways and functions regulated should pave the way for future studies into their regulation which may lead to the identification of potential drug targets for enhancement of the resolution of inflammatory disease in the CNS.

Ethics Statement
This research was approved by the Animal Ethics Committee of the University of Queensland. The approval number was UQCCR/770/08/MSRA.

Induction of EAE
EAE was induced in 11 wk old female Lewis rats by subcutaneous injection at the base of the tail of 50 mg recombinant MBP, using our standard methods (McCombe et al., 1996). Healthy age-matched unimmunized rats were used as controls. Weakness of the tail, hindlimbs and forelimbs was assessed by grading the degree of weakness of each region separately on a scale of 0 (no weakness) to 4 (total paralysis), as described previously [43]. On each day of examination, the scores from each region were added together to give a total clinical score per rat (maximum total clinical score = 12). Figure 1 shows the typical clinical course of rats with this form of EAE.

Tissue Collection
On day 13 after inoculation, which is the peak of clinical disease, rats were sacrificed and the spinal cord was excised from the spinal column by insufflation using a large bore blunt ended needle and a 60 ml syringe. Figure 1B shows the clinical course of these rats up to the day of sacrifice. Tissue was snap frozen in liquid nitrogen prior to being stored at 280uC. Four rats in each group were used for micro-array analysis and a total of 8 rats in each group were used for rt-PCR.

RNA Extraction
Total RNA was extracted from tissues using QIAGEN RNeasyH Lipid tissue Midi kits as per the manufactures instructions and treated with DNase1 (QIAGEN) to remove all traces of genomic DNA prior to storage at 280uC. The RNA quality analysis was carried out using the BioRadExperion automated gel electrophoresis system (BioRad Laboratories Inc.).

Array Hybridization
For each group, 4 biological replicates were used. cDNA synthesis and amplification was performed using the Applause WT-Amp Plus ST kit (NuGEN). Samples were enzymatically fragmented and biotinylated using the Encore Biotin Module labelling kit (NuGEN). Samples were hybridized to Affimetrix-GeneChip Rat Exon 1.0 ST Arrays as per the manufacturer's instructions. Briefly, 5 mg of fragmented biotinylated ssDNA was hybridised for 16 hrs at 45uC, 60 rpm to the array chip. After 16 hrs GeneChips were washed on a FS_450 Fluidics station using the washing script FS450_0001 with buffers and stains supplied with the GeneChip Hybridisation, Wash and Stain Kit from Affymetrix.

Data Acquisition and Analysis
Data was acquired on a 7G GeneChip Scanner 3000 and.CEL file generation performed using AGCC. Expression Console with Robust Multi-chip Average (RMA) was used initially to extract probe intensity data. This data was used to access Affimetrix Tissue samples were snap frozen in liquid nitrogen and stored at 280uC prior to total RNA preparation using the QIAGEN RNeasy Lipid tissue kit. RNA quality analysis was carried out on the BioRadExperion automated electrophoresis system. All preparations used in both assays had RNA quality indicator (RQI) values of .9.5. For RT-PCR, total RNA was reverse transcribed and amplified as described in the methods. Analysis of selected genes up or down regulated at the peak of disease in EAE. Bars represent the average fold change between expression in the spinal cord level at peak of disease compared to normal healthy animals (+/2 SEMs, Microarray n = 4, RT-PCR n = 8). Dark columns represent fold change derived from the microarray data. Similar amplification patterns were obtained from RT-PCR amplification of the same total RNA samples and a second set of 4 animals samples at an identical time point. doi:10.1371/journal.pone.0048555.g003 Table 3. Cellular functions significantly up-regulated in MBP induced EAE in the Lewis rat. Kcnj10,Ace

Gjd2
Activation of astrocytes

Gene Level Expression Analysis
For gene transcript expression analysis, core probe sets were averaged for each gene. The significance of differences in gene expression between the EAE and healthy control data sets was evaluated by analysis of variance (ANOVA) that compared the controls and normal rats and investigated the effect of different batches. There was no significant batch effect in this study.
Although we are aware that the data may not be distributed normally and that a reduced homogeneity of variance may exist, the alternative non-parametric tests such as the Kruskal Wallis analysis would not produce the ratio of means required to identify statistically significant functional categories and pathways and make predictions about the possible downstream functional effects of fold change for various groups of genes. As a result, the above violations to the assumptions made by this analysis were considered to be adequately compensated for by the robust nature of the analysis. We have produced a list of all the genes that are regulated in EAE, and also a smaller list of the genes that are most highly regulated in MBP-EAE.

Functional Association Analysis
IPA (IngenuityH systems, www.ingenuity.com) was used to identify biological functions, gene networks and pathways which were significantly regulated in the spinal cord in EAE compared to healthy controls. In IPA, right-tailed Fisher's exact test was used to calculate a p-value determining the probability that each function network or pathway assigned to that data set is due to chance alone. Molecules from the data set that met the p,0.05 cut off and were associated with biological functions in the Ingenuity knowledge base were considered for the analysis.

Downstream Effects Analysis in IPA
The Ingenuity downstream effects analytic were used to identify biological functions which were predicted to be increased or decreased based on the observed gene expression changes in the data-set. The software compared the direction of change for each gene to the causal effects information between genes and biological functions stored in the IPA Knowledge base. Statistical validation of the predictions made were carried out based on the calculation of the regulation z score ''which is designed such that data sets composed of randomly chosen perturbed genes with random sign of fold change do not lead to significant results on average'' (Ingenuity Downsteam Effects Analysis, whitepaper). Each biological function considered is associated with a number of upstream genes based on annotations derived from the literature and stored in the Ingenuity knowledgebase. The previously observed effects on function of each gene in this set are then used to predict an increase or decrease in the function under the conditions studied. Table 3 shows the biological functions found to be significantly up or down-regulated in MBP induced EAE in the Lewis rat.

Quantitative RT-PCR Validation of Micro-array Data
Validation of the results obtained from the micro-array experiments was carried out using real time PCR. Five genes (rtCasp1, rtOmg, rtScn1a, rtFas and rtSod2) found to be differentially regulated in EAE in the micro-array gene expression   QIAGEN). At the end of each PCR run a melting curve analysis was carried out to verify that homogeneous products of equal size had been amplified for each primer pair. The same set of RNA samples that was used for the micro-array experiments was also analysed by RT-PCR and in addition we tested four more independent validation samples from each group, which had been collected and prepared in an identical manner, to test for sampling error. Expression was normalized to rat peptidylprolylisomerase H (rtPpih). Fold changes were calculated using the DDCT method [44] and expressed as fold-change at peak of disease in EAE compared to the expression levels in healthy control Lewis rats ( Figure 3). Figure S1 RNA Virtual Gel Image. Tissue samples were snap frozen in liquid nitrogen and stored at 280uC prior to total RNA preparation using the QIAGEN RNeasy Lipid tissue kit. RNA quality analysis was carried out on the BioRad Experion automated electrophoresis system. This shows an example of the virtual gel image, illustrate the high quality of samples used. All preparations had RQI values of .9.5.

Supporting Information
(DOCX) Figure S2 RNA Electropherogram Profile. Tissue samples were snap frozen in liquid nitrogen and stored at 280uC prior to total RNA preparation using the QIAGEN RNeasy Lipid tissue kit. RNA quality analysis was carried out on the BioRad Experion automated electrophoresis system. This shows an example of the electropherogram profiles illustrating the high quality of samples used. All preparations had RQI values of .9.5. Table S1 (XLSX)