Gene Expression in Spontaneous Experimental Autoimmune Encephalomyelitis Is Linked to Human Multiple Sclerosis Risk Genes

Recent genome-wide association studies have identified over 230 genetic risk loci for multiple sclerosis. Current experimental autoimmune encephalomyelitis (EAE) models requiring active induction of disease may not be optimally suited for the characterization of the function of these genes. We have thus used gene expression profiling to study whether spontaneous opticospinal EAE (OSE) or MOG-induced EAE mirrors the genetic contribution to the pathogenesis of multiple sclerosis more faithfully. To this end, we compared gene expression in OSE and MOG EAE models and analyzed the relationship of both models to human multiple sclerosis risk genes and T helper cell biology. We observed stronger gene expression changes and an involvement of more pathways of the adaptive immune system in OSE than MOG EAE. Furthermore, we demonstrated a more extensive enrichment of human MS risk genes among transcripts differentially expressed in OSE than was the case for MOG EAE. Transcripts differentially expressed only in diseased OSE mice but not in MOG EAE were significantly enriched for T helper cell-specific transcripts. These transcripts are part of immune-regulatory pathways. The activation of the adaptive immune system and the enrichment of both human multiple sclerosis risk genes and T helper cell-specific transcripts were also observed in OSE mice showing only mild disease signs. These expression changes may, therefore, be indicative of processes at disease onset. In summary, more human multiple sclerosis risk genes were differentially expressed in OSE than was observed for MOG EAE, especially in TH1 cells. When studying the functional role of multiple sclerosis risk genes and pathways during disease onset and their interactions with the environment, spontaneous OSE may thus show advantages over MOG-induced EAE.


INTRODUCTION
Although animal models are widely used in human research, it is still discussed whether they can adequately mirror diseases like multiple sclerosis (MS) that only exist in humans. MS is a chronic inflammatory disease of the central nervous system (CNS), with both environmental and genetic risk factors contributing to disease susceptibility. The recent identification of more than 230 genetic risk loci for MS (1,2) requires a reassessment of the widely used experimental autoimmune encephalomyelitis (EAE) animal models. To support analyses of the primary cause and etiology of MS, animal models should ideally replicate mechanisms taking place during MS disease induction.
Most EAE models are actively induced by injection of myelinderived antigens in conjunction with potent adjuvants (3). One such antigen is myelin oligodendrocyte glycoprotein (MOG), a component of the outer surface of myelin (4). Injection of the MOG 35−55 peptide into C57BL/6 mice leads to chronic EAE (5) and thus serves as a popular animal model to date. A related model, passively-transferred EAE, is caused by bulk transfer of in vitro-activated myelin-specific T cells (6).
By contrast, transgenic models such as opticospinal EAE (OSE) spontaneously develop autoimmune disease and may, therefore, be better suited to study disease onset than induced EAE is. Spontaneous models can be used for identifying environmental triggers of MS (7,8) and might support analyses of genetic risk factors for human MS. They circumvent problems specific to induced ones, such as adjuvant inoculation, with its partially unknown effects. In OSE, ∼50% of the animals develop a spontaneous inflammatory demyelinating CNS disease, predominantly affecting optic nerves and the lumbar part of the spinal cord (9). These mice carry two transgenic modifications: they express a T cell receptor (TCR) recognizing the MOG 35−55 peptide and B cells with MOG-specific receptors. In OSE, MOGspecific B cells function as antigen-presenting cells to trigger disease onset by activating MOG-specific T cells (10). Notably, B cell-depleting treatments for MS appear to target primarily cellular and not humoral B cell responses, and, thus, result in a reduced T cell activation (11).
For a long time, T H 1 cells were considered as the predominant drivers of EAE and MS (4). This hypothesis was challenged by emerging evidence for a substantial role of T H 17 cells in the disease etiology, including the discovery that the transfer of T H 17 cells can induce EAE. In fact, both T H cell types can induce EAE, albeit with distinct pathologies (12). In humans, genome-wide association studies (GWAS) have identified many MS risk loci that support a central role of T H cells and T H cell differentiation in the pathophysiology of MS (1,2,13).
Despite their valuable contributions to our understanding of MS pathophysiology and drug development, the relationship of EAE to human MS remains controversial (14). All available EAE models are, to some degree, artificial. Therefore, knowledge of whether gene expression changes in diseased mice involve MS risk genes can support the choice of an EAE model for specific research projects. The present study had three aims: First, to characterize gene expression differences in diseased OSE and MOG 35−55 EAE mice, two widely used EAE models with markedly different forms of induction. Second, to explore which of OSE or MOG-induced EAE resembles human MS more closely. To this end, we examined to which degree genes differentially expressed in spinal cord samples of OSE and MOG EAE showed significant enrichment of human MS risk genes. Third, to analyze expression differences of T H cell-specific transcripts in both EAE models.

Mice, Animal Handling, and Scoring
All mice used in this study had a C57BL/6 background and were bred in the animal facilities of the Max Planck Institute of Biochemistry and Neurobiology, Martinsried, Germany. For the OSE model, double-transgenic 2D2 (TCR MOG ) × IgH MOG (OSE) mice were used (9). For MOG EAE, wildtype C57BL/6 mice were immunized subcutaneously with 200 µg of a MOG peptide consisting of the amino acids 35-55, emulsified in complete Freund's adjuvant supplemented with 5 mg/ml Mycobacterium tuberculosis (strain H37Ra, Thermo Fischer Scientific BD Difco), as described previously (9). Pertussis toxin (400 ng, List Biological Laboratories) was injected intraperitoneally on the day of immunization and 48 h later. Control mice (CFA) received the same treatment but without the MOG peptide. For the analysis of EAE models, only female mice were used. For the T H cell analyses, OSE mice of mixed gender were used (15).
Scores for clinical signs of EAE were assessed daily according to the standard 5-point scale (9,16): 0: healthy animal; 1: animal with a flaccid tail; 2: animal with impaired righting reflex and/or gait; 3: animal with one paralyzed hind leg; 4: animal with both hind legs paralyzed; 5: moribund animal or death of the animal after preceding clinical disease. See Supplementary Figure 1 for the disease course of MOG EAE compared to control mice. Following our ethically approved protocol, the mice were sacrificed when they reached a score of 4. The animal welfare committee of the government of Upper Bavaria (Tierschutzkommission der Regierung von Oberbayern, Munich, Germany) approved the protocol. The animal procedures were in strict accordance with the guidelines set down by the animal welfare committee of the government of Upper Bavaria.
In vitro CD4 + T Cell Differentiation T cells derived from the spleen of a mixed-gender pool of four OSE mice were used to polarize pathogenic effector T H 1 and T H 17 cells, as described previously (15). In brief, four separate batches of four mice each were used for this experiment. To generate T H 1 cells, total erythrocyte-lysed spleen cells from OSE mice were cultured in the presence of a MOG peptide (amino acids 1-125), IL-12, IL-18, and anti-IL-4. After 3 days, IL-2 was added to the culture. To generate T H 17 cells, total erythrocytelysed spleen cells from OSE mice were cultured in the presence of a MOG peptide (amino acids 1-125), TGF-β1, IL-6, IL-23, anti-IL-4, and anti-IFN-γ. After 3 days, IL-23 was added to the culture. In both cases, cells were re-stimulated after 6 days and harvested after 9 and, once more, after 12 days. Naïve T H 0 cells were harvested on day 0. The success of polarization was evaluated by flow cytometry, ELISA, and quantitative real-time PCR (Supplementary Figure 2 and Supplementary Methods).

Microarrays
Two separate microarray experiments were performed on the Illumina gene expression profiling platform: The first comprised RNA isolated from total spinal cord preparations of healthy and diseased EAE mice. The second experiment analyzed gene  17) were used for the T H cell analysis. In each experiment, all samples and chips were processed in parallel. RNA processing, array hybridization, and quantification followed the same protocols in both experiments: First, concentration and purity of total RNA were assessed by 260 nm UV absorption and by 260/280 ratios, respectively (Nanophotometer, Implen, Munich, Germany). Second, RNA integrity was evaluated using a chip-based electrophoretic assay (Agilent RNA 6000 Nano Kit used in conjunction with the Agilent 2100 bioanalyzer, Agilent Technologies, Waldbronn, Germany). Mean RNA integrity numbers were 8.4 (SD = 0.5) for the EAE and 9.0 (SD = 0.5) for the T H cell experiment. Third, RNA was amplified and labeled using the Illumina TotalPrep RNA Amplification Kit (Ambion, Houston, TX, USA) and hybridized onto Illumina gene expression arrays following the manufacturer's instructions. Fourth, fluorescence signals were scanned on an Illumina BeadStation and analyzed by inhouse software routines. The manufacturer's built-in controls were analyzed, including hybridization controls and sampledependent parameters. All microarrays fulfilled Illumina's recommendations for quality control (QC).

Quality Control of Microarrays
Raw probe intensities were exported as summary data using Illumina's GenomeStudio, and further statistical processing was carried out using R v3.3.2 (17). For the analysis of EAE models, summary data was loaded using the Bioconductor package beadarray (18), and QC was conducted with lumi (19) and vsn (20). Each probe was transformed and normalized through variance stabilization and normalization. Probes were removed if they showed a detection p-value < 0.05 in >10% of the samples or had a "no match" or "bad" probe quality in the illuminaMousev2.db package. This procedure left 21,483 transcripts from 24 samples. For the T H cell experiment, summary data was loaded using limma (21) and QC was conducted with limma and vsn. Probes were transformed, normalized, and filtered as described above, based on the illuminaMousev1p1.db package. This pipeline left 17,858 transcripts. Technical batch effects were examined by inspecting the association of the first ten principal components of expression levels with expression chip and position on the chip.

Analysis of Differential Expression
Principal component analysis (PCA) was conducted in R using the function prcomp without scaling of variables; PCs were scaled for display. K-means clustering was performed using kmeans with k = 4; the analysis was repeated 100 times and the most stable clustering solution was selected. Differential expression was assessed with limma. For the analysis of differential expression across the EAE models, six mouse types were examined (with four mice each): wild-type (WT); healthy OSE controls (OSE 0 ); OSE with disease score 1 (OSE 1 ); OSE score 4 (OSE 4 ); as a MOG EAE control, healthy control mice injected with complete Freund's adjuvant but not with a MOG peptide (CFA); as MOG 35−55 EAE, C57BL/6 wildtype mice injected with adjuvant and MOG 35−55 peptide, rated score 4 (MOG 4 ). The design matrix was constructed from the six mouse types. Each expression chip contained one sample per mouse type. The four chips were added to the model as random effects via the duplicateCorrelation function. The five contrasts MOG 4 -CFA, CFA-WT, OSE 4 -OSE 0 , OSE 1 -OSE 0 , and OSE 4 -WT were computed on the fitted linear model and moderated t-tests were calculated using the eBayes function. For the T H cell experiment, the design matrix was constructed from the three cell types (naïve T H 0, T H 1, T H 17), with the four mouse pools treated as random effects. Only T H 1 and T H 17 cells harvested on day 9 were analyzed. The two contrasts T H 1-T H 0 and T H 17-T H 0 were examined.

Overrepresentation Analyses
Overrepresentation analyses (ORA) were conducted using WebGestalt v2019 (22) in R, based on the gene ontology (GO) biological process database. Genes were submitted as unique Entrez IDs, and the reference was genome protein-coding. The significance level was determined using a hypergeometric test, followed by calculation of the Benjamini-Hochberg false discovery rate (FDR) (23).

Enrichment Tests
The enrichment of genes was calculated using a permutation test in R. For this test, the unique Entrez IDs of genes were used. First, the amount of unique differentially expressed genes was determined, and the same number of random genes was selected. Second, the number of these random genes overlapping with the test set of genes (e.g., MS susceptibility genes) was determined. These two steps were repeated 100,000 times. Third, to calculate a p-value, the number of observations where the overlap between random genes and test genes was equal to or larger than the overlap between differentially expressed genes and test genes was counted and divided by the number of permutations.
For the enrichment analysis with MS susceptibility genes, the 558 genes outside the major histocompatibility complex (MHC) region listed in Supplementary Table 18 of the MS genomic map published in 2019 by the IMSGC were used (2). From this list, CTB-50L17.10, RP11-345J4.5, JAZF1-AS1, ZEB1-AS1, GATA3-AS1, SSTR5-AS1, and RPL34-AS1 were excluded to generate the list of 551 prioritized putative MS susceptibility genes described in the IMSGC publication.

RESULTS
We compared gene expression profiles of total spinal cord preparations derived from two EAE models, OSE and MOG 35−55 EAE. Double-transgenic OSE mice developed CNS autoimmunity spontaneously, predominantly affecting the lumbar part of the spinal cord. In the MOG 35−55 EAE model, the disease was induced in C57BL/6 wildtype (WT)  Table 1). Because of the spontaneous nature of disease development in OSE mice, gene expression in diseased OSE 4 animals showed more variance than was, e.g., observed in MOG 4 mice, which exhibit a more stereotypic disease course (9). PC, principal component; SD, standard deviations. (B) Venn diagram highlighting the number of transcripts differentially expressed in the analyzed contrasts. For this plot, up-and downregulated transcripts were analyzed separately, and transcripts differentially expressed in opposing directions are therefore included in the counts. (C,D) OSE 4 mice showed greater fold changes of expression levels than (C) OSE 1 and (D) MOG 4 mice, each compared to its control condition (OSE 0 and CFA, respectively) (Supplementary Table 3). For each group, the ten most differentially expressed genes (with Entrez IDs) are labeled. If two probes per gene were present among the top differentially expressed transcripts, the gene was counted only once, but both probes are plotted. The groups are: differentially expressed in OSE 1 only (light magenta), differentially expressed in OSE 4 only (dark magenta), differentially expressed in MOG 4 only (red), differentially expressed in (C) both OSE 1 and OSE 4 or (D) both MOG 4 and OSE 4 (brown; with higher expression levels observed for OSE 4 ).
mice by immunization with a MOG 35−55 peptide. PCA of gene expression profiles separated healthy (OSE 0 , CFA, and WT) from diseased [OSE score 1 (OSE 1 ), OSE score 4 (OSE 4 ), and MOG score 4 (MOG 4 )] animals along the first component ( Figure 1A). Most variance in gene expression was thus observed between healthy and diseased mice and not between EAE models. Because of the spontaneous nature of disease development in OSE mice, gene expression in diseased OSE 4 animals showed more variance than was, e.g., observed in MOG 4 mice, which exhibit a more stereotypic disease course (9).
Unsupervised k-means clustering on PCs further supported this finding, which consistently (in 97 of 100 replications) placed healthy and diseased animals into separate clusters (Figure 1A, Supplementary Table 1). The most frequent cluster solution (34/100) placed all healthy mice together in cluster 1; additional clusters were OSE 1 only (cluster 2), OSE 4 only (cluster 3), and a mixed cluster of the remaining diseased animals (cluster 4). We could thus successfully detect disease-relevant gene expression changes in the animals.

Stronger Gene Expression Changes in the OSE Model
Next, we analyzed gene expression changes in OSE and MOG EAE mice. We examined differential expression for five contrasts: OSE 1 -OSE 0 , OSE 4 -OSE 0 , MOG 4 -CFA, and the two control contrasts CFA-WT and OSE 0 -WT (Figure 1B, Supplementary Table 2). In the control contrast CFA-WT, no transcript was differentially regulated. A single transcript was upregulated in OSE 0 -WT, T cell receptor alpha chain (Tcra), which was also upregulated in all other contrasts except CFA-WT. The number of significantly up-and downregulated transcripts was higher for OSE 4 -OSE 0 (n = 5,555) than for MOG 4 -CFA (n = 3,182). In total, the expression of 864 transcripts differed significantly between MOG 4 and OSE 4 mice (Supplementary Table 2). Interestingly, 4.88× more transcripts were differentially expressed specifically only in OSE 4 -OSE 0 than only in MOG 4 -CFA ( Figure 1B). Moreover, fold changes were higher in OSE 4 -OSE 0 than in either OSE 1 -OSE 0 (binomial test: p = 1.4 × 10 −65 for all transcripts, p = 9.9 × 10 −119 for transcripts differentially expressed in both contrasts, Figure 1C) or MOG 4 -CFA (p = 5.8 × 10 −3 for all, p = 2.7 × 10 −221 for differentially expressed transcripts, Figure 1D; Supplementary Table 3). Stronger global gene expression changes were thus triggered in OSE than in MOG EAE.

Overrepresentation of Immune System Processes Especially for OSE
To characterize the expression changes in the different EAE models further, we conducted ORA analyses of the analyzed contrasts (Supplementary Table 4 Table 4). Together with other immune-related gene sets, immune response, regulation of immune system process, and T cell activation were among the top-associated terms (adjusted p < 2 × 10 −16 ). These and other immune-associated processes remained significant in OSE 4 sp (adjusted p ≤ 3.5 × 10 −2 , Supplementary Figure 4). By contrast, no immune systemspecific process was significant for MOG 4 sp. More expression changes in the immune system were, therefore, triggered in OSE than in MOG EAE.

Activation of the Adaptive Immune System in OSE 1 Mice
While MOG EAE develops rapidly in a highly stereotypical manner, the clinical course of OSE is usually slower and shows more inter-individual variability (24). OSE thus allows for studying disease at different stages, and we analyzed mice showing a mild disease score of 1 (OSE 1 ). Compared to OSE 0 , 34 transcripts were differentially expressed specifically in OSE 1 animals and not in any other contrast [OSE 1 -specific transcripts (OSE 1 sp), Supplementary  Figure 4). Furthermore, the gene sets B cell mediated immunity and antigen processing and presentation were significantly overrepresented not only in the analysis of CDT but also for the OSE 1 ex transcripts, indicating a potential role of B cells also in mildly affected OSE mice.

Enrichment of MS Susceptibility Genes Among Transcripts Expressed in OSE
Over 230 independent genetic loci associated with MS susceptibility in humans have been identified (1,2). Based on these GWAS loci, 551 human MS susceptibility candidate genes have been proposed (2), for which expression data of 499 transcripts were available in our dataset. We conducted a PCA on these transcripts (265 genes) to analyze whether the expression of MS risk genes was increased in the EAE models. The first component, explaining 75.7% of the variance in expression of these transcripts, was significantly higher in all disease groups than in controls, indicating high expression levels of MS-associated genes in EAE, with the highest levels observed for OSE 4 (Figure 2A, Supplementary Table 6). Also individual MS risk genes, e.g., H2-Ab1, Cd52, and Cd86 (1, 2), as well as further putative MS-associated genes like Cd74, were among the transcripts showing the lowest differential expression  Table 2). In all three cases, diseased mice showed significantly higher expression levels, with the highest expression observed for OSE 4 mice. Significance levels: *** adjusted p < 0.001. p-values. They were significantly upregulated in all three diseased mouse types (Figures 1C,D, 2B-D, Supplementary Figure 5 Table 8). Significance levels: n.s.: p ≥ 0.05, * adjusted p < 0.05, *** adjusted p < 0.001.

Gene Expression in OSE Overlaps With T H Cell-Specific Transcripts
T H cell differentiation was identified as a key pathway in the etiology of MS (13). We, therefore, analyzed whether gene expression changes in EAE models were related to T H cell differentiation. To this end, gene expression profiling of in vitro polarized T H 1 and T H 17 cells was conducted, derived from OSE mice. Compared to naïve T H 0 cells, 8× more transcripts were differentially expressed specifically in T H 1 than in T H 17 cells (Figure 3A). None of the transcripts differentially expressed in both T H 1 and T H 17 were regulated in opposite directions. We examined via PCA whether the expression of T H 1-and T H 17-specific, differentially expressed probes was higher in EAE models than controls. The first component of T H 1-and T H 17specific gene expression explained 49.6 and 68.6% of the variance, respectively. For both T H cell types, the first component of cellspecific transcripts was significantly higher in all disease groups than in controls, with the highest levels for OSE 4 (Figures 3B,C,  Supplementary Table 7). Among signature molecules for T H 1 cells, Tbx21 (T-bet) was significantly upregulated in all diseased mice, Ifng only in OSE 4 ( Figure 3D, Supplementary Table 8). Of the examined T H 17 markers, only Il17f was upregulated in OSE 4 , neither Rorc nor Il17a were differentially expressed ( Figure 3E, Supplementary Table 8).
After correction for multiple testing, the CDT, OSE 4 sp, and OSE 1 ex analysis sets were significantly enriched for T H 1-and T H 17-specific transcripts ( Table 2). In the case of MOG 4 sp transcripts, the overlap was lower and only significant for T H 1specific probes. These experiments indicate a stronger overlap of known MS-associated immune responses involving T H cells with OSE than with MOG EAE.
Finally, we analyzed whether EAE-associated genes differentially expressed in T H 1 or T H 17 cells were more closely connected to human MS. To this end, we intersected the lists of EAE-specific and T H -specific transcripts. Immunerelated biological processes were overrepresented for CDT, OSE 4 sp, and OSE 1 ex genes intersected with T H 1-specific genes. (Supplementary Table 9, Supplementary Figure 6). No terms were significantly overrepresented for any T H 17-specific or MOG 4 sp genes.
CDT and OSE 1 ex genes differentially expressed in T H 1 cells were significantly enriched for the IMSGC MS risk genes (p < 7 × 10 −4 , Table 3). The enrichment for OSE4sp did not withstand correction for multiple testing. Neither any of the T H 17-specific gene sets nor the genes from the MOG 4 sp group were enriched for these risk genes. Thus, we conclude that OSE entails gene expression changes involving human MS gene risk genes, especially in T H 1 cells, which were not observed to the same degree for MOG EAE.

DISCUSSION
With the identification of over 230 MS risk loci in recent GWAS, we move closer to understanding the etiology of MS. Further research relies on adequate animal models that have to be reassessed in the context of GWAS data. Given the interplay of genetics and environment in human MS, spontaneous EAE models like OSE might be more apt for studying the genetic risk component of MS than induced EAE models that require active experimental manipulation. In the present study, we performed spinal cord gene expression profiling to, first, characterize differences between spontaneous OSE and MOG-induced EAE and, second, to analyze the relationship of both models to human MS risk genes and T H cell biology.

OSE May Reflect the Etiology of MS Better Than MOG EAE Does
In comparison to MOG EAE, gene expression changes in OSE were stronger and more closely linked to immune pathways. This might reflect a more complex mode of disease induction in OSE than is the case for MOG EAE. OSE features active B and T cell cooperation, a mechanism highly relevant for the pathophysiology of human MS, as demonstrated by the effectivity of B cell-depleting treatments (10,11). More than MOG EAE, OSE-specific transcripts were enriched for both human MS risk genes and T H cell-specific transcripts and showed    an overrepresentation of immune-specific gene sets. We thus hypothesize that OSE shows advantages over MOG EAE in studying the functional role of human MS risk genes and their associated immune pathways. Nevertheless, many of the differentially expressed genes indicate that both EAE models faithfully recapitulate critical functional pathways of MS, especially regarding the role of antigen presentation and CD4 + T cells in MS immunopathogenesis (25,26). Transcripts for the HLA genes H2-Eb1 and H2-Ab1, homologous to HLA-DRB5 and HLA-DQB1, were among the most differentially expressed probes. The alleles HLA-DRB5 * 01:01 and HLA-DQB1 * 06:02 are part of the DR15-DQ6 haplotype and are, most likely because of linkage disequilibrium with HLA-DRB1 * 15:01, strongly associated with MS risk (27). In MS, memory B cells mediate autoproliferation of brain-homing T H 1 cells in a HLA-DRB1 * 15:01-dependent manner (28). Interestingly, the antigen-presenting function of MOG-specific B cells is, in cooperation with T cells, important for the development of OSE (10). Among putative non-MHC MS risk genes (1, 2), Cd86, Cd52, and Cd74 showed very robust support for differential expression.
We could thus show that EAE, and in particular OSE, constitutes a valuable model for studying the role of human MS risk genes. Several previous studies support this finding: First, humanized EAE models successfully replicated HLA-related risk variants, including HLA-DRB * 15:01 (29). Second, knockout mice lacking the MS-associated Il7r are resistant to EAE (30). Third, shared human and EAE risk loci exist that are linked to T H cell differentiation (31). Fourth, an overlap of upregulated genes between myelin-reactive T cells from MS patients and encephalitogenic CD4 + T cells isolated from EAE was described (32). Fifth, in a passive-transfer EAE study, several MS risk genes were suggested to be implicated in the transition from in vitrogenerated MOG-specific T H 17 cells to encephalitogenic CD4 + T cells (33).
Functional pathways involving MS risk genes interact with environmental factors to trigger an autoimmune response, as demonstrated by the role of epigenetic factors for MS risk (1,34). Spontaneous EAE models might resemble gene-environment interactions more faithfully than MOG EAE does. For instance, in a spontaneous EAE model, disease onset could be prevented in mice kept under germ-free conditions (7). In this model, a higher incidence of EAE was observed following the transfer of the human gut microbiome from MS patients than when transferring the microbiome from the patient's healthy twin (8).

T H 1-Specific Transcripts Are Enriched for MS Risk Genes
Our gene set analyses point at a central role of lymphocyte activation in EAE induction and shed light on the ongoing controversy regarding the relative importance of T H 1 and T H 17 cells in mediating CNS autoimmunity (35). In accordance with previous studies (9), we observed a higher differential expression of selected T H 1-than of T H 17-specific transcripts in diseased mice. Interestingly, a high T H 1/T H 17 ratio is indicative of a lesion distribution pattern characterized by prominent spinal cord involvement, as is the case for both EAE models investigated in our study (12,15,36).
CDT and OSE 1 ex transcripts differentially expressed in T H 1 cells were significantly enriched for MS risk genes ( Table 3). We did not observe such an enrichment for transcripts differentially expressed in T H 17 cells. Albeit also OSE 4 sp genes were only enriched for risk genes in T H 1 cells at nominal significance (unadjusted p = 0.0097), T H 1-expressed MOG 4 sp transcripts showed no trend for the enrichment of MS risk genes at all (unadjusted p = 0.51). In GO overrepresentation analyses, immune-related biological processes like positive regulation of T cell proliferation were significant for OSE 4 sp-genes differentially expressed in T H 1 cells, but no GO gene sets at all were overrepresented in T H 1-specific MOG 4 sp genes. In the context of T H 1-driven immune responses, the OSE model might thus be linked more closely to human MS risk genes than MOG 35−55 EAE is. However, T H 17 cells can shift toward a T H 1 phenotype in EAE (37,38). The T H 1 markers analyzed in the EAE models may, accordingly, reflect expression in a significant proportion of former T H 17 cells. Therefore, our findings do not argue against a relevant impact of T H 17 cells in either EAE model.

Expression Patterns Across Different Disease Stages Can Be Studied Using OSE
Most genes differentially expressed in OSE 1 mice were also recapitulated in severely affected OSE 4 mice and showed the same direction of regulation in both disease stages. Many factors active in severe EAE thus also influence EAE during a mild or, potentially, early disease course. Effective immunotherapy is facilitated if the same biological pathways are continuously active throughout the entire disease. For example, the gene set response to interferon-beta was highly overrepresented in both OSE 1 ex and CDT and Cd52 was differentially expressed in all diseased mice. Studying mild OSE cases might, therefore, constitute an interesting model for defining the initial triggers of MS and the identification of novel therapeutic options.

LIMITATIONS
Our gene expression analysis of two EAE models had several limitations: First, the microarrays used covered only part of the murine transcriptome and thus, some MS risk genes could not be analyzed. Second, the statistical power of our analyses was restricted by the sample size. Third, the initial phases of EAE are hard to define since the disease develops over a short period. We thus analyzed mild OSE cases as a proxy for early disease. It is, however, unknown whether these animals would have developed more severe EAE later.

CONCLUSIONS
Although hundreds of genetic MS risk loci have been identified, their functional role in the etiology of the disease still has to be resolved. Ideally, suitable animal models recapitulate molecular and functional pathways involving these genes. They may thus move research closer to the primary cause and etiology of MS, thereby supporting the identification of effective immunotherapies. No animal model fully reflects a heterogeneous human disease like MS and each EAE model available today only replicates a part of the human disease. Researchers will thus continue to study different aspects of MS using a variety of EAE models. Our results indicate that OSE, with its closer link to MS risk genes and T H cell biology, may be better suited for studying the etiology of MS and for defining specific therapeutic targets than MOGinduced EAE is. Future studies will show whether OSE can fulfill this promise to model the human MS genetic risk landscape faithfully.

DATA AVAILABILITY STATEMENT
The datasets presented in this study can be found in online repositories. The names of the repository and accession number(s) can be found at: https://www.ebi.ac.uk/arrayexpress/, E-MTAB-9132; https://www.ebi.ac.uk/arrayexpress/, E-MTAB-9133.

ETHICS STATEMENT
The animal study was reviewed and approved by Tierschutzkommission der Regierung von Oberbayern, Munich, Germany.

AUTHOR CONTRIBUTIONS
HF, GK, PW, and FW contributed to the original conception and design of the study. HF, GK, and PW conducted experiments. DK and TA devised the statistical analyses. DK, BP, and TA conducted statistical analyses. BM-M and FW supervised the study. HF and TA drafted the manuscript. All authors contributed to manuscript revision, read, and approved the submitted version.  Table 4) for the transcript groups (A) CDT, common disease transcripts (differentially expressed for both contrasts OSE 4 -OSE 0 and MOG 4 -CFA but not in the two control contrasts OSE 0 -WT or CFA-WT), (B) OSE 4 sp, OSE 4 -specific transcripts (differentially expressed for the contrast OSE 4 -OSE 0 but not in MOG 4 -CFA, OSE 0 -WT, or CFA-WT), (C) OSE 1 ex, OSE 1 -expressed transcripts (differentially expressed in OSE 1 -OSE 0 but not in OSE 0 -WT or CFA-WT). Note that no GO terms that are descendants of the term Immune System Process were significantly overrepresented for the group MOG 4 sp, MOG 4 -specific transcripts (differentially expressed for the contrast MOG 4 -CFA but not in OSE 4 -OSE 0 , OSE 0 -WT, or CFA-WT). The -log 10 (FDR) from hypergeometric tests is shown on the x-axis and used for coloring the plots (darker colors represent lower FDRs).  Table 2). Significance levels: * adjusted p < 0.05, * * adjusted p < 0.01, * * * adjusted p < 0.001.
Supplementary Figure 6 | The top 40 overrepresented immune system pathways in the differentially expressed transcripts groups intersected with T H 1-specific genes. The plots show the top 40 overrepresented GO terms that are descendants of the term Immune System Process (Supplementary Table 9) for the transcript groups (A) CDT (common disease transcripts) intersected with T H 1-specific genes, (B) OSE 4 sp (OSE 4 -specific transcripts) intersected with T H 1-specific genes, (C) OSE 1 ex (OSE 1 -expressed transcripts) intersected with T H 1-specific genes. Note that no GO terms were significantly overrepresented for any T H 17-specific or MOG 4 sp genes. The -log 10 (FDR) from hypergeometric tests is shown on the x-axis and used for coloring the plots (darker colors represent lower FDRs).

Supplementary
Supplementary Table 4 | Results from WebGestalt gene ontology (GO) over-representation analyses (ORA) for the analysis groups CDT, OSE 4 sp, MOG 4 sp, and OSE 1 ex (see the legend of Table 1 and Supplementary Figure 5) as well as the contrasts OSE 1 -OSE 0 , OSE 4 -OSE 0 , MOG 4 -CFA, and MOG 4 -OSE 4 (Supplementary Figure 4). FDR, 5% false discovery rate. Gene sets that are labeled in bold if they are descendants of the GO term "Immune System Process".
Supplementary Table 5 | Differential expression results for the analysis groups (see the legend of Table 1) CDT, CDT, OSE 4 sp, MOG 4 sp, OSE 1 ex, and OSE 1 sp. OSE 1 sp transcripts were differentially expressed in OSE 1 -OSE 0 but not in any other contrast. FC, fold change, CI, 95% confidence interval. For CDT transcripts, separate coefficients are provided for the two contrasts MOG 4 -CFA and OSE 4 -OSE 0 .
Supplementary Table 6 | Results from the principal component analysis (PCA) of gene expression profiles of putative MS risk genes (Figure 2A). FC, fold change, CI, 95% confidence interval. P-values were adjusted for multiple testing using Holm-Bonferroni correction. Table 7 | Results from the principal component analysis (PCA) of gene expression profiles of T H cell-specific transcripts ( Figures 3B,C). FC, fold change; CI, 95% confidence interval. P-values were adjusted for multiple testing using Holm-Bonferroni correction.  Figures 3D,E). DE, differentially expressed.

Supplementary
Supplementary Table 9 | Results from Webgestalt gene ontology (GO) over-representation analyses (ORA) for the analysis groups CDT, OSE 4 sp, and OSE 1 ex (see the legend of Table 1) among T H 1-specific transcripts (Supplementary Figure 6). For MOG 4 sp or T H 17-specific transcripts, no significant gene sets were found. FDR, 5% false discovery rate. Gene sets that are labeled in bold if they are descendants of the GO term "Immune System Process".