Defining the transcriptional and cellular landscape of type 1 diabetes in the NOD mouse

Abstract

Our ability to successfully intervene in disease processes is dependent on definitive diagnosis. In the case of autoimmune disease, this is particularly challenging because progression of disease is lengthy and multifactorial. Here we show the first chronological compendium of transcriptional and cellular signatures of diabetes in the non-obese diabetic mouse. Our data relates the immunological environment of the islets of Langerhans with the transcriptional profile at discrete times. Based on these data, we have parsed diabetes into several discrete phases. First, there is a type I interferon signature that precedes T cell activation. Second, there is synchronous infiltration of all immunological cellular subsets and a period of control. Finally, there is the killing phase of the diabetogenic process that is correlated with an NF-kB signature. Our data provides a framework for future examination of autoimmune diabetes and its disease progression markers.

Figure 1. Examination of NOD islets throughout diabetogenesis.

(A) Islets of Langerhans were isolated from NOD mice at the indicated ages and stained for blood vessels (PECAM-1), intra islet myeloid cells (CD11c), and T cell (CD4). Shown are representative images obtained from a pool of 6 mice per age from two independent experiments. Insets show contacts between intraislet myeloid cells and T cells. White bars represent 50 µm. (B) Islets were isolated and stained for the indicated markers and then scored for presence or absence of staining. Bars represent mean+/−S.D. of the percentage of marker positive islets of Langerhans obtained from a pool of 6 mice per group from two independent experiments. (C) Islets of Langerhans were dispersed and cells were examined by flow cytometry for the indicated cell surface markers. Bars represent the mean+/−S.D. of the percentage of total islet cells identified in two independent experiments per age. Results were obtained from a pool of 8 to 10 mice per group.

Figure 2. Pairwise and principal component analysis of microarray data.

(A–B) Scatter-plots of the normalized probe intensity of all annotated microarray signals are shown. Each dot represents the mean of 6 independent biological replicates. Numbers in the box represent the number of features that were at least 2-fold different at a 99% confidence interval by moderated t test with Benjamini-Hochberg false discovery rate analysis. Data are plotted at a log₂ scale. Panel (A) compares NOD.Rag^−/− mice versus NOD mice at 2 wks of age and panel (B) compares them at 6 wks of age. (C) Hierarchically clustered heat map (Euclidean row distance) of the 24 differentially expressed transcripts identified in (B). (D–E) Scatter-plots as in (A–B) except the comparison is between 2 and 6 wk old NOD or 2 and 6 wk NOD.Rag^−/−. (F) GO terms and transcription factor binding analysis of the differences identified in (D–E). The two graphs on the left, labeled ‘Shared’ show signatures common to both NOD and NOD.Rag^−/−. The signatures on the right, labeled ‘NOD only’, show changes specific to the NOD strain. (G) Principal component analysis of microarray samples. Each of 57 microarrays is summarized as a point and drop-line. Samples are color coded as indicated by the title of each group. (H) Hierarchical clustered heat map (Euclidean distance) of the top 1% most variant genes identified by principal component. (i) GO and transcription binding analysis of the top 1% most variant genes identified by principal component analysis.

Figure 3. Identification of significant gene changes in different aged NOD mice.

Examination of transcriptional changes that took place between 2 wks and either (A) 4 wks, (B) 12 wks, or (C) newly diabetic NOD mice. The red portion of each Venn diagram shows the genes identified by Pearson’s correlation as following a positive correlation throughout diabetogenesis. The green portion of each Venn diagram shows statistically significant changes from 2 wks to the given wk as determined by 2-fold upregulation and 99% confidence interval by moderated t test with Benjamini-Hochberg false discovery rate analysis. The blue portion of the Venn diagrams is the intersection of all the pairwise statistically significant differences from 2 wks to the indicated ages. Hierarchically clustered heat maps show the Euclidean distance of genes identified by the yellow intersection in the Venn diagrams. Those were the genes that showed both a positive correlation by Pearson’s correlation and a pairwise fold and statistical change at the indicated time. Gene names in red are type I interferon-inducible and those in green are inducible by both type I and type II interferons.

Figure 4. Analysis of inflammatory genes changes throughout diabetogenesis.

Differentially expressed genes identified in Figure 3 and Supplemental Table 1 were interrogated for transcriptional regulation signatures and gene ontology (A) or for their immunological role (B). (A) Cell enrichment, transcription factor binding, and gene ontology analysis were performed. Numbers in parentheses indicate the number of cell type specific genes identified as statistically significantly changed at the given age. (B) Hierarchically clustered heat maps of cell type specific gene changes throughout diabetes. Values were adjusted to a per row color scale so all changes were relative to 2 wk NOD mice.

Figure 5. Quantitative RT-PCR validation of microarray data.

(A) Quantitative RT-PCR was performed using SYBR green detection for the indicated genes. Bars show the mean (log₂) +/− S.E.M. of at least three independent experimental replicates from 3–6 biological replicates per group. All data is represented relative to the expression of actin (ΔC_t). In order to facilitate visualization on a log₂ scale, values were transformed as indicated on the y-axis label. (B) Microarray results for the same genes interrogated in (A). (C) Taqman qPCR quantification of pan-IFNα, Ifnb1, or Ifng throughout diabetogenesis. Bars represent the mean of the normalized probe intensity +/− S.E.M. of 3–6 biological replicates per group. Asterisks indicate statistical significance (P<0.05) from 2 wk NOD sample by one-tailed Mann-Whitney test.

Figure 6. Data modeling of transcriptional networks.

The most significant interaction networks were calculated at (A) 6 wks (FDR 0.00006), (B) 8 wks (FDR = 0.0000001) and (C) 18 wks of age (FDR = 0.00012). Interaction significance was based on the p values calculated for each age group. Notably, at 6-wk, strong developmental signature persists along with AP-1 module, generally attributed to myeloid cells. At 8-wks, strongest interacting subnetwork is T-cell specific, while 18-wk additional B-cell and cytoskeleton specific modules appear.