1Department of Surgery, The University of British Columbia, and Child & Family Research Institute, Vancouver, British Columbia, Canada

2Department of Medicine and Centre for Heart Lung Innovation, The University of British Columbia, and Prevention of Organ Failure (PROOF) Centre of Excellence, St. Paul’s Hospital, Vancouver, British Columbia, Canada

Abstract

Type 1 diabetes (T1D) is caused by immune-mediated destruction of insulin-producing β-cells. Insufficient control of autoreactive T cells by regulatory T cells (Tregs) is believed to contribute to disease pathogenesis, but changes in Treg function are difficult to quantify because of the lack of Treg-exclusive markers in humans and the complexity of functional experiments. We established a new way to track Tregs by using a gene signature that discriminates between Tregs and conventional T cells regardless of their activation states. The resulting 31-gene panel was validated with the NanoString nCounter platform and then measured in sorted CD4+CD25hiCD127lo Tregs from children with T1D and age-matched control subjects. By using biomarker discovery analysis, we found that expression of a combination of six genes, including TNFRSF1B (CD120b) and FOXP3, was significantly different between Tregs from subjects with new-onset T1D and control subjects, resulting in a sensitive (mean ± SD 0.86 ± 0.14) and specific (0.78 ± 0.18) biomarker algorithm. Thus, although the proportion of Tregs in peripheral blood is similar between children with T1D and control subjects, significant changes in gene expression can be detected early in disease process. These findings provide new insight into the mechanisms underlying the failure to control autoimmunity in T1D and might lead to a biomarker test to monitor Tregs throughout disease progression.

Introduction

The precise mechanisms underlying the failure to control autoreactive T-cell expansion in type 1 diabetes (T1D) are unclear, but evidence suggests that changes in FOXP3+ regulatory T cells (Tregs) may be involved. People lacking functional Tregs due to genetic mutation are born with T1D (1), and in NOD mice, Treg deletion accelerates T1D progression, whereas transfer of Tregs prevents disease (2). In patients with T1D, Treg proportions in blood appear to be normal (3), although accurate enumeration is difficult because activated conventional T cells (Tconvs) also express FOXP3 and CD25. Some evidence points to impaired Treg function in T1D (4–7), possibly related to differentiation into cells producing interleukin (IL)-17 (6) or interferon-γ (7), but in vitro suppression assays may not reflect in vivo function because cells are removed from their inflammatory environment. T1D Tregs may also have diminished IL-2–mediated activation of STAT5 (8), resulting in destabilized FOXP3 expression and diminished Treg stability.

Although many clinical trials in T1D have aimed to enhance Tregs, the lack of a simple method to accurately monitor Tregs in blood is a major limitation that hinders identification of subjects with dysfunctional Tregs and, thus, the tracking of changes in Tregs throughout disease progression. We report a new gene signature of human Tregs created by comparing the gene expression profiles of memory and naive Tregs with those of memory and naive Tconvs ex vivo and after in vitro activation. Measurement of the Treg signature in T1D revealed a subset of genes that accurately classifies children with new-onset T1D versus healthy children. These data provide new evidence that T1D Tregs have intrinsic changes and provide a clinically applicable way to monitor for changes in Tregs in T1D.

Research Design and Methods

Peripheral blood from 12 adult and 24 pediatric control subjects, 29 patients with new-onset T1D, and 27 patients with established T1D was obtained upon written consent and age-appropriate assent as approved by The University of British Columbia Clinical Research Ethics Board. Samples were collected on average 1.5 (range 0.9–3.4) and 28.9 (9.4–50.2) months after diagnosis of new-onset and established T1D, respectively.

Affymetrix Gene Expression Analysis and Gene Signature Determination

mRNA was extracted, reverse transcribed, and amplified, and cDNA was hybridized to Affymetrix human genome U133A 2.0 microarrays (GEO accession GSE76598). After preprocessing, normalization, and quality control, four samples were excluded. The data were analyzed with the Marker Selection module (GENE-E microarray analysis software) with a false discovery rate cutoff of <25%. Two rounds of Gene Set Enrichment Analysis were performed with and without a cutoff (750-count difference between Tregs and Tconvs) and including publicly available data sets of mouse (9,10) or human (11–15) Treg/Tconv gene expression profiles.

Hierarchical Clustering

The hierarchical clustering module (public GenePattern server, Broad Institute [18]) was used to generate heat maps and cluster trees by pairwise average linkage according to Pearson correlation coefficient.

Statistical Analysis

Analyses between two groups were performed with one-tailed Student t test or Mann-Whitney t test, and for more than two groups, the one-way ANOVA with Šidák or Bonferroni post hoc test was performed. P < 0.05 was considered significant. Error bars in figures represent the SD.

Results

Development of a Human Treg Gene Signature

To generate Tregs and Tconvs in multiple states of activation, naive (CD4+CD25hiCD45RA+) and memory (CD4+CD25hiCD45RA−) Tregs and naive (CD4+CD25-CD45RA+) and memory (CD4+CD25−CD45RA−) Tconvs were sorted from blood of seven healthy adults (Fig. 1A). RNA was isolated ex vivo or after stimulation for 40 h, promoting activation-induced FOXP3 in Tconvs (Fig. 1B). Tregs, but not Tconvs, were suppressive in vitro (Fig. 1C).

Development and validation of a gene signature that discriminates between Tregs and Tconvs regardless of activation state. A: Sorting strategy. CD4+ T cells from seven healthy adults were enriched as CD25+ or CD25− cells then sorted into four populations: naive Treg (nTreg) (CD4+CD25hiCD45RA+), memory Treg (mTreg) (CD4+CD25hiCD45RA−), naive Tconv (nTconv) (CD4+CD25−CD45RA+), and memory Tconv (mTconv) (CD4+CD25−CD45RA−). B: Expression of FOXP3 and CD25 immediately postsort or after activation with α-CD3/α-CD28–coated beads for 40 h. C: In vitro suppressive activity of sorted T-cell subsets. Proliferation of CD4+ responders (8,000 cells/well) stimulated with α-CD3/α-CD28–coated beads for 6 days in the presence of the indicated ratios of sorted T-cell subsets. [3H]-thymidine was added to the cultures for the final 16 h before harvesting; means with SD are shown. Statistical analysis by two-way ANOVA with Bonferroni multiple comparisons test between Treg and Tconv subsets: *P ≤ 0.05, **P ≤ 0.01. D: mRNA isolated after sorting or after 40-h stimulation with α-CD3/α-CD28–coated beads was analyzed with Affymetrix U133A microarrays. Thirty-one genes were selected (see research design and methods) that effectively discriminated between Tregs and Tconvs regardless of their state of activation by hierarchical clustering. E: Gene expression in replicate samples was measured by NanoString. F: Gene expression in naive, memory Tregs/Tconvs sorted as in A as well as in total Tregs/Tconvs sorted as CD25hiCD127lo or CD25−CD127+, respectively, was measured by NanoString. A and B depict representative data (n = 7); C shows mean with SD (n = 7); and in D–F, each column represents a single donor, with biological replicates of five to seven for each cell type in the Affymetrix data and three to four in NanoString data. Each row is a gene, with the gene symbol on the right. stim, stimulated with α-CD3/α-CD28–coated beads for 40 h; us, unstimulated.

The gene expression profile of the eight cell subsets was analyzed by Marker Selection module and Gene Set Enrichment Analysis. Hierarchical unsupervised clustering revealed that the composite expression of 31 genes distinguished Tregs from Tconvs regardless of activation state (Fig. 1D). This differential gene signature was validated using NanoString (Fig. 1E) and separated Tregs from Tconvs by using three independent Treg sorting strategies: naive CD25hiCD45RA+, memory CD25hiCD45RA−, or total CD25hiCD127lo (Fig. 1F).

Analysis of the Treg Signature in T1D

Because many of the 31 signature genes are associated with Treg function (e.g., FOXP3, CTLA-4), we hypothesized that a change in Treg signature may be detected in diseases associated with Treg dysfunction, such as T1D. Therefore, age-, sex-, and BMI-matched (Fig. 2A) (data not shown) subjects with new-onset (n = 29) or established (n = 27) T1D and pediatric control subjects (n = 24) were recruited. Subjects with T1D had stable glycemic control as indicated by comparable fed-state levels of glycated hemoglobin and daily insulin dose (Fig. 2B).

CD4+CD25hiCD127lo Tregs and CD4+CD25loCD127hi Tconvs were sorted (Fig. 2C), with comparable FOXP3 expression in all three cohorts (Fig. 2C) and demethylated TSDR in Tregs (Fig. 2D). The Treg gene signature expression was measured by NanoString. PCA of the raw log2-transformed data showed that the composite 31-gene signature effectively discriminated between Tregs and Tconvs from pediatric samples regardless of disease state (Fig. 2E and Supplementary Fig. 1). IL1R1 was below background and excluded from further analysis, leaving a 30-gene panel.

A Six-Gene Panel in Tregs Is a Sensitive and Specific Biomarker of New-Onset T1D

The composite 30-gene signature was then compared among subject groups using PCA, and the algorithm performance was tested (α = 0) (Fig. 3 and Supplementary Fig. 2A–C). Control versus new-onset T1D Tregs had biomarker potential, with an AUC of 0.73 ± 0.15 compared with an AUC of 0.59 ± 0.16 for control subjects versus subjects with established T1D and 0.50 ± 0.16 for subjects with new-onset versus established T1D. No significant correlations were found in Tconvs gene expression or the ratio of Tregs/Tconvs (data not shown).

To optimize classification power of a T1D-specific algorithm, the biomarker discovery analysis was repeated with α = 0.65 (median 10 genes/panel, error rate 0.31 ± 0.12). By using this approach, we identified six genes present in >65% of panels (Fig. 4A). Considering each of these six genes individually, three were significantly downregulated (TNFRSF1B [CD120b], TMEM23 [SGMS1], and FOXP3), ANK3 and ZNF532 were not significantly changed, and one gene was slightly upregulated (LRRC32 [GARP]) in new-onset T1D compared with control Tregs (Fig. 4B). Despite the difference in FOXP3 mRNA, FOXP3 protein levels were not different (Fig. 2C), and indeed, due to multiple layers of posttranscriptional regulation, differences in mRNA expression are not necessarily linked to simultaneous differences in protein expression (20). PCA depicting expression of these six genes showed distinct clustering of new-onset T1D versus control Tregs (Fig. 4C). Tregs from established T1D had intermediate gene expression, with only FOXP3 and ANK3 expression significantly different from that of control subjects.

Biomarker discovery and refinement. A: Top 10 most frequent genes in biomarker discovery pipeline using elastic net cross-validation with α = 0.65; 65% is marked by the dashed line. B: Relative expression levels of the six most frequent genes used in the refined six-gene panel (normalized and log2-transformed NanoString data); common gene name in square brackets. C: PCA representing expression of six genes by sample group. D and E: Biomarker performance of the six-gene panel. Elastic net calculations with α = 0 to include all six genes, and AUC, error rate, PPV, NPV, sensitivity, and specificity were assessed (D). Classification accuracy of the six-gene panel, cutoff for control (ped control) vs. new-onset T1D (new T1D) classification at 0.5, and means from 1,000 panels are shown for ped control and new T1D. Established T1D (est T1D) were run through the algorithm once (E). Throughout, ped control n = 24, new T1D n = 29, and est T1D n = 27. Horizontal lines represent means, with SD as error bars; statistical analysis among three groups by one-way ANOVA with Bonferroni multiple comparisons test. *P ≤ 0.05, **P ≤ 0.01, ****P ≤ 0.0001. PC, principal component.

To explore the utility of these six genes as a biomarker algorithm in T1D, the performance test was repeated with α = 0, restricting panels to these six genes. The six-gene panel improved biomarker performance, resulting in an AUC of 0.84 ± 0.12 and sensitivity of 0.86 ± 0.14 (Fig. 4D and Supplementary Fig. 2A, D, and E), correctly classifying 26 of 29 subjects with new-onset T1D and 20 of 24 control subjects across 1,000 panels (Fig. 4E). In contrast, the gene expression from established T1D Tregs was variable but significantly different from control Tregs, possibly reflecting the different stages of autoimmune activity in these subjects.

Discussion

Changes in Treg function are believed to contribute to T1D, but tools to accurately monitor Tregs throughout the disease course are lacking. We describe a novel gene signature of human Tregs and a refined six-gene algorithm as a sensitive and specific biomarker in new-onset T1D. The data show that monitoring Treg gene expression could potentially be used as a clinically applicable biomarker to track the fitness of Tregs, overcoming many of the current limitations of protein-based markers and functional assays.

Previous attempts to define specific markers of human Tregs focused on comparing CD25+ and CD25− cells (11,12,21,22), thereby introducing bias toward identifying T-cell activation/memory markers rather than Tregs. We developed a new 31-gene signature of human Tregs that can, when used in composite, effectively discriminate between Tregs and Tconvs regardless of their activation state. Although some of the Treg-specific genes in this signature overlap with previously identified ones (21,22), some do not, illustrating the importance of integrating data from activated Tconvs. This signature has broad applicability to T1D and many other diseases because it represents the first tool to define Tregs not limited by similar expression patterns in activated Tconvs.

Although Tregs are believed to be dysfunctional in T1D, defining whether this dysfunction relates to changes in number, intrinsic function, and/or the susceptibility of other immune cells to Treg suppression has been difficult. The proportion of CD25hiCD127loFOXP3+ Tregs in blood of subjects with T1D is normal (3); we confirmed this in a pediatric cohort (including similar levels of FOXP3 protein expression) and further showed that ex vivo Tregs from subjects with T1D have a demethylated TSDR. Whether peripheral T1D Tregs are able to suppress effector cell proliferation is unclear (4,5,23). Subtle changes in function may only be revealed under specific conditions (4), in subfractions of Tregs (6,7), and/or in subjects with specific genetic risk factors (24). For example, single nucleotide polymorphisms affecting IL-2R signaling (CD25 rs12722495, rs3118470, rs2104286) or the function of PTPN2 (rs478582, rs45450798, rs1893217) are associated with changes in Treg function (24–29), so future studies to correlate changes in the Treg signature with genetic risk factors are warranted.

Traditional assays to test Treg function are not well suited for routine clinical implementation to monitor Tregs over time because they are labor intensive, variable, and may not accurately represent in vivo function. In the current study, we show that measuring the expression of a few genes in Tregs sorted from blood represents a simple and sensitive method to detect Treg changes in children with new-onset T1D. Tregs from subjects with new-onset T1D had significantly lower expression of FOXP3, TNFRSF1B (CD120b), and TMEM23 (SGMS1), a small decrease in ZNF532 and ANK3, and slightly increased LRRC32 (GARP). These changes in gene expression are consistent with the notion that functional changes in T1D Tregs are subtle and related to FOXP3 expression, signaling, and stability (3). Preliminary experiments revealed that the six-gene algorithm may also be applicable to unfractionated PBMCs. We found an AUC of 0.64 for control subjects (n = 8) versus subjects with new-onset T1D (n = 6). Future studies with a larger sample size will be needed to fully test the possibility of using the Treg gene signature in PBMCs, which would be ideal for clinical application.

Both TNFRSF1B (CD120b) and LRRC32 (GARP) have been previously shown to be associated with Treg function. There is some evidence that CD120b may mediate some of the pro-Treg functions of tumor necrosis factor-α. For example, tumor necrosis factor-α stimulates Treg proliferation through CD120b (30,31), CD120b signaling stabilizes Tregs in inflammatory environments (31), and expression of CD120b on Tregs is necessary for suppression of colitis (32). GARP is a receptor for latent transforming growth factor-β and can transfer the mature transforming growth factor-β to its receptor on the same cell or another target cell (33). GARP is upregulated upon T-cell receptor activation (34) and is part of a positive feedback loop with FOXP3 (35). Thus, three of the six signature genes have been directly linked to Treg function and stability. The other three might also be involved in Treg signaling because they have been associated with processes important for signal transduction and immunological synapse formation (TMEM23 [SGMS1]), cellular arrangement (ANK3 is integral to the cytoskeleton), and transcriptional control (ZNF532) (36).

The Treg gene signature in subjects with established T1D was variable, with ∼40% of subjects possessing a signature similar to control subjects and the other ∼60% being similar to cells from subjects with new-onset T1D. Patients with established T1D with a disease-classified Treg signature may be suited for Treg-targeted therapy, whereas patients with a control-classified Treg signature might have quiescent autoimmunity. More research in prospective longitudinal samples (e.g., from clinical trials studying the effects of new immunomodulatory treatments) is required to test the possibility that defects in self-tolerance may be more readily detected early in disease when fulminant autoimmunity exists. Additionally, testing of prediabetic samples would be of use to determine whether the Treg signature algorithm predicts progression to diabetes. Alterations in the Treg gene signature also may be detectable more broadly in autoimmunity and/or other immune-mediated diseases, with a correlation to the presence of genetic susceptibility to autoimmunity. Whether similar or distinct changes in Treg gene expression occur in other autoimmune diseases and, if so, whether they are linked to genetic risk and/or change with disease progression/treatment are of significant interest.

In conclusion, we developed a tool to discriminate between Tregs and Tconvs regardless of their activation state, which will have broad utility in monitoring and identifying human Tregs. Because changes in expression of a few genes in the Treg signature accurately classify subjects with new-onset T1D, this assay could be used as a biomarker to track Treg changes over the disease course and to identify individuals who would be best suited for Treg-targeted therapies.

Article Information

Acknowledgments. The authors thank Susan Farrell (British Columbia Children’s Hospital) for subject recruitment, Rusung Tan and Huilian Qin (Child & Family Research Institute [CFRI] and The University of British Columbia) for contributions to the diabetes biobank, and Paul C. Orban and C. Bruce Verchere (CFRI and The University of British Columbia) for ongoing support and discussions.

Funding. This study was supported by grants from the Canadian Institutes of Health Research (CIHR) (MOP-115199), the Canucks for Kids Fund, Genome British Columbia, Genome Québec, and STEMCELL Technologies Inc. A.M.P. holds fellowships from the CIHR Strategic Training Initiative in Health Research Program in Transplantation and the 4 What Matters Foundation and a JDRF postdoctoral fellowship. A.Y.W. was funded by a CFRI Research Methodology Training Grant. A.S. is supported by a CIHR Banting and Best Doctoral Award. M.K.L. receives a Scientist Salary Award from CFRI.

Duality of Interest. No potential conflicts of interest relevant to this article were reported.

Author Contributions. A.M.P. and A.Y.W. contributed to the study design, execution of the experiments, and writing of the manuscript. A.S. performed the bioinformatic and statistical analysis. J.G., Y.K., and D.N. contributed to the sample processing and analysis. C.A.P. contributed to the study design. W.N.H. supervised the gene expression analysis and contributed to the study design. S.J.T. supervised the bioinformatic and statistical analysis and contributed to the study design and writing of the manuscript. C.P. oversaw subject recruitment and sample biobanking and contributed to the writing of the manuscript. M.K.L. contributed to the overall design and execution of the research and writing of the manuscript. M.K.L. is the guarantor of this work and, as such, had full access to all the data in the study and takes responsibility for the integrity of the data and the accuracy of the data analysis.