Natural genetic variation profoundly regulates gene expression in immune cells and dictates susceptibility to CNS autoimmunity

Abstract

Regulation of gene expression in immune cells is known to be under genetic control, and likely contributes to susceptibility to autoimmune diseases such as multiple sclerosis (MS). How this occurs in concert across multiple immune cell types is poorly understood. Using a mouse model that harnesses the genetic diversity of wild-derived mice, more accurately reflecting genetically diverse human populations, we provide an extensive characterization of the genetic regulation of gene expression in five different naive immune cell types relevant to MS. The immune cell transcriptome is shown to be under profound genetic control, exhibiting diverse patterns: global, cell-specific and sex-specific. Bioinformatic analysis of the genetically controlled transcript networks reveals reduced cell type specificity and inflammatory activity in wild-derived PWD/PhJ mice, compared with the conventional laboratory strain C57BL/6J. Additionally, candidate MS-GWAS (genome-wide association study candidate genes for MS susceptibility) genes were significantly enriched among transcripts overrepresented in C57BL/6J cells compared with PWD. These expression level differences correlate with robust differences in susceptibility to experimental autoimmune encephalomyelitis, the principal model of MS, and skewing of the encephalitogenic T-cell responses. Taken together, our results provide functional insights into the genetic regulation of the immune transcriptome, and shed light on how this in turn contributes to susceptibility to autoimmune disease.

Figure 1. Genetic control of gene expression in immune cells

(A) Top significantly DE probes between B6 and PWD (for each strain biological replicates, male, n=3; female n=3) immune cells exhibiting |FC| > 4, at an FDR < 0.05 in any of the five cell types are shown. FC indicates expression in PWD relative to B6. Mean indicates the FC average across all five cell types, which was used to sort the order of the probes. The |FC| > 4 cutoff was chosen to highlight the significantly DE genes with the largest effect size, and to facilitate visualization of the data. (B) Top 40 DE genes in PWD vs. B6 cells, identified by top 20 maximum FC (upregulated) and top 20 minimum negative FC (down-regulated) in any cell type. (C) and (D) Distribution and overlap of DE genes between selected cell types is shown. The number of genes indicates DE genes passing the filter of |FC| > 2 and FDR < 0.05, darker shading indicates a higher number of genes.

Figure 2. Predicted activation state of PWD immune cells

(A) Canonical pathway analysis of differential gene expression between PWD and B6 immune cells. The top significantly enriched canonical pathways (p < 0.01) are shown. The heat map indicates the Z-score, indicative of predicted direction of change (orange – upregulated; blue, downregulated). (B) Upstream analysis of differential gene expression between PWD and B6 immune cells. Top predicted activators (genes and proteins only) (p < 0.01) are shown. Enrichment analysis of MS-GWAS (C) or MS-CIS (D) genes within the DE transcripts between B6 and PWD was performed as described in the Materials and Methods. Enrichment p values [displayed as negative log(p)] for each cell type are shown. (E) Gene set enrichment analysis of DE genes in CD4 T cells and APCs was performed as described in Materials and Methods. Top panels show genes upregulated in PWD, the bottom show downregulated genes. (F) Quantitation of cell lineage-specific gene expression (as a measure of enrichment) across different cell lineages was performed as described in the Materials and Methods. Significance of differences between genes upregulated in PWD relative to B6 vs. those downregulated is indicated as follows: *, p<0.05; $, p<0.01, #, p<0.001. Abbreviations, Stem – stem cells, GN – granulocytes, TC – T cells, Str – Stroma, ILC = innate-like lymphoid cells.

Figure 3. Genes exhibiting significant sex-by-strain interaction across five cell types

(A) Genes exhibiting sex-by-strain interactions were identified as outlined in the Materials and Methods. Genes passing the filter of |FC| > 1.5 and FDR < 0.05 are shown. FC represents the change in male:female ratio between PWD and B6, calculated as (FC PWD^Male – PWD^Female) – (FC B6^Male – B6^Female), see Materials and Methods. Thus, a positive value is indicative of more male biased expression of a gene in PWD compared to B6. (C) and (D) Relative expression values of the indicated genes exhibiting significant sex-by-strain interaction in CD8 cells. Relative expression values represent log2-scaled normalized raw expression values. Error bars indicate standard error of the mean.

Figure 4. PWD mice are resistant to EAE

EAE was induced and evaluated B6 (female, n=6; male n=10) and PWD (female, n=5; male n=15) mice as described in the Materials and Methods. The following EAE quantitative traits were calculated: (A) incidence, (B) cumulative disease score (CDS), (C) severity index, (D) peak score, (E) days affected, and (F) day of onset. Significance of differences in (A) was determined by Fisher’s exact test. Significance of differences in (B–F) was determined by two-tailed Student’s t-test. Significance of differences between B6 and PWD is indicated using asterisks as follows: *, p < 0.05; **, p<0.01, ***, p<0.001. Error bars indicate standard error of the mean. The data are pooled from two independent experiments, both of which yielded similar results.

Figure 5. PWD mice display skewed peripheral immune responses

B6 (n=6) and PWD (n=9) mice were immunized as in Fig. 3. At day 10 post-immunization, cells were isolated from the spleen (A, B, E) and draining lymph nodes (C, D, E), restimulated with PMA/Ionomycin (except in E), stained for surface markers, followed by fixation and intracellular staining for the indicated cytokines (A–D) or FoxP3 (E), and flow cytometric analysis. Percentages of cytokine positive cells among the live CD19-TCRβ+CD4+ (A, C, E), or CD19-TCRβ+CD8+ (B, D) populations are shown. Significance of differences was determined by two-way ANOVA with Bonferroni’s multiple comparison test. The data represent one independent experiment.

Figure 6. PWD mice display reduced immune responses in the CNS during EAE

B6 (n=4) and PWD (n=6) mice were immunized as in Fig. 3. At day 30 post-immunization, mononuclear cells were isolated from the CNS by Percoll gradient, counted, and enumerated by flow cytometry. Numbers of cells positive for the indicated markers were calculated by multiplying the total number of isolated mononuclear cells by the percentage of CD45+ cells (A) or by the percentage of CD45+TCRβ+ cells (B). In (C), mononuclear cells were restimulated with PMA/Ionomycin and analyzed by intracellular staining and flow cytometry, as in Fig. 4. Percentages of CD45+TCRβ+CD4+ cells positive for the indicated cytokines are shown. Significance of differences was determined using the Student’s t-test. The data represent one independent experiment.