Epigenetic Changes in Islets of Langerhans Preceding the Onset of Diabetes

Alarming incident rates worldwide are projected to increase the current prevalence of type 2 diabetes (T2D) from 422 million to 592 million in 2035 (1). T2D is a progressive, chronic disorder with a long asymptomatic phase, averting detection for many years (2,3). Better disease management might be possible with earlier detection through robust, sensitive, and easily accessible biomarkers of T2D.

T2D is characterized by chronic hyperglycemia, which is caused by an impaired insulin secretion from pancreatic β-cells and an insulin resistance of target tissues. Aging, a sedentary lifestyle, and obesity contribute to insulin resistance. After long-term exposure to elevated lipid and glucose levels, pancreatic islet function decreases (4,5), which leads to insufficient compensation and a loss of β-cells.

The involvement of epigenetic mechanisms in T2D development emerged as a promising research area (6). One epigenetic modification is DNA methylation, which mainly occurs at the 5′ carbon of cytosine-phosphate-guanine (CpG) sites (7,8). DNA methylation marks are established during prenatal and early postnatal development and function throughout life to maintain the diverse gene expression patterns of different cell types (9,10) but can also arise later in somatic cells either by random events or under the influence of the environment (11,12). Thus, tissue specificity and flexibility of DNA methylation in response to the environment are two major problems to face in order to identify stable epigenetic biomarkers of disease risk. The aim of our translational study was to identify early epigenetic marks related to T2D by uncovering methylome alterations in pancreatic islets of mice that occur before the onset of severe hyperglycemia and assessing prospective T2D risk information conveyed by congruent differential methylation in human blood.

A full description of animals and diets was detailed previously (13,14). At 5 weeks of age, female New Zealand Obese (NZO) mice were placed on a high-fat diet (HFD) (20% kcal protein, 20% kcal carbohydrate, 60% kcal fat; D12492; Research Diets) for 5 weeks. Five weeks after switching the diet, mice were killed during midlight cycle with acute exposure to isoflurane (Fig. 2A). Animal studies were approved by the animal welfare committees of the German Institute of Human Nutrition Potsdam-Rehbruecke and local authorities (Landesamt für Umwelt, Gesundheit und Verbraucherschutz, Brandenburg, Germany).

Islet isolation was performed as described (15). Islets isolated from two to four mice (30–110 islets/mouse) were pooled per sample for RNA sequencing (RNA-seq) and whole-genome bisulfite sequencing (WGBS). The total number of islets used for nucleic acid extraction was 900 and 1,500 for RNA and DNA, respectively. For RNA-seq of diabetes-resistant (DR) mice, four individual islet pools were used that contained islets from 2 mice/pool; for diabetes-prone (DP)mice, five pools were used that contained 2–3 mice/pool. WGBS was performed with five pools per group from 4 animals/pool (Supplementary Table 1). Thus, each sample of islet pools comprised islets from different mice.

Body weight and blood glucose were measured from 7:00 to 9:00 a.m. on a weekly basis by using a Contour blood glucose meter (Bayer). At 5, 7, and 10 weeks of age, body composition and liver fat content were analyzed using nuclear magnetic resonance and computed tomography as described (14).

Plasma adiponectin and leptin levels were measured by Mouse Adiponectin/Acrp30 (DY1119; R&D Systems) and Mouse/Rat Leptin (MOB00; R&D Systems) ELISA kits, respectively. Plasma triglycerides (T2449, F6428, G7793; Sigma), free fatty acids (91096, 91898, 91696; Wako), cholesterol (10017, HUMAN), ALT (12212, HUMAN), AST (12211, HUMAN), and γ-glutamyl transferase (GGT) (12213, HUMAN) levels were measured according to the manufacturer’s protocol.

Total RNA was extracted by using miRNeasy Micro Kit (QIAGEN, Hilden, Germany) according to the manufacturer’s protocol, with additional DNase treatment. RNA samples with RNA integrity number ≥8 (Agilent Bioanalyzer) were selected for RNA-seq. Transcriptome sequencing was carried out by GATC biotech (Konstanz, Germany) on an Illumina HiSeq platform. Adapters were trimmed and reads filtered for quality by using the wrapper Trim Galore! v0.4.2 and Cutadapt 1.9.1 with option phred33. FastQC v0.11.5 was used to check sample quality. Alignment of reads to reference genome was performed with HISAT2 v2.1.0, and fragments per kilobase of exon model per million reads mapped values for transcripts was determined by Cufflinks 2.2.1, both with default options for paired reads. We considered only transcripts with fragments per kilobase of exon model per million reads mapped mean values >1/group. Kyoto Encyclopedia of Genes and Genomes (KEGG) analysis was performed by the using DAVID 7 tool (16), with cutoff enrichment score set to >1.7 and enriched P < 0.05. Network analysis was obtained with Ingenuity Pathway Analysis (IPA) (QIAGEN) (Supplementary Table 1).

Genomic DNA from NZO islets was isolated using Invisorb Genomic DNA Kit II. One microgram of genomic DNA from each pool was bisulfite converted (Zymo Research Corporation, Irvine, CA), and library preparation and sequencing steps were carried out by GATC. WGBS data in fastq format were generated using an Illumina HiSeq platform for further analyses. Raw data have been quality controlled and processed using Trim Galore! v0.4.2, FastQC v0.11.5, Bismark v0.17.08 (17), and MethPipe v3.4.2 (18) (Supplementary Table 1).

A reference genome file was generated by combining a GRCm38.68 B6 reference and a GRCm38p4 single nucleotide polymorphism file in order to exchange all B6 with NZO high-quality single nucleotide polymorphisms. Methylation counting was carried out with MethPipe v3.4.2 default options.

Nonsymmetric CpG sites have been withdrawn, P values have been calculated using log-likelihood ratio test. For the final analysis, CpG sites that fit to the following criteria were used: 1) at least four of five samples with read counts in both groups, 2) average of read counts per group >20, and 3) SD per group >0 in both groups.

The study sample is a nested case-control study derived from the prospective European Prospective Investigation into Cancer and Nutrition (EPIC)-Potsdam cohort study (n = 27,548) designed to estimate the association of baseline measurements to incident T2D (19). Baseline recruitment for EPIC-Potsdam was conducted in Potsdam, Germany, and surrounding municipalities between 1994 and 1998. Included participants’ ages ranged from 35 to 64 years. The study was approved by the ethics committee of the Medical Society of the State of Brandenburg, Germany, and participants provided informed consent. Potential incident cases were systematically and continually identified (self-report, biennial questionnaire, death certificates, tumor registries, clinical records linkage) and verified (by last treating physician and study physician). The final sample for analysis comprised 270 case-control pairs matched on age, sex, fasting time before blood draw, time of day, and season at blood sampling. DNA methylation was measured using the Illumina EPIC 850K array (20). Raw DNA methylation data were processed and normalized using the R package meffil (21). Associations of baseline DNA methylation levels within each CpG site with incident T2D were evaluated in conditional logistic regression models accounting for the prospective nested case-control design using z-standardized and winsorized (95%) β-values adjusted for age, waist circumference, smoking, alcohol intake, leisure time, and physical activity as well as for estimated cell composition (22) and unaccounted batch effects (R package smartSVA [23]). T2D risk information on the basis of DNA methylation of a specific gene was summarized through least absolute shrinkage and selection operator (LASSO) regression computed on all available CpG sites annotated to a respective gene (R package clogitL1). β-Values were adjusted (i.e., use of residuals after regressing each CpG site on the adjustment variables) for the aforementioned variables to ensure that prediction of the outcome was independent of those variables. Area under the receiver operating characteristic curve (ROC-AUC) served as the ranking of the conveyed risk information per gene (Supplementary Table 1). For each CpG site, we further tested whether the association with T2D differed by sex by including interaction terms in the pooled models.

Human pancreatic islets were provided by the Nordic Network for Islet Transplantation, Uppsala University, Sweden. Genomic DNA was extracted from 15 participants with T2D and 34 control participants without T2D. Diabetes status was diagnosed before death of participants; control participants with glycated hemoglobin (HbA_1c)<6% (42 mmol/mol) were selected (24). DNA methylation analysis was reconsidered for the current study (Supplementary Table 1).

Human RNA-seq data were extracted from Gene Expression Omnibus repository GSE50244, and participants without data on HbA_1c (55 without T2D and 26 with T2D) were excluded from the statistical analysis. The National Human Genome Research Institute-European Bioinformatics Institute Genome-Wide Association Study (GWAS) Catalog (25) was used to screen for differentially methylated ortholog genes. The islet expression quantitative trait loci data set from the two published studies (26,27) was downloaded and compared with the human orthologs identified in the current study.

For animal experiments, Welch t test was performed for the comparison of two groups. All statistical tests were conducted using R 3.5.0 software (23 April 2018). Significance levels were set for P < 0.05, 0.01, and 0.001. After Benjamini-Hochberg correction, we used a threshold of P < 0.1.

All mouse data sets are available in the Gene Expression Omnibus repository (GSE143875) and will be publicly available upon acceptance. All human results are available in the Supplementary Tables. EPIC data sets are not publicly available because of data protection regulations. In accordance with German federal and state data protection regulations, epidemiological data analyses of EPIC-Potsdam may be initiated upon an informal inquiry addressed to the secretariat of the Human Study Center ([email protected]). Each request will have to pass a formal process of application and review by the respective principal investigator and a scientific board.

Figure 1 summarizes the strategy of a semiexplorative approach used in the current study to identify early changes in DNA methylation marks related to T2D. We first analyzed the patterns of gene expression by RNA-seq and DNA methylation by WGBS in islets of genetically identical prediabetic mice that are supposed to differ in their later development of hyperglycemia (see below) (Fig. 2). To later translate the findings to human, we focused on genes with conserved CpG sites and further examined them in two different studies: 1) in blood cells of baseline-healthy donors to be related to incident T2D and 2) in islets of donors with and without diabetes. DNA methylation of human samples was measured with different methods: 850K array in blood cells and 450K array in human islets.

The NZO mouse is a model of polygenic obesity and shows a T2D-like phenotype (28). The diabetes susceptibility of NZO female mice depends on the diet and heterogeneity. At the age of 5 weeks when fed an HFD (60% fat), ∼64% of the mice displayed a rapid increase in blood glucose starting at the age of 8 weeks, whereas the other mice were protected from severe hyperglycemia and showed only a moderate increase in blood glucose concentration to maximally 11.5 mmol/L by the age of 18 weeks (Fig. 2B). DP prediabetic mice exhibited a slightly higher body weight than diabetes-resistant (DR) mice (Fig. 2C), which was mainly caused by elevated fat mass (Supplementary Fig. 1). Plasma lipids were not different between DP and DR mice. Only the concentration of leptin was slightly lower and that of adiponectin to some extent higher in DR mice, resulting in a significantly lower leptin-to-adiponectin ratio (Table 1). The immunohistochemical costaining of pancreatic sections for insulin, glucagon, and somatostatin revealed that the proportion of β-, α-, and δ-cells were similar between 10-week-old DR and DP mice (Supplementary Fig. 2). In a new cohort of mice, liver fat content was quantified by computed tomography at weeks 5, 7, 8, and 10, showing significant differences at week 10 (Fig. 2D). The ability to predict later development of severe hyperglycemia was evaluated by using data of early blood glucose and liver density alone or in combination. As shown in the ROC-AUC curves of Fig. 2E–G, prediction of the T2D-like phenotype was most accurate by combining blood glucose and liver fat content.

To study the impact of altered DNA methylation on the islet transcriptome in DP compared with DR mice, islets were isolated at the age of 10 weeks for genome-wide transcriptome and methylome analyses. Comparative analysis of samples from DR and DP mice identified 3,546 differentially expressed genes, and among these, 3,120 mRNAs displayed a fold change >1.5 (unadjusted P < 0.05) (Fig. 3A). To evaluate the molecular events involved in the transition from mild to severe hyperglycemia in DP mice, KEGG pathway enrichment analysis was performed. Transcripts that are lower abundant in islets of DP mice are linked to lysosome, fatty acid metabolism, tricarboxylic acid cycle, and others (Fig. 3B); those that were higher expressed in DP islets are involved in processes of ribosome, oxidative phosphorylation, proteasome function, and DNA replication (enrichment score >2; P < 10⁻³) (Fig. 3C). IPA resulted in five significantly enriched networks mainly related to cell death (Supplementary Table 2) and carbohydrate metabolism. The latter links 14 differentially expressed genes to the transcription factor pancreatic and duodenal homeobox 1 (PDX1), which is crucial for islet function (29,30) (Fig. 3D).

Apart from the alteration in metabolic pathways, RNA-seq analysis revealed changes in the gene expression of 39 chromatin modifiers, 75 transcription factors, 20 RNA-binding proteins, and 11 enzymes involved in the transfer and the maintenance of DNA methylation (unadjusted P < 0.05) (Supplementary Table 3). The differential expression of transcription factors and epigenetic modifiers supports our assumption that epigenetic mechanisms participate in differences in diabetes susceptibility in NZO female mice.

To clarify to what extent early alterations of gene expression in NZO female mice resemble those detected in human islets of donors with and without diabetes, we compared our transcriptome results with published RNA-seq data (27). Interestingly, 1,374 genes displayed differential expression in islets of both DP prediabetic animals and human donors with diabetes (Fig. 3E). This number of overlapping genes was more frequent than statistically expected by chance (χ² test P = 0.02), and several gene products are involved in oxidative phosphorylation, ribosome, and cell cycle pathways (Supplementary Table 4). This strong overlap between human and mouse islet expression data further suggests that the NZO mouse is a suitable model to study molecular alterations that precede onset of T2D (28).

WGBS was assessed in islets of DR and DP mice. After quality control, DNA methylation data were obtained for a 21,862,896 CpG sites. There were no differences in the average of β-methylation between DR (76.67%) and DP (76.11%) islets. The degree of DNA methylation throughout the genome was highly intercorrelated, and global DNA methylation displayed the same pattern in islets of DR and DP mice (Supplementary Fig. 3A and B). A total of 37,628 CpG sites exhibited a different degree of methylation between islets of DP and DR mice; 61% were hypermethylated and 39% hypomethylated in DP islets. The majority of DNA methylation changes were observed in intergenic regions (∼80%), and the remaining changes were distributed within promoters (1.3%), gene bodies (12.6%), first introns (5.4%), and first exons and 3′ untranslated regions (UTRs) (0.7%) (Supplementary Fig. 3C–F). To further investigate genomic regions with consistent and extended differential methylation, we calculated differentially methylated regions (DMRs), genomic regions with at least two differentially methylated CpG sites consistently hypermethylated or hypomethylated within a maximum width of 1,000 base pairs (bp). We detected 223 DMRs in proximity to 211 genes and 1,107 DMRs in intergenic regions (Supplementary Table 5).

To translate our findings on changes in DNA methylation to humans, a fully automated method was developed to identify conserved CpG sites in mouse and human. Using pairwise alignment of mouse and human DNA (University of California, Santa Cruz, Genome Browser mm10/hg19), 4,750 CpG sites revealed complete conservation at the specific position, and 8,711 sites presented a CpG in the surrounding 10 bp (Fig. 4A). Among these conserved 13,461 CpG sites, only 1,519 CpG sites were differentially methylated in islets of DR and DP mice. Most of these CpG sites (900) were intergenic (59%), and 41% were located in the vicinity of genes (619 CpG sites) (Supplementary Table 6).

To relate DNA methylation and gene expression, we compared all upregulated genes with hypomethylated CpG sites in promoter regions (including 5′ UTR, exon 1, intron 1/2) and hypermethylated CpG sites in gene bodies and vice versa for the downregulated genes (31). All CpG sites not located within a promoter (2 kilobases upstream transcription start site) or within a gene were designated as intergenic and excluded from our analysis. As it is unknown which genes are regulated by intergenic DNA methylation sites and as we later compare data with results obtained from the 450K and 850K arrays (see below), we focused on alterations within or in proximity to genes. On the basis of this, 497 differentially expressed genes exhibited changes in DNA methylation (unadjusted P value), of which 176 exhibited at least one conserved CpG site in the human genome (Fig. 4B), of which 39 genes were located in DMRs (Supplementary Tables 5 and 6). The hypergeometric test demonstrated that the observed intersection of the 497 differentially expressed and methylated genes is not due to randomness (P = 6.42e-09). The top five KEGG pathways of the 497 genes, ranked by fold enrichment, are extracellular matrix (ECM)-receptor interaction, long-term potentiation, insulin secretion, and amphetamine and cocaine addiction. Consideration of the affected genes within these pathways (x-axis of Fig. 4C) indicates an overlap for those enriched in insulin secretion and amphetamine and cocaine addiction. This overlap is even stronger with other pathways, such as pancreas secretion, cAMP, and calcium signaling (Supplementary Table 7). The scatterplot in Fig. 4D illustrates the relationship between gene expression and DNA methylation changes of the 497 identified genes. Atp2b1, for instance, is higher expressed in islets of DP mice and exhibits an elevated DNA methylation within the gene body. Network analysis by IPA revealed 14 differentially expressed and methylated genes that can also be connected to PDX1 and to other key islet transcription factors (e.g., NKX2-2 and MAFB) (Fig. 4E).

In a translational approach, we assessed whether DNA methylation of genes identified in mice was related to incident T2D in the EPIC-Potsdam cohort (Fig. 1). Of the 176 genes identified to be differentially expressed in islets of DP and DR mice and that contain at least one conserved CpG site (Fig. 5B), 120 genes with 8,276 annotated CpG sites are covered in the Illumina EPIC 850K array and could thus be analyzed in human samples that were collected a median of 3.8 years before the diabetes diagnosis. In this approach, we contemplated that not only the exact conserved CpG site but also every change in DNA methylation in the ortholog gene could be an epigenetic signature associated with T2D (20) (Table 2). Overall, there was evidence for an association with incident T2D (unadjusted P < 0.05) for 605 CpG sites located in 105 genes (Fig. 5A). We did not observe systematic interaction with sex and, therefore, present results from pooled analyses (data not shown). The most statistically significant associations for single CpG sites were observed for cg25381383 (annotated to MEIS2, odds ratio per Z score 0.33 [95% CI 0.20, 0.55], P = 1.9 × 10⁻⁵), cg11995041 (WWOX, 2.36 [1.5, 3.71], 2.2 × 10⁻⁴), cg19746591 (TAOK3, 3.18 [1.7, 5.94], 1.3 × 10⁻⁴), cg09587151 (ATF7, 0.52 [0.36, 0.74], 1.4 × 10⁻⁴), and cg17429772 (PTPRN2, 2.04 [1.38, 3.02], 2.1 × 10⁻⁴). In addition, the majority of CpG sites with unadjusted P < 0.001 were located in the gene body (Supplementary Fig. 4 and Supplementary Table 8).

To rank the genes by conveyed T2D risk information of the respectively annotated CpG sites, we derived the best model for each gene on the basis of LASSO regression and calculated the discrimination of the respective model (ROC-AUC). The methylation values were adjusted for the above-mentioned adjustment variables by linear regression (subsequent use of the derived residuals). AKAP13 (ROC-AUC 0.73 [95% CI 0.69, 0.77]), TENM2 (0.70 [0.66, 0.74]), CTDSPL (0.68 [0.64, 0.73]), PTPRN2 (0.68 [0.64, 0.72]), and PTPRS (0.67 [0.62, 0.71]) emerged as the top five genes in this study sample (Fig. 5B and Supplementary Table 9). Of note, the majority of genes associated with diabetes incidence have not been described in GWASs for T2D. All GWASs with 1,503 genes associated with diabetes-related traits, including 403 association signals that were identified in an expanded GWAS discovery performed by Mahajan et al. (32), were considered for this comparison, but only 14 genes overlapped with the 105 differentially methylated genes. The comparison with gene expression quantification loci from islets revealed an overlap of 9 genes with the 2,339 genes detected in the Fadista et al. (27) study and of 7 genes that overlapped with 2,853 genes described in the van de Bunt et al. (26) study (Supplementary Table 10). Thus, the current analysis put a specific focus on novel genes, unknown in previous genetic studies, as a putative contributor in the onset of T2D.

To evaluate which of the 120 genes detected in blood cells were also affected by DNA methylation changes in pancreatic islets of donors with diabetes versus control donors without diabetes, we reevaluated 450K array methylation data previously published by Dayeh et al. (24). In total, 99 of 120 genes exhibited 655 differentially methylated CpG sites in islets of donors with diabetes compared with control donors (unadjusted P < 0.05) (Supplementary Table 11). After applying multiple testing corrections, 51 CpG sites located in 23 genes (of which 4 are shown in Fig. 5C) were significantly different between donors with and without diabetes (q < 0.05) and showed the same effect in mouse islets. Those encompass AKAP13, CTDSPL, and PTPRN2, which emerged within the list of the top five predictive genes in blood cells. Hence, we demonstrate that the majority of genes highly associated with future diabetes in blood cells exhibit altered DNA methylation levels in islets after diabetes diagnosis. Among the 99 differentially methylated genes in human islets, 34 showed altered mRNA levels in donors with diabetes according to RNA-seq data published by Fadista et al. (27) (marked in red in Fig. 5B).

Finally, we used more stringent criteria by applying a correction for multiple testing for the initial exploratory transcriptome analysis that was achieved in islets of DP and DR mice and calculated all subsequent analysis again. Benjamini-Hochberg corrections (P < 0.1) resulted in 789 differentially expressed genes (Supplementary Table 11). To observe which of these genes might be affected by DNA methylation, we analyzed those CpG sites within or in proximity that fulfill the following criteria: high sequencing coverage, methylation averages in the range of 10–90%, and methylation differences >10%. By this, we found 479 differentially methylated CpG sites in 274 differentially expressed islet genes. Of these, 253 are covered by the EPIC-850K array, and the majority (236) exhibited an altered DNA methylation on 1,481 CpG sites associated with incident T2D (Supplementary Table 13). Interestingly, 201 of 236 genes exhibited 785 differentially methylated CpG sites in islets of donors with diabetes compared with control donors (unadjusted P < 0.05) (Fig. 6A). After multiple testing, the number of altered CpG sites was reduced to 31 in 23 genes (q < 0.05) (Supplementary Table 14).

Finally, the comparison of gene lists generated by the two approaches used in the current study (Fig. 6A) resulted in an overlap of 19 genes detected in blood cells (Fig. 6A and Supplementary Table 14) of which 17 have an ROC-AUC ≥0.58. Sixteen genes appear to be relevant in islets (e.g., CUX2, DACH1, MEIS2; hypergeometric test P = 1 × 10⁻²⁰) (Supplementary Table 14).

To our knowledge, this is the first translational study to examine changes in DNA methylation in pancreatic islets before the onset of severe hyperglycemia in mice, showing a significant overlap of DNA methylation marks visible in blood cells of people with a diabetes risk. By integrating methylome and transcriptome analysis, we identified ∼500 differentially expressed and methylated genes of which several were linked to insulin secretion and ECM-receptor pathways. From 176 genes with conserved CpG sites in humans, 120 were detected by the EPIC-850K array, and of those, 105 exhibited an association for at least one differentially methylated CpG site in blood cells and incident T2D diagnosis. Furthermore, 99 genes with strong association with T2D incidence in blood cells exhibited altered DNA methylation profiles in islets of donors with diabetes. Accordingly, >80% of the selected candidates for human translation showed significant results in two independent data sets (blood cells and islets).

In pancreatic islets, early molecular alterations occur and mediate islet dysfunction before the onset of diabetes (3). The clarification of the pathomechanisms remains challenging because of the difficulty to investigate them in humans. Therefore, an appropriate mouse model that resembles the metabolic syndrome and T2D with impaired β-cell function, such as the NZO mouse, is of particular interest (28). It was used to screen for methylome alterations arising at a very early stage of the disease.

At 10 weeks of age, in a prediabetic state when blood glucose and liver fat were slightly elevated, the islet gene expression profile reflected the metabolic observations and the expected β-cell failure at later time points. Although DP animals were not diabetic (blood glucose <16.6 mmol/L) at the stage of examination, they exhibited higher blood glucose levels than DR mice, which probably affected the methylome and thereby accelerated diabetes development.

The methylome analysis of the current study provides the first catalog of epigenetic alterations in pancreatic islets before the development of severe hyperglycemia. It identified ∼500 differentially methylated and expressed genes mainly involved in ECM-receptor interaction and insulin secretion. Several candidates are well known in T2D, such as Insr (insulin receptor), Gcg (glucagon), and others described to be affected by DNA methylation changes in islets of patients with diabetes (e.g., Park2, Adcy5). By silencing or overexpressing Park2 and Socs2 in INS-1 cells, their role in insulin secretion was shown (33,34).

To translate results to humans, a list of candidate genes that contain conserved CpG sites in the human sequence was generated and analyzed in blood cells of apparently healthy individuals. By this, we uncovered evidence for an association of DNA methylation marks in the selected genes and incident T2D. The most significant association for single CpG sites was observed for the Meis homeobox 2 (MEIS2) gene. Rat experiments demonstrated that Meis2 protein forms a multimeric complex with Pbx1 and modulates the transcriptional activity of Pdx1 (35).

In an attempt to summarize the risk information per gene through LASSO regression, we identified that nearly all genes showed at least some degree of discriminatory capacity. The top predictive genes even reached relatively high ROC-AUC values of ∼0.7, which is noteworthy given that information of other risk factors was regressed out before analysis. However, overfitting still plays a role in this setting, and the interpretation of the absolute predictive capacity is optimistic and requires replication in other human cohorts. The highest predictive gene was AKAP13, a scaffold protein that plays an essential role in assembling signaling complexes downstream of several G-protein–coupled receptors. Such a downstream protein is the small GTPase RhoA, which is involved in cytoskeleton organization, cell migration, and cell cycle (36). Furthermore, PTPRN2 (receptor-type tyrosine-protein phosphatase N2), one of the top five predictive genes in our study, is required for the accumulation of normal levels of insulin-containing vesicles, preventing their degradation, and plays a role in glucose-stimulated insulin secretion (37). None of our 105 candidates were genome-wide statistically significant in the two largest prospective epigenetic-wide association studies (38,39). In a nested case-control study, Chambers et al. (39) identified an association between differential methylation at five genetic loci (ABCG1, PHOSPHO1, SOCS3, SREBF1, and TXNIP) with T2D incidence among Indian Asians and Europeans. Recently, Cardona et al. (38) confirmed the association of ABCG1, SREBF1, and TXNIP with T2D incidence in a British population (EPIC-Norfolk study) and identified 15 additional CpG sites. Similar to all genome-wide approaches, a conservative threshold for statistical significance to account for multiple testing was used in these studies. Although such corrections are mandatory in fully explorative analyses, associations with lower precision or effect size might be overlooked. To our knowledge, both mentioned publications did not publish full results of all analyzed CpG sites, thereby not allowing a comparison with our results.

One limitation of our broad screening and translational approach is that the results somehow differ, depending on the type of analysis. We mainly focused on a semiexploratory approach, 1) starting from a selected number of differentially expressed and methylated genes determined in islets of a suitable mouse model, 2) focusing on CpG sites conserved between mouse and human, and 3) considering the direction of expression changes. By this, we detected 105 genes with an altered DNA methylation associated with incident T2D, and most of them (94%) showed an altered methylation in islets of patients with T2D. Correcting the transcriptome data for multiple comparisons resulted in a smaller list of differentially expressed and methylated genes in mouse islets (274) of which the majority (236 genes) revealed an altered DNA methylation in blood cells to be associated with incident T2D diagnosis, of which 85% carried an altered DNA methylation in islets of donors with diabetes (Fig. 6). The number of the genes that overlapped in both analyses is relatively low (19 genes in blood cells and 16 in islets) but with a high degree of significance (hypergeometric test P = 1 × 10⁻²⁰). Thus, these overlapping genes appear to be the most relevant for the prediction of the disease.

A second limitation of our study was discarding all CpG sites located in intergenic regions because earlier epigenetic studies have shown that T2D-associated DMRs are located in these regions (33) and that environmental and nutritional effects induce most DNA methylation changes in intergenic regions (40–42). Recent discoveries on the three-dimensional structure of the epigenome (43,44) have provided a better understanding on genomic folding and looping. For example, enhancers and superenhancers can be located ∼50–100 kilobases up- or downstream of the transcription starting site. Therefore, it is difficult to link intergenic methylation alterations to a specific differentially expressed gene.

We and others identified some changes in DNA methylation in tissues to be mirrored in blood cells (20,45,46), whereas other alterations can appear exclusively in one tissue (47). Although, our results clearly showed that several DNA methylation marks are visible in islets and blood cells, it is difficult to speculate which molecular mechanisms are responsible for these observations. One explanation is that these CpG sites (Fig. 5A and B) are established during early development (e.g., metastable epiallele) (48) and remain stably affected, and the second explanation is that the slightly elevated ectopic fat content in liver and elevated blood glucose levels modulate DNA methylation of specific loci in different tissues in a similar manner.

In summary, we identified novel epigenetic alterations in murine pancreatic islets that occur in the prediabetic state with a mild hyperglycemia before β-cell failure. Alterations in 105 genes were mirrored in blood cells of participants with prediabetes, and 655 CpG sites within 99 genes were also affected in pancreatic islets of participants with diabetes. Thereby, our study provides novel DNA methylation marks as putative biomarkers for T2D.

Acknowledgments. The authors thank Christine Gumz, Andrea Teichmann, and Anett Helms (German Institute for Human Nutrition Potsdam-Rehbruecke) for excellent technical assistance. Regarding the EPIC-Potsdam Study, the authors thank the Human Study Centre of the German Institute of Human Nutrition Potsdam-Rehbruecke, namely the trustee and data hub for data processing and the biobank for the processing of biological samples, and the head of the Human Study Centre, Manuela Bergmann, for contribution to the study design and leading the underlying processes of data generation. Furthermore, the authors thank all EPIC-Potsdam participants for invaluable contributions to the study.

Funding. This work was supported by the Federal Ministry of Science, Germany (BMBF:DZD grant 82DZD00302) and the European Union (grant SOC 95201408 05 F02) for the recruitment phase of the EPIC-Potsdam Study, by the German Cancer Aid (grant 70-2488-Ha I) and the European Community (grant SOC 98200769 05 F02) for the follow-up of the EPIC-Potsdam Study, and by a grant from the German Federal Ministry of Education and Research (BMBF) to the German Center for Diabetes Research (DZD) (82DZD00302) and the State of Brandenburg. For human islets study, the work was supported by Novo Nordisk Foundation, Swedish Research Council, Region Skåne, Medical Training and Research Agreement (ALF); ERC-Co Grant (PAINTBOX, No. 725840); European Foundation for the Study of Diabetes; EXODIAB; Swedish Foundation for Strategic Research (IRC15-0067) and Swedish Diabetes Foundation grants (DIA2018-328).

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

Author Contributions. M.O., S.S., and M.S. performed data acquisition for the animal study. M.O., S.S., F.E., M.J., C.W., T.R., L.Z., and P.G. performed the data analysis. F.E. and M.S. designed and provided data on DNA methylation in blood cells from the EPIC-Potsdam study. T.R. and C.L. provided and analyzed data on human DNA methylation in pancreatic islet samples and critically edited and revised the manuscript. M.O. and A.S. contributed to the study conception and design and wrote the manuscript. All authors read and approved the final manuscript. A.S. 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.