Novel Subgroups of Type 2 Diabetes Display Different Epigenetic Patterns That Associate With Future Diabetic Complications.

OBJECTIVE: Type 2 diabetes (T2D) was recently reclassified into severe insulin-deficient diabetes (SIDD), severe insulin-resistant diabetes (SIRD), mild obesity-related diabetes (MOD), and mild age-related diabetes (MARD), which have different risk of complications. We explored whether DNA methylation differs between these subgroups and whether subgroup-unique methylation risk scores (MRSs) predict diabetic complications. RESEARCH DESIGN AND METHODS: Genome-wide DNA methylation was analyzed in blood from subjects with newly diagnosed T2D in discovery and replication cohorts. Subgroup-unique MRSs were built, including top subgroup-unique DNA methylation sites. Regression models examined whether MRSs associated with subgroups and future complications.
RESULTS: We found epigenetic differences between the T2D subgroups. Subgroup-unique MRSs were significantly different in those patients allocated to each respective subgroup compared with the combined group of all other subgroups. These associations were validated in an independent replication cohort, showing that subgroup-unique MRSs associate with individual subgroups (odds ratios 1.6-6.1 per 1-SD increase, P < 0.01). Subgroup-unique MRSs were also associated with future complications. Higher MOD-MRS was associated with lower risk of cardiovascular (hazard ratio [HR] 0.65, P = 0.001) and renal (HR 0.50, P < 0.001) disease, whereas higher SIRD-MRS and MARD-MRS were associated with an increased risk of these complications (HR 1.4-1.9 per 1-SD increase, P < 0.01). Of 95 methylation sites included in subgroup-unique MRSs, 39 were annotated to genes previously linked to diabetes-related traits, including TXNIP and ELOVL2. Methylation in the blood of 18 subgroup-unique sites mirrors epigenetic patterns in tissues relevant for T2D, muscle and adipose tissue.
CONCLUSIONS: We identified differential epigenetic patterns between T2D subgroups that associated with future diabetic complications. These data support a reclassification of diabetes and the need for precision medicine in T2D subgroups.

Introduction

Diabetes is responsible for >1.5 million deaths per year (World Health Organization 2021) (1). A better prediction, prevention, and targeted treatment of diabetes and its complications may decrease mortality rates and reduce the burden of this disease. Diabetes is mainly classified into type 1 diabetes and type 2 diabetes (T2D). With this traditional classification, T2D includes >85% of all patients with diabetes. T2D is, however, a complex and heterogeneous disease, influenced by genetic, epigenetic, and environmental factors and characterized by several pathological conditions, including insulin resistance, β-cell dysfunction, and elevated hepatic glucose production (2). Classifying T2D as one group has been shown to be insufficient to adequately treat diabetes and predict related complications (3). New reclassifications of T2D were therefore performed in the All New Diabetics in Scania (ANDIS) cohort and in several other cohorts (4–9). In the original report, four different subgroups of T2D were identified based on unsupervised data-driven clustering analysis of six phenotypes: age at onset of diabetes, BMI, HbA1c at diagnosis, HOMA2-B (measure of β-cell function), HOMA2-IR (measure of insulin resistance), and GAD autoantibodies. These novel subgroups were labeled as severe insulin-deficient diabetes (SIDD), severe insulin-resistant diabetes (SIRD), mild obesity-related diabetes (MOD), and mild age-related diabetes (MARD) (4). The subgroups have different patient characteristics and risk of diabetic complications (4). Differences in their genetic, metabolomic, and proteomic signatures further support that diverse etiologies exist between the subgroups (10,11). This reclassification may hence provide a better basis for understanding differences in patients with T2D, representing an important step toward precision medicine in diabetes.

Our group and others have found epigenetic differences in tissues from patients with T2D versus control subjects (11–15), demonstrating that epigenetic mechanisms contribute to the pathogenesis of T2D. Moreover, there has been an increasing interest in identifying blood-based epigenetic biomarkers for risk assessment in patients with diabetes. For example, DNA methylation in blood was associated with future T2D, insulin secretion, and response to therapy (12,17–19). However, it remains unknown whether the epigenetic patterns differ between the novel subgroups of T2D and whether these epigenetic differences may predict complications in patients with newly diagnosed diabetes.

Therefore, we analyzed the methylome in the blood of patients with newly diagnosed T2D from the prospective ANDIS cohort, and our first goal was to investigate whether DNA methylation differs between the four T2D subgroups identified by Ahlqvist et al. (4) and to find “subgroup-unique” methylation sites (i.e., sites that show different methylation levels in one subgroup compared with each of the other subgroups) (Supplementary Fig. 1). We then tested whether combined subgroup-unique methylation risk scores (MRSs) generated from identified top subgroup-unique methylation sites in the ANDIS discovery cohort 1) associated with T2D subgroups in replication cohorts and 2) associated with future diabetic complications such as cardiovascular disease (CVD), chronic kidney disease (CKD), and retinopathy (Supplementary Fig. 1).

Research Design and Methods

Study Populations

ANDIS is an ongoing prospective study of patients with newly diagnosed diabetes that aims to document all new incidences of diabetes within the Scania region in Southern Sweden (https://andis.ludc.med.lu.se) (4,10,19). Blood samples for DNA extraction are taken at registration (i.e., within 1 year from diagnosis of diabetes). The ANDIS protocol was approved by the Regional Ethical Review Board in Lund, Sweden (584/2006, 2011/354, 2014/198).

All New Diabetics in Uppsala County (ANDiU) is a similar study to ANDIS but includes patients with newly diagnosed diabetes living in the Uppsala region (https://www.andiu.se/) (19). The ANDiU study protocol was approved by the Regional Ethics Review Committee in Uppsala, Sweden (2011/155).

Discovery and Replication Cohorts

To study the association between DNA methylation in blood and the recently defined novel subgroups of T2D (4), we included 280 patients with T2D from ANDIS who were previously assigned to the four novel T2D subgroups based on unsupervised clustering (4) and who had available DNA methylation data from blood at diagnosis in the ANDIS discovery cohort. We selected these 280 patients while blind to their subgroup information based on a power calculation from a previous study showing 80% power to detect differences in DNA methylation of 4,000 sites with false discovery rate (FDR) of <5% (13).

We are only aware of a few cohorts, other than ANDIS (e.g., ANDiU), with available blood samples at T2D diagnosis and available phenotypes required for clustering of subgroups in newly diagnosed subjects with T2D. We therefore selected one replication cohort from ANDIS, the ANDIS replication cohort, including 76 additional patients with T2D, and one replication cohort from ANDiU, the ANDiU replication cohort, including 197 patients with T2D.

Clinical characteristics of these cohorts are shown in Supplementary Table 1 and Fig. 1. The inclusion criteria and the flowchart of the selection of patients in these cohorts are displayed in Supplementary Fig. 2.

Phenotype Measurements

Age and HbA1c were considered at diagnosis, and BMI, HOMA2-B, and HOMA2-IR were measured at registration in ANDIS and ANDiU. Standard protocols were applied for measuring weight and height to calculate BMI (kg/m2). HbA1c was measured using the Variant II Turbo HbA1c Kit 2.0 (Bio-Rad Laboratories, Copenhagen, Denmark) (4). C-peptide concentrations for HOMA2 were measured using an electrochemiluminescence immunoassay on cobas e 411 (Roche Diagnostics, Mannheim, Germany) or a radioimmunoassay (Human C-peptide RIA; Lincom, St Charles, MO; or Peninsula Laboratories, Belmont, CA) and used with the HOMA calculator (20).

Diabetic Complications

CVD was defined as having had either coronary events (defined by ICD-10 codes I20–I21, I24, I251 and I253–I259) or stroke (defined by ICD-10 codes I60, I61, I63, and I64). CKD was defined as a minimum of two measurements of estimated glomerular filtration rate (eGFR) <60 mL/min/1.73 m2 for >90 days or a single measurement of eGFR <15 mL/min/1.73 m2 (kidney failure) (4). eGFR was calculated with the MDRD-4 study equation (21). Diagnosis of diabetic retinopathy was based on ICD-10 codes E113 and H36.0. Patients with complications before DNA methylation samples were excluded for the respective analyses. Analyses related to complications were done in combined ANDIS and ANDiU cohorts to improve statistical power due to the modest number of individuals with complications in each individual cohort. Phenotypes of individuals with each complication are presented in Supplementary Tables 2–4.

Genome-Wide DNA Methylation Analysis

Participants’ whole blood samples were taken at ANDIS and ANDiU registration, and DNA was then extracted using the Gentra Puregene Blood kit (Qiagen, Hilden, Germany). DNA methylation analysis of the ANDIS discovery and replication cohorts was performed at two different times and places, the Swegene Centre for Integrative Biology at Lund University (SCIBLU) genomics center and at Lund University Diabetes Centre, respectively. ANDiU samples were analyzed at Lund University Diabetes Centre. Bisulfite was used to treat 500–1000 ng of genomic DNA with the EZ DNA Methylation Kit (Zymo Research, Irvine, CA). DNA methylation was analyzed for all participants using Illumina MethylationEPIC BeadChip microarrays (Illumina, San Diego, CA) according to the Infinium HD assay methylation protocol. Detailed information about quality control and bioinformatic analyses of the genome-wide DNA methylation data are available in Supplementary Fig. 3 and elsewhere (19).

DNA Methylation in Other Tissues

To test whether DNA methylation of sites included in the subgroup-unique MRSs in blood mirror DNA methylation levels in other tissues, we used Illumina 450K array DNA methylation data from blood, skeletal muscle, and adipose tissue taken from the same subject in the Monozygotic Twin cohort (14). Here, methylation data were extracted if methylation sites in MRSs were also covered by the 450K array. Twins with available methylation data were included (Supplementary Table 5). Characteristics of the full twin cohort used for these analyses and additional information has been previously published (14,18).

Statistical Analysis

Statistical analyses were performed using R software. Clinical patient data are presented as means (SD) or percentages. Differences between the four subgroups regarding continuous clinical variables were analyzed with a Kruskal-Wallis one-way ANOVA and a Dunn post hoc analysis corrected for multiple testing using Benjamini-Hochberg procedure. Differences in categorical variables were assessed using a Pearson χ2 test.

To find differences in DNA methylation between the four T2D subgroups, we first performed an ANCOVA adjusting for sex in the ANDIS discovery cohort (Supplementary Fig. 1). Here, a BenjaminiHochberg was applied to correct for multiple testing, and methylation sites with FDR <5% (q < 0.05) were included in further analyses. The X-chromosome was then removed to mitigate the effect of sex on DNA methylation data. Pairwise comparisons were then used to identify subgroup-unique methylation sites, which were defined as sites with differences in methylation levels in one subgroup compared with the methylation levels in all of the other subgroups based on q < 0.05. Here, we did six Bonferroni corrected pairwise comparisons using linear regression models adjusted for sex for each of the sites with q < 0.05 in the ANCOVA.

To integrate epigenetic information across the identified subgroup-unique methylation sites, we calculated a weighted MRS for each subgroup, so-called subgroup-unique MRS. To calculate these weighted MRSs, the sum of the standardized methylation values for each of the included subgroup-unique sites was multiplied by the β-coefficient for the respective site (19,22) (Eq. 1).

These β-coefficients were obtained from sex-adjusted linear regression models for each of the included subgroup-unique sites in the ANDIS discovery cohort. Here, the methylation value for each site was the dependent variable, while the respective subgroup versus the combined group of all other subgroups was the binary independent variable. To select the best combination of subgroup-unique methylation sites to be included in the MRSs, the subgroup-unique sites were 1) rank-ordered based on their significance using q-values in the ANDIS discovery cohort; and 2) subgroup-unique sites included in each subgroup-unique MRS were then selected starting with the highest rank (lowest q-value) and going down in rank until the best possible combination of sites were included in the MRS based on its ability to discriminate between subjects with a particular diabetes subgroup and those without, performing separate analyses in the ANDIS discovery cohort and the ANDIS replication cohort. The ability to discriminate subjects was based on the best area under the curve using C-statistics in both ANDIS cohorts (Supplementary Table 6). Subsequently, four different subgroup-unique MRSs (i.e., SIDD-MRS, SIRD-MRS, MOD-MRS, and MARD-MRS) were generated for each person independently of which subgroup they belonged to. MRSs were adjusted for cell composition using a reference-based method (23). MRSs were generated in the same way in the independent ANDiU replication cohort, and their ability to discriminate between subjects with a particular diabetes subgroup and those without was examined using linear and logistic regression in crude models and when adjusting for the clinical variables defining the subgroups.

To evaluate whether the subgroup-unique MRSs associate with future diabetic complications, sex-adjusted weighted Cox regression models were applied. Hazard ratios (HRs) are presented with 95% CIs, and P < 0.05 was considered statistically significant. A statistical power of 85% (α = 0.05) was achieved with a sample size of 500 (probability of events: 0.15) and assuming a HR of 0.5 or 2 and an SD of 0.5. Regression models were not adjusted for age at onset, BMI, HbA1c, HOMA2-B, or HOMA2-IR due to multicollinearity with the subgroup-unique MRSs.

To examine the correlation between methylation in blood and methylation in other tissues of sites included in subgroup-unique MRSs, Pearson correlation tests were performed. Benjamini-Hochberg was used to correct for multiple testing, and q < 0.05 was considered significant.

Results

Different DNA Methylation Patterns in T2D Subgroups

We explored whether DNA methylation in blood is associated with the four novel T2D subgroups (SIDD, SIRD, MOD, and MARD) using the ANDIS discovery cohort, the ANDIS replication cohort, and the ANDiU replication cohort. These cohorts include newly diagnosed subjects with T2D who had previously been assigned to a subgroup using data-driven clustering and who had DNA methylation data available. In line with our previous study (4), significant differences in age at onset of diabetes, BMI, HbA1c, HOMA2-B, and HOMA2-IR were found between the four T2D subgroups in all three cohorts. Subjects with SIDD had higher HbA1c levels, subjects with SIRD had higher HOMA2-IR and HOMA2-B, subjects with MOD had higher BMI and lower age, whereas subjects with MARD were the oldest of all subgroups (Fig. 1 and Supplementary Table 1).

Figure 1

Patient distribution and phenotype characteristics by T2D subgroups in the discovery and replication cohorts. Phenotypes were measured in the ANDIS and ANDiU cohorts. Included patients were previously defined as SIDD, SIRD, MOD, or MARD. Pie charts show the subgroup distribution in the ANDIS discovery cohort (n = 280) (A), the ANDIS replication cohort (n = 76) (B), and the ANDiU replication cohort (n = 197) (C). Box plots show the distribution of age at diagnosis, BMI, HbA1c, HOMA2-B, and HOMA2-IR, and bar charts show the prevalence of male sex for each T2D subgroup in the respective cohort. Statistical differences between the subgroups were evaluated using Kruskal-Wallis for continuous variables and χ2 test for categorical variables. P < 0.05 was considered significant. Post hoc pairwise comparisons for continuous variables were done using the Dunn test, including correction for multiple testing based on Benjamini-Hochberg. Significance is indicated as *q < 0.05, **q < 0.01, and ***q < 0.001. For detailed characteristics see Supplementary Table 1.

Supplementary Fig. 1 presents our study design. First, we assessed whether DNA methylation of any individual sites associated with the T2D subgroups in the ANDIS discovery cohort. Here, 22,034 sites showed differences in methylation between any of the four subgroups based on an ANCOVA adjusted for sex (FDR <5%, q < 0.05) (Supplementary Table 7). We then performed post hoc pairwise comparisons to identify “subgroup-unique” sites among these sites (i.e., sites that showed different methylation levels in one subgroup compared with the level in each of the other subgroups in the ANDIS discovery cohort). We identified 4,465 subgroup-unique methylation sites, including 56 sites unique for SIDD, 74 sites unique for SIRD, 4,135 sites unique for MOD, and 200 sites unique for MARD in the ANDIS discovery cohort (Supplementary Table 8).

Next, we selected top-ranked subgroup-unique methylation sites from Supplementary Table 8 to build subgroup-unique MRSs (SIDD-MRS, SIRD-MRS, MOD-MRS, and MARD-MRS) that best discriminate T2D subgroups in not only the ANDIS discovery cohort but also in the ANDIS replication cohort. Based on this, we included 54 SIDD-unique sites to generate SIDD-MRS, 2 SIRD-unique sites for SIRD-MRS, 31 MOD-unique sites for MOD-MRS, and 8 MARD-unique sites for MARD-MRS (Supplementary Table 9). The subgroup-unique MRSs were significantly different in patients allocated to each respective subgroup compared with the combined group of all other subgroups in both the ANDIS discovery and ANDIS replication cohorts (Fig. 2). The MRSs also differed between the four subgroups, showing a higher SIDD-MRS in individuals with SIDD versus all other subgroups, a higher SIRD-MRS in individuals with SIRD versus all other subgroups, a higher MOD-MRS in individuals with MOD versus all other subgroups, and a higher MARD-MRS in individuals with MARD versus all other subgroups (Supplementary Fig. 4). Moreover, these subgroup-unique MRSs remained significant after adjusting for cell composition (Supplementary Table 10). We then visualized the methylation level of the sites included in the subgroup-unique MRSs in relation to all sites analyzed with the MethylationEPIC BeadChip. Supplementary Fig. 5 shows that a large proportion of SIDD-MRS sites are hypermethylated, whereas overall SIRD-MRS sites, MOD-MRS sites, and MARD-MRS sites have intermediate levels of methylation.

Figure 2

Subgroup-unique MRSs associate with T2D subgroups and future diabetic complications and play a biological function in the pathogenesis of T2D. The respective subgroup-unique MRSs differ statistically significantly between patients with SIDD, SIRD, MOD, and MARD and patients without the respective T2D subgroup in the ANDIS discovery (A), ANDIS replication (B), and in the independent ANDiU replication cohort (C). Patients within each subgroup had statistically significantly (P < 0.05) higher subgroup-unique MRSs compared with the combined group of all other subgroups. Differences in MRSs were compared using the Mann-Whitney U test. D: Subgroup-unique MRSs associate with the T2D subgroups in the independent ANDiU replication cohort (n = 197). ORs are shown per 1-SD increase in MRSs. In the logistic regression model, the dependent variable is the corresponding subgroup for each MRS vs. the combined group of all other subgroups, so for SIDD-MRS it is SIDD vs. non-SIDD individuals, for SIRD-MRS it is SIRD vs. non-SIRD individuals, for MOD-MRS it is MOD vs. non-MOD individuals, and for MARD-MRS it is MARD vs. non-MARD individuals. E: Associations between subgroup-unique MRSs and the risk of developing diabetic complications during 8 years of follow-up (mean ∼4.5 years) in the combined ANDIS discovery, ANDIS replication, and ANDiU replication cohorts. P < 0.05 was considered significant. The results for the sex-adjusted weighted Cox regression are presented as HRs and 95% CIs. For CVD, there are 410 control subjects and 76 case subjects (n = 486); for CKD, there are 444 control subjects and 73 case subjects (n = 517); and for diabetic retinopathy, there are 490 control subjects and 54 case subjects (n = 544). CVD was defined as having had either stroke (ICD-10 codes I60, I61, I63, and I64) or coronary events (ICD-10 codes I20-I21, I24, I251, and I253-I259). CKD was defined as having had an eGFR <60 mL/min/1.73 m2 for a minimum period of 90 days or a single measurement of eGFR <15 mL/min/1.73 m2. Diagnosis of diabetic retinopathy was based on ICD-10 codes E113 and H36.0. MRSs were normalized to show the risk per 1-SD increase. Patients with the respective complication before DNA methylation samples were excluded for the respective analyses. F and G: Relevant genes annotated to the 95 sites included in the subgroup-unique MRSs associated with diabetes, NAFLD, and/or with some subgroup-defining phenotypes and might therefore be important in the pathogenesis of T2D. We performed a systematic literature search using each gene symbol and the following terms: diabetes, insulin secretion/β-cell function, insulin resistance, obesity, age, and NAFLD. Of the 72 genes, 39 (54%) have been associated with diabetes and/or with some characteristics which defined the subgroups or NAFLD (F), and when looking at individual subgroup-unique MRSs, 23 of 44 genes (52%) included in SIDD-MRS, both genes (100%) included in SIRD-MRS, 12 of 21 genes (57%) included in MOD-MRS, and 2 of 5 genes (40%) included in MARD-MRS were associated with any of the terms representing the subgroup traits (G).

We proceeded to validate whether these subgroup-unique MRSs could discriminate between the four subgroups in the independent ANDiU replication cohort. Importantly, the respective MRSs were statistically significantly different between one subgroup and the combined group of all other subgroups (Fig. 2) and also differed between the four subgroups (SIDD-MRS, P = 3.8e−04; SIRD-MRS, P = 5.4e−04; MOD-MRS, P = 9.8e−14; and MARD-MRS, P = 1.5e−09) in a similar pattern observed in the ANDIS cohorts (Supplementary Fig. 4). When a single subgroup corresponding to the respective MRS was taken as the reference category (i.e., SIDD for SIDD-MRS, SIRD for SIRD-MRS, MOD for MOD-MRS, and MARD for MARD-MRS), all other subgroups had a significantly lower MRS in sex-adjusted models, which shows the ability of each MRS to differentiate the reference subgroup from the others in an independent replication cohort (Table 1). The only one that did not reach significance was the comparison between SIDD versus SIRD for the SIRD-MRS, although SIDD shows a lower MRS than SIRD. Moreover, higher values of all MRSs were associated with a higher probability of clustering to a particular subgroup; that is, higher SIDD-MRS was associated with SIDD (odds ratio [OR] 2.08, P = 2e−04), higher SIRD-MRS was associated with SIRD (OR 1.61, P = 0.011), higher MOD-MRS was associated with MOD (OR 6.06, P = 1.85e−09), and higher MARD-MRS was associated with MARD (OR 2.52, P = 6.72e−07) (Fig. 2). These associations remained with similar effect sizes after adjusting for the primary variable defining each respective subgroup, except for SIDD-MRS, where the effect was lost after adjusting for HbA1c (Supplementary Table 11). Overall, these results show that the four novel subgroups of T2D display different epigenetic patterns in discovery and validation cohorts.

Table 1

Associations between the four subgroups and subgroup-unique MRSs in the ANDiU replication cohort

	ANDiU replication cohort adjusted for sex
	β-Coefficient (SE)	P
SIDD-MRS
SIDD	0 (Ref.)
SIRD	−0.12 (0.05)	0.021
MOD	−0.13 (0.05)	0.008
MARD	−0.19 (0.04)	2.0e−05
SIRD-MRS
SIDD	−0.02 (0.01)	0.168
SIRD	0 (Ref.)
MOD	−0.05 (0.01)	3.8e−05
MARD	−0.02 (0.01)	0.037
MOD-MRS
SIDD	−0.79 (0.15)	3.7e−07
SIRD	−0.91 (0.16)	2.4e−08
MOD	0 (Ref.)
MARD	−1.14 (0.13)	6.3e−15
MARD-MRS
SIDD	−0.08 (0.02)	0.002
SIRD	−0.09 (0.02)	2.9e−04
MOD	−0.17 (0.02)	2.0e−11
MARD	0 (Ref.)

Linear regression coefficients for the associations between the four subgroups and the MRSs, taking the corresponding subgroup for each MRS as the reference group.

Subgroup-Unique Methylation Risk Scores Associate With Future Diabetic Complications

A previous study found that patients allocated to certain T2D subgroups show a higher risk of developing diabetic complications, including CKD or diabetic retinopathy (4). We therefore tested whether subgroup-unique MRSs were associated with future complications. To increase statistical power and due to the modest number of patients who develop complications during follow-up, discovery and replication cohorts were combined for complication-related analyses. Their characteristics are presented in Supplementary Tables 2–4. Sex-adjusted weighted Cox regression models showed that all MRSs were associated with development of future CVD, all but SIDD-MRS were associated with future CKD, and none of them were associated with future retinopathy during a mean follow-up of 4.5 years (Fig. 2). Higher SIDD-MRS (HR 0.72, P = 0.032) and MOD-MRS (HR 0.65, P = 0.001) were associated with a lower risk of developing CVD, whereas higher SIRD-MRS (HR 1.47, P = 0.002) and MARD-MRS (HR 1.41, P = 0.007) were associated with a higher risk of future CVD. Regarding CKD, higher MOD-MRS (HR 0.50, P = 3.11e−07) was associated with a lower risk, whereas higher SIRD-MRS (HR 1.55, P = 0.007) and MARD-MRS (HR 1.90, P = 1.72e−06) were associated with a higher risk of developing renal disease. These associations remained significant after further adjustment for blood cell types, except for SIRD-MRS and the risk for CKD (P = 0.09) (Supplementary Table 12). Owing to multicollinearity, these associations could not be adjusted by subgroup-defining phenotypes. However, none of these phenotypes, except age and HOMA2-B, were associated with future vascular and kidney complications (Supplementary Table 13). We also found associations between some T2D subgroups and complications in the current study (Supplementary Table 14). However, none of the subgroup-unique MRS, subgroup-defining phenotypes, or T2D subgroups generated area under the curves >0.75 (Supplementary Table 15).

Biological Function of Subgroup-Unique Sites

To better understand the biological function of the 72 genes annotated to the 95 differentially methylated sites included in any of the subgroup-unique MRSs, we performed a systematic literature search using each gene symbol and the following terms: diabetes, insulin secretion/β-cell function, insulin resistance, obesity, age, and nonalcoholic fatty liver disease (NAFLD). Any study with any of these search terms and any of the 72 gene symbols was considered. This showed that 39 of the 72 genes (54%) have been associated with diabetes and/or with some characteristics defining the subgroups or NAFLD (Supplementary Table 16). Furthermore, among the 72 genes, 26 genes (36%) have been associated with diabetes (e.g., LMNB2, NBPF20, RREB1, IFIH1), 15 (21%) with insulin secretion/β-cell function (e.g., TXNIP, TFEB), 9 (12%) with insulin resistance (e.g., GRK5, SOD3), 20 (28%) with obesity (e.g., GSN, MOGAT1, RREB1, STK3, SLC6A4), 16 (22%) with age (e.g., ELOVL2, SOD3, TRIM59, TFEB), and 11 (15%) with NAFLD, highlighting the relevance of identified methylation sites for stratification of the subgroups (Fig. 2). When looking at individual subgroup-unique MRSs, 23 of 44 genes (52%) included in SIDD-MRS, both genes (100%) included in SIRD-MRS, 12 of 21 genes (57%) included in MOD-MRS, and 2 of 5 genes (40%) included in MARD-MRS were associated with any of the terms representing the subgroup traits (Fig. 2 and Supplementary Table 16). Genes annotated to SIDD-MRS sites include AATK, CPLX1, CTDSPL, GRK5, LMNB2, RREB1, SMARCA4, SOD3, SYT2, and TXNIP, which play a role in insulin secretion/β-cell function; genes annotated to SIRD-MRS sites include RAB27B and RBL2, previously associated with diabetes; genes annotated to MOD-MRS sites include ELOVL2, PDGFC, SCN9A, SLC6A4, and TFEB, previously associated with obesity; and genes annotated to MARD-MRS sites include CRMP1 and RNF170, previously associated with age (Supplementary Table 16).

Cross-Tissue Methylation in Different Tissues

Finally, we examined whether blood-based DNA methylation of the 95 sites included in any of the subgroup-unique MRSs mirror methylation in other central tissues for T2D. Methylation of these sites was compared between blood, skeletal muscle, and adipose tissue. We used 450K methylation array data from blood, muscle, and adipose tissue of the Monozygotic Twin Cohort because we had access to data from the same individuals for these cell types (Supplementary Table 5) (14). Methylation data for 57 of 95 sites were available in the 450K array. Among these, blood methylation correlated positively with methylation of 18 sites in adipose tissue and 3 sites in muscle after correcting for multiple testing (Table 2), including cg14013597 and sites annotated to AATK, CRMP1, ELOVL2, KCNQ2, MOGAT1, PGAM2, and SLC6A4, suggesting that methylation of some subgroup-unique sites may play a role in relevant tissues for T2D (Supplementary Fig. 6).

Table 2

Cross-tissue DNA methylation of sites included in subgroup-unique MRSs in different human tissues

			Blood—adipose tissue (n = 32)		Blood—skeletal muscle (n = 28)
CpG site	Subgroup	Annotated gene	r	q	r	q
cg05963087	SIDD	ENOX1	0.96	5.53e−14
cg14013597	SIDD		0.81	2.18e−06	0.78	9.82e−06
cg23616741	SIDD	PGAM2	0.82	2.18e−06
cg25356393	SIDD	KCNQ2	0.73	0.0001
cg16867657	MOD	ELOVL2	0.74	0.0001
cg13379325	SIDD	KCNQ2	0.72	0.0001
cg22891868	SIDD	MOGAT1	0.69	0.0003
cg13930790	SIDD	NCLN	0.68	0.0005	0.56	0.016
cg15081033	SIDD	SPSB4	0.64	0.001
cg26161329	MOD	PPM1E	0.62	0.002
cg16276209	SIDD	AATK	0.59	0.004
cg14578612	SIDD	KCNAB2	0.59	0.004
cg02789526	SIDD	A2BP1	0.57	0.007
cg01542019	MOD	TECR	0.56	0.007
cg07963234	MARD	CRMP1	0.52	0.016
cg14692377	MOD	SLC6A4	0.49	0.025
cg06933824	MOD	NEURL1B	0.49	0.025
cg15225267	SIDD		0.49	0.027	0.72	8.31e−05

Correlations between DNA methylation of sites included in subgroup-unique MRSs in blood and DNA methylation of these sites in adipose tissue and skeletal muscle taken from the same subjects for these cell types from the Monozygotic Twin Cohort based on FDR <5% (q < 0.05). Pearson correlation tests show significant correlations between DNA methylation of sites in blood and skeletal muscle and adipose tissue, respectively, for the subgroup-unique sites included in the subgroup-unique MRSs. A FDR analysis based on Benjamini-Hochberg was performed, and FDR <5% (q < 0.05) was considered significant. DNA methylation of 57 of the 95 subgroup-unique sites included in any of the MRSs was available from the 450K array and used to analyze DNA methylation in blood, muscle, and adipose tissue in subjects from the Monozygotic Twin Cohort. For 32 subjects, methylation data were available for blood and adipose tissue, and for 28 subjects, methylation data were available for blood and skeletal muscle. Here, DNA methylation in blood correlated positively with DNA methylation in adipose tissue of 18 sites and in skeletal muscle of 3 sites (q < 0.05).

Conclusions

We demonstrate for the first time that there are epigenetic differences between the novel T2D subgroups SIDD, SIRD, MOD, and MARD already at diagnosis. Importantly, epigenetic markers differed between the four subgroups in an independent validation cohort, further establishing a clear heterogeneity of these T2D subgroups. The subgroup-unique epigenetic markers did also associate with future diabetic complications, supporting development of blood-based epigenetic biomarkers for precision medicine of diabetes.

T2D is a heterogeneous disease with individual variation in obesity, insulin resistance, insulin secretory defects, and/or age between different patients. There is a need to improve treatment strategies for T2D and to better identify individuals with increased risk for complications. We recently reclassified T2D into four novel subgroups based on age at onset of diabetes, BMI, HbA1c at diagnosis, HOMA2-B, HOMA2-IR, and GAD autoantibodies (4). We also demonstrated that genetic risk scores for diabetes-related traits associate with these subgroups (10). Interestingly, here, we identified subgroup-unique epigenetic modifications in patients with newly diagnosed T2D, further supporting different underlying etiopathological processes for each subgroup. The fact that ∼50% of genes annotated to the subgroup-unique sites included in the MRSs have been previously associated with diabetes, insulin secretion, insulin resistance, obesity, and/or age, suggests that these epigenetic modifications have important functions in the pathogenesis of diabetes as well. This is further supported by methylation in blood of some subgroup-unique sites mirroring the methylation pattern in adipose tissue and muscle, two relevant tissues for T2D. One needs to consider that these are “only” correlations and based on 450K arrays that miss methylation data for several sites included in the subgroup-unique MRS. Future studies should further examine whether there are epigenetic differences between subgroups also in adipose tissue, muscle, liver, and islets. Nevertheless, our data point to adipose tissue being important, as 18 subgroup-unique methylation sites showed positive correlations between blood and adipose tissue, and interestingly, 12 of these are SIDD-unique sites, while 5 sites are MOD-unique sites, suggesting a potential role for methylation in adipose tissue in SIDD- and MOD-specific pathogenesis.

To better understand the biology of our findings, we performed a systematic literature search using each gene symbol annotated to the differentially methylated sites included in the subgroup-unique MRSs, and diabetes, insulin secretion/β-cell function, insulin resistance, obesity, or age. For example, SIDD is characterized by β-cell dysfunction, poor metabolic control, and higher HbA1c. We found two SIDD-MRS sites annotated to TXNIP, encoding a thioredoxin-binding protein associated with oxidative stress and regulation of β-cell function, and TXNIP methylation is associated with HbA1c and incident T2D, and is altered in muscle and islets from individuals with T2D (13,18,24,25). Moreover, methylation of CTDSPL, a gene annotated to SIDD-MRS sites, was associated with incident T2D and islet dysfunction (26). Interestingly, we found that higher SIRD-MRS was associated with an increased risk of CKD with certain SIRD-MRS sites annotated to RAB27B, whose downregulation contributes to exosome dysfunction in diabetic kidney disease (27). MOD-MRS sites were annotated to genes previously associated with obesity, including SLC6A4, TFEB, and ELOVL2 (28–30). SLC6A4 encodes a serotonin transporter regulating energy balance, and methylation of SLC6A4 has been linked to obesity (28). The transcription factor encoded by TFEB regulates lysosomal biogenesis and autophagy and is upregulated in adipose tissue during obesity (29), while ELOVL2 regulates synthesis of very long polyunsaturated fatty acids (30). MOD is characterized by higher BMI, but also by younger age, and notably, genes annotated to MOD-MRS sites are associated with aging, including ELOVL2, IFIH1, NAV2, TFEB, and TRIM59 (31–35). Two genes annotated to MARD-MRS sites have been linked to aging, including CRMP1, whose expression decreased with age (36), and RNF170 associated with age-dependent gait abnormalities (37). These results clearly show that the novel sites included in subgroup-unique MRSs identified in this study are linked to diabetes and subgroup-defining phenotypes.

Patients with T2D have a higher risk of CVD, CKD, and/or retinopathy (38). T2D subgroups also have different risks of these complications, confirming variability in patients with T2D and justifying reclassification of diabetes (4). CVD is responsible for the majority of deaths among patients with diabetes, while CKD is the leading cause of hospitalization in individuals with diabetes (38), reflecting the need for reliable prediction tools in this at-risk population. Interestingly, the identified subgroup-unique MRSs were associated with future risks of developing diabetes complications. Higher SIRD-MRS and MARD-MRS were associated with an increased risk of developing CVD and CKD, whereas a higher MOD-MRS was associated with lower risk for developing these two complications. Our results are supported by previous studies, where patients with SIRD had a higher risk of CKD and patients with MOD had a lower risk of CKD and coronary events (4,11). In the current study, we also found that the associations between patients with MOD and those with MARD and the risk for CVD and CKD were similar to the associations observed for MOD-MRS and MARD-MRS. On the other hand, SIRD-MRS and SIDD-MRS were able to predict CVD and/or CKD, whereas patients with SIRD and those with SIDD were not associated with the risk of these diseases, suggesting a potentially better ability of these MRSs on predicting complications compared with the T2D subgroups. Of note, our study compared one subgroup versus all remaining patients, while MARD has previously been used as a reference group (4). Overall, a new classification of T2D, now supported by epigenetic markers, could identify patients at high risk of developing complications already at diagnosis. This tool may be further developed to decrease suffering for patients and costs for the society.

This study has some limitations. The first discovery step identified a large number of methylation sites, which may have resulted in an inflated type I error rate, and was also supported by a λ of 1.7. It should, however, be noted that λ is not suitable to measure inflation in epigenome-wide associations studies (39). Nevertheless, we mitigated this potential issue by studying MRSs in three different cohorts including external validation, thus reducing possible bias and showing the reliability of the methylation sites included in the MRSs. At this point, we cannot fully conclude that the subgroup-unique sites presented here are the optimal combination of sites to differentiate between the individual subgroups and/or predict complications. Nevertheless, these MRSs seem to be robust since they remain significantly different between subgroups after adjusting for confounding factors, such as cell composition, sex, and clinical variables defining the subgroups, and could discriminate between the subgroups in an independent validation cohort. Moreover, taking one subgroup at a time as a reference group showed that the MRSs robustly differ between that specific subgroup from all the other subgroups, demonstrating the capability of the MRSs at characterizing each respective subgroup.

Importantly, all subgroup-unique MRSs, except for SIDD-MRS, were strongly associated with their respective subgroup, independently of adjustment for the corresponding subgroup-defining variable. SIDD-MRS was associated with a twofold higher probability of being SIDD, but this result was partially driven by HbA1c, which is the identifying variable of this subgroup. In contrast, similar effect sizes were observed for an association between MARD-MRS and MARD after adjusting for age and between SIRD-MRS and SIRD after adjusting for HOMA-IR. Notably, the effect size was even greater when the clinical variable was included in the model for the association between SIRD-MRS and SIRD after adjustment for HOMA-B, and the association between MOD-MRS and MOD after adjusting for BMI. This suggests that our epigenetic markers are equally robust at discriminating these subgroups as the identifying trait of the subgroup, further supporting that epigenetics play a key role in the etiology of these subgroups.

Owing to multicollinearity, analyses for complications were not adjusted for the subgroup-defining phenotypes. Therefore, we cannot rule out that the associations between MRSs and complications are influenced by these phenotypes. Hence, these associations may be both due to altered DNA methylation and clinical phenotypes, a conundrum difficult to disentangle. However, only age and HOMA2-B were associated with future CVD and CKD, suggesting a minor effect for the rest of subgroup-defining phenotypes. Nevertheless, age could have influenced some of the associations between MRSs and diabetic complications, especially for MOD-MRS and SIDD-MRS, associated with a low risk for CVD and/or CKD. Individuals with MOD and SIDD were younger and therefore less likely of developing events after 8 years of follow-up. Certain risk factors for diabetic complications (e.g., smoking, LDL, albumin-to-creatinine ratio) were completely or partially missing at baseline for patients included in this study, and we could therefore not assess whether MRSs affected the predictive value of such factors. Limited sample size and the short period of follow-up may affect associations with the studied complications. Moreover, future epigenetic studies focusing only on diabetic complications may identify additional markers to be included in scores for prediction of CVD, CKD, and retinopathy. However, such analyses are not suitable in this study since they are unlikely to support reclassification of diabetes, the main goal of this study. Since this study was based on predominantly northern European subjects, the applicability to other ethnicities needs to be evaluated to establish the generalizability of associations between the identified methylation patterns and SIDD, SIRD, MOD, and MARD.

Conclusion

A reclassification of T2D might help tailor prevention strategies to individual subgroups and personalize care for those affected by diabetes. The observed associations between DNA methylation of certain sites and the novel T2D subgroups support the hypothesis that subgroups have epigenetic differences. Our study also suggests that epigenetic mechanisms may be more important for some of the subgroups, which should be further dissected in future studies. Identified epigenetic markers could successfully discriminate between diabetes subgroups and associated with future incidence of CVD and CKD, thus further validating the reclassification of T2D on an epigenetic level. These epigenetic markers may be developed for precision medicine to improve treatment of T2D subgroups and prevent their complications.