An Early Islet Transcriptional Signature Is Associated With Local Inflammation in Autoimmune Diabetes Free

Type 1 diabetes (T1D) is driven by the destruction of insulin-producing β-cells in pancreatic islets and affects the lives of millions worldwide. Both genetic and environmental factors, such as virus infection, contribute to disease onset. Viruses, after all, may drive the production of cytokines and chemokines that recruit immune cells such as autoreactive T cells to islets in a process called “insulitis” that results in β-cell destruction. Much of our understanding of inflammatory factors involved in the development of insulitis derives from assessments of peripheral blood, which reveal an increased in immune response genes, including “virus exposure signature genes” in the period preceding T1D onset (1–3). Additional information has been obtained by studying pancreas material from deceased donors with T1D, in which β-cell mass is reduced and inflammatory cells, including CD8⁺ T cells, CD20⁺ B cells, and macrophages, are preferentially found in islets with insulin immunoreactivity (4–7). Human islet heterogeneity is common, as some islets from new-onset T1D appear free of insulitis and insulin remains present.

Examining and defining the processes associated with the onset of human insulitis remains challenging, given that environmental triggers occur long before clinical diagnosis of T1D. Defining features of insulitis prior to diagnosis is crucial for understanding T1D and for identifying intervention strategies. Samples from autoantibody-positive individuals with prediabetes (8) have become increasingly available through resources such the Network for Pancreatic Organ Donors with Diabetes and molecular technologies continue to advance; thus, our understanding of the development of human insulitis continues to expand (9–13).

Animal models are valuable in evaluating the immunopathological mechanisms driving T1D. Such models allow for the direct study of the pancreas and to define the temporal development of insulitis before the onset of diabetes. The NOD mouse is often used for T1D investigation, but mice exhibit peri-insulitis and immune cell infiltration that is distinct from that found in humans with T1D (14,15). Insulitis and autoimmune diabetes can be induced in LEW.1WR1 rats, which, like humans, have a strong genetic MHC class II association with diabetes, with rat RT1-B homologous to human HLA-DQ and RT1-D to HLA-DR. LEW.1WR1 rats develop spontaneous diabetes at a low frequency (˜2%) (16), but infection of weanling rats with the parvovirus Kilham rat virus (KRV) triggers autoimmune diabetes at a rate of 30–40% (17). Furthermore, autoimmune diabetes can be induced with high penetrance (∼90%) in a reproducible time course in weanling LEW.1WR1 rats if KRV infection is preceded by polyinosinic-polycytidylic acid (pIC) treatment on days −3, −2, and −1 (18). This occurs at equal incidence in males and females. Notably, we and others (19–21) have found that KRV does not directly infect pancreatic β-cells, but rather infects splenic cells and regional lymph nodes.

We analyzed transcriptome profiles of islets over a time course to capture temporal dynamics of the prediabetic stages of disease using single-cell RNA-sequencing (scRNA-Seq) to define changes in cells within complex, heterogeneous pancreatic islets during the development of insulitis and before the onset of overt diabetes. For scRNA-Seq, we harvested and dissociated islets from pIC plus KRV–treated LEW.1WR1 rats at three time points when insulitis is detected before diabetes development and compared them to islets from control LEW.1WR1 rats. We identified the emergence of a specific transcriptional signature in subsets of β- and α-cells and determined that such signature cells increase in frequency over time. We also detected a temporal increase in inflammatory markers within mononuclear phagocytic cells. We identified T cells in islets only following pIC plus KRV treatment. We performed an RNA-in situ hybridization (RNA-ISH) on pancreatic sections of pIC plus KRV–treated LEW.1WR1 rats at these time points to spatially map signature β- and α-cells and to correlate with local islet T cell recruitment.

LEW.1WR1 rats were bred at the University of Massachusetts Chan Medical School (UMass Chan). Animals were housed in a viral antibody–free facility in accordance with the UMass Chan Institutional Animal Care and Use Committee, which approved all animal work in this study. Diabetes induction was performed by pretreating rats (either sex, 21–24 days old) with intraperitoneal (i.p.) injection of high-molecular-weight pIC (InvivoGen, San Diego, CA) at 1–2 μg/g body weight on 3 consecutive days (days −3, −2, and −1) followed by a single i.p. dose of 1 × 10⁷ plaque-forming units of KRV on day 0. Control rats received either an i.p. injection of PBS or pIC on days −3, −2, and −1 without the virus. Diabetes was defined when daily blood glucose measurements were >250 mg/dL (hyperglycemia) for 2 consecutive days.

Islets were handpicked from at least two animals per group on days 7, 9, and 11 to capture islets in the prediabetic stages. Islets were placed in PBS, washed once with PBS and once with PBS/10 mmol/L EDTA, then dissociated in TrypLE (Invitrogen, Waltham, MA) for 12 min, followed by trituration until ˜80% dissociated. Digestion was stopped with 2% BSA in PBS, and islets were filtered once with a cell strainer (Falcon #352235).

Samples were processed for 10X Genomics per the manufacturer’s protocol. The cDNA libraries were prepared using the 10X Library prep kit (Next GEM Single Cell 3′ Library Kit v3) and sequenced with the Illumina NextSeq 500/550.

The paired-end fastq files contain the 16-base cell barcode and the12-base unique molecular identifier (UMI) in the R1 read and the 55-base transcript mRNA sequence in the R2 read. These paired-end reads were preprocessed using a custom Python script to extract the cell barcode and UMI from each R1 read and append them as a colon-delimited pair to the corresponding R2 read name. Reads with cell barcodes that were not in the 10X cell barcode dictionary, or with ‘N’s in either the cell barcode or UMI were discarded. The resulting fastq files were then processed through our DolphinNext analysis pipeline (22) as single-ended reads, removing reads from any cell barcode with <500 reads. The mRNA sequences were aligned to the rat genome (rn6) using TopHat2 (v2.0.12) with default settings and RefSeq transcript annotations. Gene transcripts were quantified using the End Sequencing Analysis Toolkit (ESAT; https://github.com/garber-lab/ESAT), again using RefSeq transcript annotations. ESAT ignores reads that result from PCR duplication during the library preparation process using the UMI. If reads from the same cell barcode map to the same gene and have the same UMI, they are considered PCR duplicates, and only one is counted. The output of ESAT is an array containing the observed transcript counts for each gene for each cell.

Cell type identification was a multistep process using custom R scripts (R v3.5.0) based on the Signaling Single Cell package (https://github.com/garber-lab/SignallingSingleCell). All functions referenced refer to that package. The data for each sample were first processed to remove any cells with <1,000 transcripts. All samples were combined into a single ExpressionSet data object for further analysis (R/Bioconductor Biobase package). The 2,000 genes were selected using the variance-stabilizing transformation gene selection method to reduce noise introduced by low variance and low expression genes, followed by t-distributed stochastic neighbor embedding (t-SNE) mapping [dim_reduce() with default parameters] and density clustering [cluster_sc() with method = ‘density’ and num_clust = 12] using the selected genes to produce an initial segmentation of the cells. We then used the cluster and mapping boundaries to classify the cells into major groups using the expression levels of known marker genes to identify the groups (e.g., Ins1 for β-cells, Sst for ∂-cells). The expression of each endocrine marker gene had a bimodal distribution, so the cell type–specific expression threshold was set as the local minimum. Because endocrine cell types (α, β, δ, pancreatic polypeptide [PP], and ε) accounted for 75% of the total cells, we focused on removing endocrine–endocrine doublets (e.g., α/β, β/δ, etc.) from further analysis. Endocrine–endocrine doublets were defined as cells that had high expression of marker genes from more than one endocrine cell type.

For each identified cell type, we selected the cells from each condition (day 7, 9, or 11 post–pIC plus KRV treatment versus PBS control or pIC control) and performed differential expression (DE) analysis (edgeR) to identify genes that were most differentially expressed.

Pathway enrichment analysis was performed on differentially expressed genes between disease states and visualized as a list of enriched gene ontology (GO) terms with Benjamini-Hochberg adjustment of P values for multiple hypothesis testing. GO enrichment was performed with the online tool DAVID (https://david.ncifcrf.gov/tools.jsp) using lists of significantly upregulated or downregulated genes using the Rattus norvegicus background and the “GOTERM_BP_DIRECT” biological process database. GO analysis was performed for α-cells, β-cells, and macrophages.

At designated time points, pancreata were harvested from euthanized rats, fixed immediately in 10% buffered formalin, and then embedded in paraffin. Paraffin-embedded pancreas was sectioned and prepared for light microscopy in the UMass Chan Morphology Core laboratory (https://www.umassmed.edu/morphology/protocols). Sections were stained with hematoxylin-eosin (H-E) and scored for insulitis as previously described (23): 0, no inflammatory mononuclear cell infiltration; 1+, small numbers of infiltrating cells with preservation of islet architecture; 2+, moderate infiltrating cells with preservation of architecture; 3+, many cells with most islets affected and distortion of islet architecture; and 4+, florid infiltration and distorted islet architecture or end-stage islets with or without residual inflammation.

The RNAscope Multiplex fluorescence V2 Kit was used according to the manufacturer’s instructions (Advanced Cell Diagnostics, Newark, CA) on formalin-fixed, paraffin-embedded tissues. Briefly, fixed cells were first pretreated by rehydrating with ethanol and 1× PBS wash, followed by treatment with RNAscope hydrogen peroxide and protease III. Formalin-fixed, paraffin-embedded tissue sections were baked at 60°C for 1 h. After deparaffinization with xylene and hydration with ethanol, tissue sections were treated with hydrogen peroxide, then heated in antigen retrieval buffer, followed by digestion with proteinase. Target RNA probes for rats were available from Advanced Cell Diagnostics as follows: Ins1 (no. 413411), Gcg (no. 315471), and Cxcl10 (no. 553191). Probes were hybridized to pretreated sections for 2 h at 40°C, followed by a series of signal amplification with the kit-provided Preamplifier and Amplifier. Fluorescence signal was generated using Opal 520, 570, and 690 dyes (Akoya Biosciences, Marlborough, MA). Sections or fixed cells were counterstained with kit provided DAPI, mounted, and stored at 4°C until image acquisition. For select samples, RNAscope was combined with immunofluorescence using mouse anti-rat CD8a antibody OX-8 at 1:100 (no. 201701; BioLegend, San Diego, CA) and goat anti-mouse IgG Alexa Fluor 647 secondary antibody at 1:1,000 (no. A21235; Invitrogen) according to the protocol specified by Advanced Cell Diagnostics.

Imaging was performed using a Nikon Eclipse Ti series microscope (×20 objective), Leica SP8 Laser Scanning Confocal Fluorescence microscope (×40 objective), or the Sanderson Center for Optical Experimentation facility TissueFAXS SL tissue cytometer (×20 and ×40 objectives).

Scanned tissue sections were visualized using TissueFAXS SL Viewer version 7.1 software. Automated counts were performed in the StrataQuest analysis package, version 7.1. Only islets formed by ≥30 cells were included to avoid possible detection errors derived from small artifacts. To confirm that findings between manual and automated analyses were consistent, the numbers of Ins1⁺ and Gcg⁺ islets, as well as the presence of Cxcl10 and CD8, were assessed manually. Manual analysis was performed by researchers blinded to the study groups. Two different operators evaluated the expression of every marker and islet independently and blinded to the experimental conditions. A third operator combined and revised the analysis.

Comparisons between independent groups were performed using one-way ANOVA with multiple comparisons or simple linear regression analysis. The significance level was 0.05 for all statistical tests. Statistical analyses were performed using GraphPad Prism version 9.

The scRNA-Seq raw fastq files, gene by cell transcript counts matrix, and metadata file, including the viral transcript count for each cell, are deposited in the National Center for Biotechnology Information Gene Expression Omnibus database under accession number GSE203412.

Weanling LEW.1WR1 rats develop progressive insulitis following pIC plus KRV administration, whereas those receiving either pIC or PBS control remain insulitis free (Fig. 1A). Diabetes occurs in ∼100% of rats by day 25 post–pIC plus KRV treatment (Fig. 1B). With the aim of capturing cellular transcripts during the prediabetic state, pancreatic islets were isolated from euthanized animals over a time course on days 7, 9, and 11 following pIC plus KRV treatment (KRV7, KRV9, and KRV11) during the development of insulitis in the prediabetic stage (Fig. 1C). Islets were collected at day 11 following administration of pIC alone (pIC11) or PBS alone (PBS11) as control samples. Dissociated islet cells were processed for scRNA-Seq (10X Genomics), and each library was sequenced to a depth of 20,000 to 25,000 reads/cell with an average of 12,174 transcripts/cell and 1,452 nonzero genes/cell.

We used a combination of clustering and marker gene expression to determine major cell groups. Endocrine cells were identified by high expression of specific markers as follows: α, Gcg; β, Ins1; δ, Sst; PP, Ppy; and ε, Ghrl. Endocrine doublets had high expression of multiple endocrine cell markers. Nonendocrine cells had high expression of the following markers: T cell, Cd3g/d/e; macrophage, Cd68; stellate, Pdgfrb, Col1a1, and Acta2; vascular, Pecam1, Esm1, and Flt1; red blood cell (RBC), Hba1 and Hbb; acinar, Cpa1; ductal, Krt19; and natural killer (NK) cells, Nkg7 (see Supplementary Fig. 1 for t-SNE mapping of select markers). Insignificant numbers of acinar and ductal cells were identified in any of the samples. Comparisons between pIC plus KRV–treated rat islets over time (KRV7, KRV9, and KRV11) were performed to assess temporal transcriptional changes associated with the development of insulitis. Comparisons were also made between day 11 pIC plus KRV–treated rat islets and control samples (KRV11, PBS11, and pIC11). The distribution of cells from all 5 times/conditions is shown in Fig. 1D. The distribution of specific cell types for all treatment times/conditions is displayed in Fig. 1E. Distinct cell types are represented based on cluster analysis, including endocrine (β, α, δ, PP, and ε), immune (T lymphocytes and macrophages), and other cell populations (vascular, stellate, and RBC). MHC class II genes RT1-Ba, Bb, Da, Db1, Db2, DOa, and Dob are expressed only in cells in the macrophage cluster, indicating that macrophages are the predominant antigen-presenting cells (Supplementary Fig. 1). While the percentage of each endocrine cell type remains constant between times/conditions, T cells are notably absent (i.e., 0% frequency) in PBS11 and pIC11 control samples, then emerge at increasing frequency following pIC plus KRV treatment, comprising 0.3%, 0.8%, and 2.3% of all cells in KRV7, KRV9, and KRV11 samples, respectively (Fig. 1F). Absolute cell counts for each cell type, as well as the percentage of all cells for each cell type for the five treatment times/conditions, are provided in the worksheet in the Supplemental Material. Thus, scRNA-Seq analyses of islets at defined time points revealed the specific cell types and transcripts engaged in the development of insulitis before diabetes onset.

Next, we examined the DE of transcripts from endocrine cells. For the initial analysis (i.e., examining the changes during the development of insulitis), data from endocrine cells were extracted from KRV7, KRV9, and KRV11 and then combined, mapped, and clustered. To analyze the differences between KRV11, PBS11, and pIC11 samples, data from endocrine cells were separately extracted and then combined, mapped, and clustered. We found select transcripts to be highly expressed in separate subclusters of β-cells as well as subclusters of α-cells from KRV7, KRV9, and KRV11 samples. To distinguish cells in these subclusters from other cells, we refer to them as inflamed or “hot.” Hot β-cells increase in frequency over time, comprising 0.075% of all β-cells in KRV7, 1.0% of all β-cells in KRV9, and 15.2% of all β-cells in KRV11. Hot β-cells are absent in pIC11 and PBS11 control samples. Complete DE results comparing different times/conditions are provided in the worksheet in the Supplementary Material.

To further analyze hot subclusters, we performed DE analysis on β-cells from the KRV11 sample, comparing hot β-cells versus all other β-cells (considered “cool”). Cells in hot subclusters express high levels of transcripts for the defense-response chemokine Cxcl10, diubiquitin (Ubd), and the long noncoding RNA (lncRNA) LOC100910973. Comparisons between hot β-cells from KRV11 versus all β-cells from either pIC11 or PBS11 control samples yield similar results. Transcriptional profiles in cool β-cells are akin to those from pIC11 or PBS11 control samples (see worksheet in Supplementary Material). Because Cxcl10, Ubd, and LOC100910973 are the most highly induced transcripts in hot β-cells, their expression in all cell types is displayed (Fig. 1G). t-SNE maps for β-cell expression of these transcripts over the time course are shown in Fig. 2A.

DE of transcripts for hot versus cool KRV11 β-cell comparisons are shown by rank in Fig. 2B, along with GO terms. They fall into three main categories, including interferon (IFN)–stimulated genes (ISGs), MHC class I–associated genes, and ubiquitin-proteasome pathway genes. ISGs include Stat1, a transcription factor essential in the IFN response pathway, as well as downstream genes such as Isg15, Mx2, and Oasl2 and the chemotactic factors Cxcl9, Cxcl10, and Cxcl11 (24). RT1-CE10 (MHC class 1b) and RT1-S3, the functional analog of HLA-E, are both increased, as are MHC class I–related transport proteins B2m, Tap1, and Tapbp. Indeed, MHC class I overexpression has been shown to be mediated by type I IFN in human islets (25) and is a feature of T1D (26). The MHC class II–associated gene CD74 is also increased. Specific genes in the ubiquitin-proteasome system are also upregulated in β-cells from prediabetic rats from day 7 onwards. Ubd expression is highly induced in pIC plus KRV–treated rats compared with controls. Ubd encodes the ubiquitin D protein, which attaches covalently to a target protein and induces its degradation through the ubiquitin-proteasome system. Likewise, β-cells from the pIC plus KRV group exhibit significant elevations in transcripts for subunits of the 26S proteasome, including Psmb8, Psmb9, and Psmb10, which is consistent with studies that show that PSMBs can be induced by type I IFN (27,28). Another highly induced transcript is variant 2 of the lncRNA LOC100910973 (National Center for Biotechnology Information reference NR_102352.1). Given that this lncRNA has not been characterized previously, we sought to determine if it is an ISG and stimulated normal rat kidney cells with recombinant rat IFN-β and measured LOC100910973 expression by quantitative RT-PCR. We found that LOC100910973 was induced by IFN-β compared with control cells, much like Isg15 under these same conditions (Supplementary Fig. 2).

Next, we wanted to exclude the possibility that hot β-cells represent β-cells that have undergone endocytosis by professional phagocytes of the immune system. We confirmed that hot β-cells cluster with other β-cells, express high levels of Ins1, Ins2, and Iapp, and do not express markers for phagocytic mononuclear cells. A macrophage engulfing a β-cell might be expected to express both Cd68 and Ins1, yet Cd68 is not expressed in any of the cells labeled β (Supplementary Fig. 3).

Likewise, a distinct subcluster of hot α-cells also emerges; such cells comprise 0.15% of all α-cells in KRV7, 0.06% of all α-cells in KRV9, and 20.2% of all α-cells in KRV11 but are absent in pIC11 and PBS11. As for β-cells, Cxcl10, Ubd, and LOC100910973 are the three most highly increased transcripts, and t-SNE maps for these transcripts are shown in Fig. 2C. However, fewer hot transcripts are induced by at least fourfold compared with cool for α-cells than for β-cells. The transcripts with the greatest increase by differential expression (DE) are shown by rank in Fig. 2D, along with associated GO terms. All transcripts highly induced in hot α-cells are also induced in hot β-cells and therefore are not unique to α-cells. Comparisons between all groups and a complete list of all DE transcripts are provided in the worksheet in the Supplementary Material.

To further examine islet macrophage transcripts, cells labeled as macrophages (i.e., Cd68⁺) were extracted from the data set, and a separate t-SNE map was generated. Notably present is a separate subcluster of macrophages, mostly comprised of cells from KRV9 and KRV11 samples, with high expression of Cxcl10, Ubd, and LOC100910973 (Fig. 3A, black circles in top panels). DE analysis was performed for the hot macrophage group compared with all other macrophages (minus macrophage–endocrine cell doublets, dendritic cells, and B cells; see below) regardless of condition/day, and a ranked list and associated GO terms are shown in Fig. 3B. Along with Cxcl10, Ubd, and ISGs, a macrophage activation program appears to be induced in the hot macrophage group, given that the transcripts Cd40 and Ly6c are greatly increased (Fig. 3A, black circles in middle panels). Furthermore, increases in Socs1 encoding suppression of cytokine signaling as well as Sod2 encoding superoxide dismutase 2, a marker of oxidative stress, are prominent. Altogether, the data support activation of a proinflammatory state in prediabetic rat islet myeloid cells.

In addition, three small groups of cells in the top left quadrant (Fig. 3A, green circles in bottom panels) likely represent endocrine/macrophage doublets on dual Cd68 plus Ins1, Gcg, or Ppy expression. These may be macrophages that have engulfed β-, α-, or PP cells, respectively; alternatively, they are macrophage/endocrine cell doublets that persist following single-cell processing. A small cluster in the top right quadrant (Fig. 3A, blue circle in middle panels) appears to be dendritic cells based on high expression of Batf3, Irf8, and Flt3 (29) and are primarily found in KRV9 and KRV11 samples. A final cluster in the lower right quadrant is Cd68, and expresses B-cell markers, including Cd79b (Fig. 3A, orange circle in middle panels). DE analysis for dendritic cells and B cells is provided in the worksheet in the Supplementary Material.

As noted above, T cells are only detected in islets collected from pIC plus KRV–treated rats, with T cells increasing in frequency from day 7 to 11 (0.3%, 0.8%, and 2.3% of all cells from KRV7, KRV9, and KRV11, respectively), similar to the increased frequency of hot β- and α-cells between days 7 and 11. These T cells are Cd8⁺ and Cd4⁻, as shown in Fig. 4A. DE analysis could only be performed on pIC plus KRV–treated rat islets with T cells since T cells are absent in PBS and pIC control rat islets; increased transcripts are listed in Fig. 4B. While no significantly enriched GO pathways are found when comparing days 11 versus 9, 11 versus 7, or 9 versus 7, expression of several immune-related genes, including Ifng, Gzmb, and Icos, is modestly increased in T cells at day 11 compared with day 9.

Given prominent Cxcl9/10/11 expression in hot β- and α-cells, we sought to define the spatial distribution of such cells in the islets. Pancreas sections from pIC plus KRV–treated LEW.1WR1 rats were stained with RNAscope probes targeting Cxcl10, Gcg, and Ins1. Cxcl10 expression is found in both Ins1⁺ and Gcg⁺ cells within Cxcl10⁺ islets (representative images, Fig. 5A). The Cxcl10⁺ and Cxcl10⁻ islet distribution appears to be random, as shown in Fig. 5B. H-E staining of the section is shown for reference (Fig. 5C). A few Cxcl10⁺ islets have decreased Ins1 expression relative to Gcg, suggesting that expression may be associated with the specific loss of β-cells (Supplementary Fig. 4), but in the LEW.1WR1 rat model, nearly all islets are eventually destroyed following pIC plus KRV treatment. Quantification of the frequency of Cxcl10⁺ islets in pancreatic sections confirmed a significant increase from day 9 to 11 following pIC plus KRV treatment and absence in PBS controls (Fig. 5F). Scanned high-resolution images of several representative rat pancreatic sections, stained with either H-E or RNAscope probes, are available in Supplementary Fig. 5.

A primary function of chemotaxis is to attract immune cells with the appropriate receptors, so we analyzed T cells in the pancreatic islets in more detail. The CXCL10-CXCR3 axis contributes to inflammation and pathogenesis of several autoimmune diseases, including human T1D (30–32). Given the high induction of Cxcl10 in hot β- and α-cells following pIC plus KRV treatment, we sought to establish if Cxcl10 expression is associated with the recruitment of CD8⁺ T cells that bear CXCR3. To determine if islets with high Cxcl10 expression are associated with CD8⁺ T cell infiltration, we stained pancreata for CD8 using an anti-rat CD8 antibody in combination with Cxcl10 and Ins1 RNAscope probes. Abundant CD8⁺ cells are recruited to islets over time following pIC plus KRV treatment and are associated with Cxcl10⁺ islets, whereas few CD8⁺ cells are associated with Cxcl10⁻ islets (Fig. 5D, E, G, and H). We confirmed the presence of Cxcr3 transcript in CD8a-expressing cells (Fig. 4C) to support the paired receptor–ligand interaction between CXCR3 and CXCL10 during insulitis. Of note, increased Cxcl10 expression in the pancreas is present in non-β- and non-α-cells, which is expected given increased expression of Cxcl10 in nonendocrine cells, such as macrophages in scRNA-Seq analyses. Given that preproinsulin-specific CD8⁺ T cells are found in healthy human exocrine pancreas (33), we inspected whole pancreatic sections from PBS-treated rats to determine if such cells are present in this model. CD8⁺ cells are abundant in pIC plus KRV–treated rat pancreas but are not found in either endocrine or exocrine cells of PBS-treated animals (Supplementary Fig. 6).

We successfully obtained high-quality islet single-cell transcriptomes to explore signature profiles associated with the development of autoimmune diabetes. Differential gene expression analysis of prediabetic rat islets identified β- and α-cell subpopulations with a highly reproducible signature transcriptional profile comprised of ISGs, including the CXCL9–11/CXCR3 axis, MHC class I–associated genes, and ubiquitin-proteasome system genes. RNA-ISH and immunostaining of pancreatic tissues provided for visualization of the heterogeneity of Cxcl10 expression in islets and expression within β- and α-cells in individual islets and defined interactions between islet cells and immune cell recruitment. Using RNA-ISH and immunostaining analyses, we confirmed that β-cells and α-cells that highly express Cxcl10 colocalize within specific islets, sounding the alarm for CD8⁺ T cell recruitment around and within Cxcl10⁺ islets and further escalating the inflammatory process. Together, our data indicate that a subset of type I IFN–driven β-cells undergo a series of transcriptional changes that make them prone to attack by immune cells.

The crucial role of type I IFNs in autoimmune diabetes pathogenesis is evident from our previous study showing that type I IFN receptor knockout rats are protected from diabetes (23). Genome-wide association studies identify polymorphisms in IFN pathway genes to be associated with T1D. Both clinical and epidemiological evidence suggest a link between viral infection and human T1D pathogenesis (1). IFNs induce proteasome modification, with proteasomes becoming immunoproteasomes that process MHC class I peptide (34). Notably, in previous work, Ubd was found to be a susceptibility gene for virus-triggered autoimmune diabetes in LEW.1WR1 rats, and UBD-deficient rats have substantially reduced diabetes after viral infection (35). Our findings from the rat model are remarkably consistent with reports that HLA class I is highly upregulated in human T1D islets (36), which has been confirmed in β-cells from T1D donors in our own scRNA-Seq studies (37). Type I IFN may drive the signature transcript in β- and α-cells, but the timing, identification, and localization of IFN-producing cells in the rat model remain undefined. Defining the exact nature of the trigger for islet inflammation will be important for preventing the subsequent autoimmune processes targeting pancreatic islets.

CXCL10 plays a dominant role in the attraction of effector T cells bearing CXCR3 in islets from donors with T1D, and it remains an active candidate for therapeutic intervention to prevent T1D (38). The CXCL10-CXCR3 axis has been studied extensively for its role in T1D. Serum CXCL10 is found at higher levels from newly diagnosed T1D donors compared with controls (39,40) as well as in long-standing T1D (41), and both CXCL10 and CXCR3 are present in the microenvironment of T1D islets (42–44). Furthermore, animal models of T1D, including NOD mice (42,45) and RIP-LCMV mice (46,47), feature increases in islet CXCL10. Our study demonstrates that Cxcl10 is highly upregulated in subsets of β- and α-cells from pIC plus KRV–treated rat islets; such islets attract chemokine receptor–expressing immune cells during autoimmunity. As we and others previously demonstrated (19,21), this occurs in the absence of KRV infection in pancreatic β-cells and is likely an effect of an IFN-rich environment. The production of IFN by lymphoid tissue (spleen and pancreatic lymph nodes) may be the trigger for islet inflammatory signatures and remains an area of active investigation. Such studies coupled with network interaction analysis of the rat model could provide novel insights on factors crucial for autoimmunity in human T1D.

We found increased expression of several additional transcripts in endocrine cells that may participate in T1D pathogenesis. PDL1 and HLA-E have been identified in insulin-containing islets from patients with T1D (48,49). PD-L1 is a ligand of PD-1, a protein expressed in T cells, and PD-L1–PD-1 binding inhibits T cell function; inhibition of PD-1–PD-L1 in the context of immune checkpoint therapy for cancer triggers immune-mediated endocrine diseases, including T1D (50,51). HLA-E allows for immune surveillance and inhibition of NK cells (52), and selective retention of HLA-E protects stem-cell–derived islets from immune rejection (53). Notably, transcripts for both PD-L1 (Cd274) and the functional analog of HLA-E (RT1-S3) are increased in hot rat β-cells, likely because of type I IFN induction. Nonetheless, the specific impact of their induction on immune cell recruitment and disease can be further studied in rats through their blockade.

lncRNAs in immune cells and in β-cells are likely to contribute to autoimmune pathology in human T1D (54). LOC100910973 is highly induced in hot β-cells and α-cells across all prediabetic rat samples, and our preliminary analyses indicate that it is an ISG. lncRNA editing could be performed in cultured rat cells and in the LEW.1WR1 rat model of diabetes to establish if and how LOC100910973 participates in pathogenesis. If knocking out LOC100910973 protects against the development of autoimmunity, it would make a strong case for identifying and targeting lncRNAs that may be associated with human T1D.

We previously reported that both MHC class II invariant chain CD74 and CIITA, the transcriptional regulator controlling class II gene expression, are upregulated in β-cells in T1D islets (37). In the rat islets, Cd74 transcripts are significantly induced in hot β- and α-cells. However, unlike in human samples, Ciita is not induced in such β-cells. Using immunostaining in future studies, we can confirm whether MHC class II protein is present or absent in β-cells in the LEW.1WR1 rat model, including at time points beyond day 11 and at the onset of full-blown diabetes, to see if we can capture such findings. The rat model could then be further used to dissect the role of β-cell MHC class II in recruiting islet-infiltrating CD4 T cells.

Following the development of islet inflammatory signatures, both rat islet α- and β-cells are ultimately destroyed. The spatial distribution of Cxcl10 expression is not confined to specific regions of the pancreas but appears scattered across the entire pancreas. This is distinct from human T1D, in which α-cell numbers remain preserved or even increase (55). During the prediabetic state, human β-cells are selectively destroyed via autoimmune mechanisms, and the α-cell component begins to expand (56). Immunohistochemical analyses of pancreata of organ donors with T1D reveal that insulitis, the hallmark of autoimmunity, evolves over time with varying degrees of insulin positivity within donors (57). Interestingly, both α- and β-cells can have elevated CXCL10 expression in cases of T1D (58), and a subset of transcriptionally activated islets are found in humans before the onset of T1D, with some markers being present in both α- and β-cells in Network for Pancreatic Organ Donors with Diabetes–obtained samples (11). Rat islets are characterized by the localization of β-cells predominantly in the central islet, and α-cells are in the periphery (59); this is in sharp contrast to human islets, in which α-, β-, and δ-cells are distributed throughout islets. Additional posttranscriptional or donor-specific factors may contribute to species-level differential findings between rat and human, and selective destruction of β- but not α-cells in human T1D is likely mediated by multiple factors. Whether heightened Cxcl10 expression might be evolutionarily beneficial or detrimental for β- and α-cells has not been determined.

The finding of activated macrophages in islets from prediabetic rats is consistent with that reported in mouse diabetes studies (60). The activation state of macrophages may orchestrate T-cell recruitment and pathogenesis. Since the number of myeloid cells in the pooled islets is somewhat low, separately isolating such cells from samples to obtain larger quantities would permit a more powerful transcriptional analysis. We are further exploring the role of inflammasome activation and interleukin-1 effects on islet function in the rat model of autoimmune diabetes. The role of inflammasomes in T1D is of high interest (see Pearson et al. [61] for review) and has preventive and therapeutic implications. Likewise, analyzing larger pools of T cells would be helpful, given the low yields in these studies.

Precise identification of factors driving immune cell activation and recruitment in the rat pancreas will aid in the long-term goal of discovering critical elements in the autoimmune destruction of β-cells. Assessing interspecies differences between rat and human autoimmune diabetes with comprehensive analyses and comparisons of the diabetic and prediabetic states can provide new mechanistic insights. As findings emerge from the study of human samples, including spatial transcriptomics of islets from pre-T1D cases, the rat data provided in this study can be mined for correlations. The rat model can be exploited through targeted editing via CRISPR technologies to understand the significance of a specific transcript on pathogenesis, which will assist in future therapeutic and prevention strategies for human T1D. In sum, transcripts from rat islets prior to onset of diabetes are markedly similar to those in islets from humans with T1D. Future investigations will elucidate when and how these might impact the development of autoimmune diabetes.

Funding. This study was supported by the National Institutes of Health/the National Institute of Allergy and Infectious Diseases grant R01 AI139095 (to J.P.W.) and Massachusetts Life Sciences Center Bits to Bytes award (to C.E.B. for SCOPE TissueFAXS SL).

Duality of Interest. A.G.D. owns shares of 10X Genomics and Illumina common stock. No other potential conflicts of interest relevant to this article were reported.

Author Contributions. M.G. and J.P.W. conceived and designed the research. A.G.D., A.A., B.S., S.D.R., N.Q., Z.G., and E.V. performed experiments and analyzed and interpreted data. M.I.T. reviewed and edited the manuscript. C.E.B., D.M.H., D.L.G., and M.G. provided scientific input and reviewed the manuscript. J.P.W. wrote the manuscript. J.P.W. is the guarantor of this work and, as such, had full access to all of the data in the study and takes responsibility for the integrity of the data and the accuracy of the data analysis.