DNA methylation age calculators reveal association with diabetic neuropathy in type 1 diabetes

Background Many CpGs become hyper or hypo-methylated with age. Multiple methods have been developed by Horvath et al. to estimate DNA methylation (DNAm) age including Pan-tissue, Skin & Blood, PhenoAge, and GrimAge. Pan-tissue and Skin & Blood try to estimate chronological age in the normal population whereas PhenoAge and GrimAge use surrogate markers associated with mortality to estimate biological age and its departure from chronological age. Here, we applied Horvath’s four methods to calculate and compare DNAm age in 499 subjects with type 1 diabetes (T1D) from the Diabetes Control and Complications Trial/Epidemiology of Diabetes Interventions and Complications (DCCT/EDIC) study using DNAm data measured by Illumina EPIC array in the whole blood. Association of the four DNAm ages with development of diabetic complications including cardiovascular diseases (CVD), nephropathy, retinopathy, and neuropathy, and their risk factors were investigated. Results Pan-tissue and GrimAge were higher whereas Skin & Blood and PhenoAge were lower than chronological age (p < 0.0001). DNAm age was not associated with the risk of CVD or retinopathy over 18–20 years after DNAm measurement. However, higher PhenoAge (β = 0.023, p = 0.007) and GrimAge (β = 0.029, p = 0.002) were associated with higher albumin excretion rate (AER), an indicator of diabetic renal disease, measured over time. GrimAge was also associated with development of both diabetic peripheral neuropathy (OR = 1.07, p = 9.24E−3) and cardiovascular autonomic neuropathy (OR = 1.06, p = 0.011). Both HbA1c (β = 0.38, p = 0.026) and T1D duration (β = 0.01, p = 0.043) were associated with higher PhenoAge. Employment (β = − 1.99, p = 0.045) and leisure time (β = − 0.81, p = 0.022) physical activity were associated with lower Pan-tissue and Skin & Blood, respectively. BMI (β = 0.09, p = 0.048) and current smoking (β = 7.13, p = 9.03E−50) were positively associated with Skin & Blood and GrimAge, respectively. Blood pressure, lipid levels, pulse rate, and alcohol consumption were not associated with DNAm age regardless of the method used. Conclusions Various methods of measuring DNAm age are sub-optimal in detecting people at higher risk of developing diabetic complications although some work better than the others.


Background
CpGs are regions of DNA where a cytosine is followed by a guanine nucleotide. Cytosines within CpGs can be methylated, and CpG methylation levels affect gene expression. Many CpGs become hyper or hypo-methylated with age [1][2][3][4]. In 2013, Horvath used publicly available DNA methylation (DNAm) data to define and evaluate a DNAm age predictor, Pan-tissue, which is accurate across most tissues and cell types. Chronological age was regressed on CpG methylation levels using a penalized regression model (elastic net) which selected 353 CpGs [1] (Supplementary Table 1). Pan-tissue has been widely used and has shown that faster epigenetic aging is associated with multiple age-related diseases and conditions (e.g., Alzheimer's, cancer, cardiovascular diseases (CVD)) indicating that epigenetic age is an indicator of health status. Some risk factors for type 2 diabetes (T2D) including BMI, waist circumference, and fasting glucose have been associated with higher Pan-tissue epigenetic age acceleration (EAA = epigenetic age − chronological age) [5,6]. On average liver Pan-tissue EAA increased significantly by 0.33 years per BMI unit [7].
Pan-tissue performed sub-optimally estimating fibroblast age in in vitro studies. Therefore, Horvath et al. described Skin & Blood DNAm age using a similar method to Pan-tissue which performed remarkably well across a wide spectrum of cells that are most frequently used in in vitro studies, including the blood (Supplementary Table 1). Skin & Blood EAA was highly predictive of time to all-cause mortality. It was positively correlated with waist/hip ratio (WHR), blood insulin, glucose, triglyceride, systolic blood pressure, and BMI and negatively correlated with HDL and physical exercise. However, the respective correlation coefficients were weak (|r| < 0.11) [2].
Pan-tissue and Skin & Blood were both developed using chronological age as a surrogate for biological age. Therefore, they may not capture CpGs that signal departure of biological age from chronological age. In a newer method called PhenoAge, chronological age was replaced with a surrogate measure of phenotypic age developed using clinical data. A Cox penalized regression model was applied where the hazard of mortality was regressed on 42 clinical markers and chronological age to select variables for inclusion in the phenotypic age score. Nine clinical markers and chronological age were selected and used to estimate the 10-year mortality risk score which was then converted into units of years. Finally, the resulting phenotypic age estimate was regressed on CpG methylation levels using an elastic net regression model (Supplementary Table 1). A 1-year increase in PhenoAge was associated with a 4.5% increase in the risk of all-cause mortality in independent populations without diabetes. PhenoAge predicted mortality significantly better than Pan-tissue. It was also associated with increased risk of CVD and differed significantly between never, current, and former smokers. It was also positively correlated with blood insulin, glucose, triglyceride, and WHR and negatively correlated with HDL and physical exercise [8].
Most recently, Horvath et al. defined another DNAm age calculator called GrimAge. To develop this method, 88 plasma proteins and smoking pack-years were individually regressed on chronological age, sex, and the CpGs methylation levels using an elastic net regression model. Twelve plasma proteins and packyears had high correlations between their DNAm estimation and the corresponding measured levels. These twelve DNAm estimated plasma proteins and smoking pack-years as well as chronological age and sex were regressed on the hazard of aging-related mortality using a Cox penalized regression model. This selected DNAm pack-years, age, sex, and the predicted DNAm for 7 of the 12 plasma proteins. The combined estimate of these factors was then transformed into GrimAge which has the same mean and variance as chronological age (Supplementary Table 1). GrimAge was highly predictive of lifespan and time-to-CVD even after adjustment for known risk factors and outperformed Pan-tissue and PhenoAge. It was also associated with hypertension and T2D. In addition, GrimAge was correlated with BMI, WHR, and physical exercise; blood insulin, glucose, HbA1c, triglyceride, and HDL; and albuminuria (all correlations were in the expected directions) [9].
Only a small proportion of CpGs are common among three of the four epigenetic ages (Supplementary Figure  1), and there is only weak/moderate correlation among them [2,8,9].
To our knowledge, there has been no study of DNAm age in type 1 diabetes (T1D). However, telomere length, another indicator of aging, has been investigated in T1D, and it was found to be significantly shorter in T1D compared to non-diabetic subjects [10]. Shorter telomere length has also been associated with T1D duration, pulse pressure [10], BMI [11], systolic blood pressure, allcause mortality [12], and diabetic nephropathy [13] in subjects with T1D.
Here, we investigated DNAm age calculated by all four methods in 499 subjects with T1D from the Diabetes Control and Complications Trial/Epidemiology of Diabetes Interventions and Complications (DCCT/EDIC) study and its association with diabetic complications (CVD, nephropathy, retinopathy, and neuropathy) and their risk factors [14][15][16]. We also examined DNAm age in a smaller subset at two time points, 16-17 years apart, to investigate changes in DNAm age over time.

Results
Illumina whole blood EPIC data Comparison of the four DNAm ages with chronological age and with each other Characteristics of the subjects with EPIC DNAm data are summarized in Table 1. All four epigenetic ages were highly correlated with chronological age and with each other. However, there were significant difference among them: GrimAge was higher than Pan-tissue, and both were higher than chronological age whereas Skin & Blood and PhenoAge were both lower than chronological age, and PhenoAge was lower than Skin & Blood (GrimAge > Pan-tissue > chronological age > Skin & Blood > PhenoAge) (all p < 0.0001) ( Table 2, Fig. 1, Supplementary Figure 2).

Changes in EAA by chronological age
The differences between chronological age and each of Pan-tissue, Skin & Blood, and GrimAge decreased by 0.13, 0.07, and 0.14 years per 1-year increase in chronological age (p = 1.18E−4, 6.49E−3, and 1.03E−4, respectively). The difference between chronological age and PhenoAge did not differ significantly by chronological age (Fig. 2, Supplementary Table 2).

Association of the four DNAm ages with development of diabetic complications
Epigenetic ages were not significantly associated with the development of CVD or retinopathy. Although DNAm ages were not associated with estimated glomerular filtration (eGFR), both PhenoAge (β = 0.023, p = 0.007) and GrimAge (β = 0.029, p = 0.002) were positively associated with repeated measures of albumin excretion rate (AER, natural log transformed) which remained significant after adjustment of HbA1c levels. The effect of GrimAge on AER increased over time (GrimAge × EDIC year interaction β (SE) = 0.0013 (0.0005), p = 5.10E−3) whereas effect of PhenoAge was not significantly different over time (PhenoAge × EDIC year interaction p = 0.85) ( Table 3, Supplementary Figure 3-4). Pan-tissue, Skin & Blood, and PhenoAge were not associated with neuropathy, but GrimAge was positively associated with both diabetic peripheral neuropathy (DPN: OR = 1.07, p = 9.24E−3) and cardiovascular autonomic neuropathy (CAN: OR = 1.06, p = 0.011). These associations also remained significant after further adjustment for time-weighted HbA1c ( Table 3).

Association of risk factors of diabetic complication with the four DNAm ages
The univariable associations of different factors with EAAs are shown in Supplementary In the multivariable analysis when all factors were included, males had on average 1.5 years higher Pan-tissue (p = 8.00E−4) and GrimAge (p = 9.99E−5) compared to females whereas females had on average 1.5 years higher PhenoAge compared to males (p = 0.005) ( Table 4). PhenoAge increased 0.4 years per 1% increase in time-weighted HbA1c (p = 0.026) and 0.01 years per 1-month increase in T1D duration (p = 0.043) ( Table 4). The effect of time-weighted HbA1c on Pheno-Age was not significantly different in either conventional (β (SD) = 0.45 (0.29), p = 0.12) or intensive (β (SD) = 0.15 (0.48), p = 0.76) therapy group (interaction p = 0.56). Skin & Blood increased 0.09 years per one-unit increase in BMI (p = 0.048) ( Table 4). Those with strenuous physical activity at work on average had 1.99 years lower Pan-tissue compared to those with sedentary jobs (p = 0.045), and those who achieved one to two times the recommended level of physical activity during leisure time on average had 0.8 years lower Skin & Blood compared to those who did not achieve the recommended level (p = 0.022) ( Table 4). Current smokers had on average 7.1 years higher GrimAge compared to non-smokers (p = 9.03E−50). The other factors were not significantly associated with epigenetic ages (Table 4).

Illumina 450K whole blood data
Characteristics of the subjects with 450K data are summarized in Supplementary Table 4. Chronological age and all four epigenetic ages were highly correlated. However, there were significant differences among them: Pan-tissue > (greater than) GrimAge > Skin & Blood > chronological age > PhenoAge (all p < 0.0001) (Supplementary Table 5 Table 6).
The difference between chronological age and GrimAge decreased by 0.13 years per 1-year increase in chronological age (p = 0.01) (Supplementary Table 6 and Figure 9). Sex, cohort, treatment group, T1D duration, stimulated C-peptide, and time-weighted HbA1c were not significantly associated with epigenetic ages (Supplementary Table 6, Table 5).

Illumina 450K monocyte data
Chronological age and all four epigenetic ages were highly correlated (Supplementary Table 7, Supplementary Figure 8-9). However, except for Pan-tissue and GrimAge, there were significant differences among them: PhenoAge > Pan-tissue ≈ GrimAge > Skin & Blood > chronological age (p < 0.001). Epigenetic age acceleration (years) DNAm age (years) 40.3 (6.9) 40.3 (6.6) 40.7 (6.5) 42.2 (7.0) 40.8 (6.8) Epigenetic age acceleration (years) 5.7 (4.7) 5.5 (4.3) 5.7 (4.8) 6.6 (5.0) 5.9 (4.7) All factors were obtained at DNAm measurement except for stimulated C-peptide which is measured at DCCT eligibility MET metabolic equivalent of task *Time-weighted HbA1c since DCCT baseline †Level of activity on the job, at school, or in home making: sedentary such as office work with occasional inter-office walking; moderate activity requires considerable but not constant lifting, walking, bending, pulling, etc. such as homemaker with family and without domestic assistance; and strenuous activity requires almost constant lifting, bending, pulling, scrubbing, etc. such as furniture mover ‡According to the international classification by Ainsworth used by American College of Sports Medicine (ACSM), light, moderate, hard, and very hard activity was allocated 3, 4, 6, and 9 METs, respectively. For each participant, these allocated MET values were multiplied by the time (min) spent in that activity to obtain the MET for that level of activity. The sum of METs from all activities was recorded as the total leisure time activity for each participant. Subjects then were categorized into three groups based on the ACSM recommendation for METs min/week [17,18] Although Pan-tissue EAA was on average 2.9 years lower in the former DCCT intensive versus conventional treatment group, the difference was not significant in the multivariable analysis. The other EAAs were also not significantly different between the two treatment groups ( Table 5, Supplementary Table 8).
The differences between chronological age and each of Pan-tissue, Skin & Blood, and GrimAge decreased by 0.3, 0.2, and 0.2 years per 1-year increase in chronological age (p = 0.005, 0.003, and 0.004, respectively) (Supplementary Table 8 and Figure 10).
All four EAAs were highly correlated between whole blood and monocyte (Fig. 3).
In the multivariable analysis, Pan-tissue on average decreased 13.05 years per 1 pmol/ml increase in stimulated C-peptide at DCCT eligibility (p = 0.016). Association of C-peptide with Pan-tissue was not significantly different in the two treatment groups (p = 0.61). Sex, cohort, treatment group, T1D duration, and time-weighted HbA1c were not significantly associated with epigenetic ages (Table 5).

Discussion
We used DNAm data in whole blood measured by Illumina EPIC array and four different methods to estimate epigenetic ages in 499 subjects with T1D. Subsequently, we compared estimated epigenetic ages with chronological age and investigated if the epigenetic ages were associated with development of T1D complications (CVD, nephropathy/decreased renal function, retinopathy, and neuropathy) and their risk factors. All four epigenetic ages were correlated with chronological age and with each other, but there were significant differences between them. Pan-tissue and Skin & Blood were developed to predict chronological age in healthy individuals whereas PhenoAge and GrimAge used biological biomarkers associated with time-to-death to predict differences in life expectancy of individuals. In addition, they were developed and tested in datasets measured in different tissues by different arrays (Illumina 27K, 450K, and EPIC) and used different statistical methods. As a result, there have been only low to moderate correlations between them [2,8,9]. Pan-tissue and GrimAge were significantly higher than chronological age, but Skin & Blood and PhenoAge were significantly lower than chronological age. T1D usually has a negative impact on general health [19,20]. Therefore, we expected DNAm age to be higher than chronological age in subjects with T1D [21]. However, DCCT subjects were a relatively healthy group of subjects with T1D at baseline due to the extensive inclusion/exclusion criteria applied. Out of~7000 individuals who made initial contact, only 1441 subjects aged 13-39 years with 1-15 years of T1D and no serious long-term complications of diabetes were included in DCCT. Subjects were excluded if they were at risk for adverse effects (e.g., history of frequent ketoacidosis, hypoglycemic coma, or seizure) had known risk factors for vascular complications, were unlikely to comply with the demands of treatment protocols, did not demonstrate an adequate understanding of the DCCT's purpose, or had drug addiction, chronic alcoholism, or major mental illness [22]. In addition, individuals with EPIC DNAm data are not an entirely random sample of DCCT/EDIC subjects (Supplementary Tables 9-13).
With the exception of PhenoAge, the difference between the epigenetic ages and chronological age decreased in older subjects. This is consistent with previous findings where longitudinal changes in Pantissue EAA were tracked using linear mixed models (LMMs) within five different cohorts and showed that Fig. 2 Epigenetic age acceleration vs. chronological age in EPIC data. The dash line is epigenetic age acceleration = 0. The solid line is the line fitting linear regression model with chronological age as predictor and epigenetic age acceleration as outcome and is present when chronological age is significantly associated with epigenetic age acceleration (p < 0.05) We did not find significant association between epigenetic ages and development of CVD. Studies have investigated the association of Pan-tissue and risk of developing CVD in non-diabetic subjects. In 2543 African Americans, the hazard ratio of fatal coronary heart disease increased by 1.03 per year increase in Pan-tissue over a ≈ 17-year follow-up period (144 events) [21]. Similarly, in a cohort of 1863 subjects aged 50-75 years from Germany, the risk of CVD mortality increased by 1.04 per year increase in Pan-tissue over 13 years of follow-up (194 events) [24]. However, our study had only 22% power to detect an effect with this size, due to the relatively small number (N = 58) of CVD events. Another study did not find any significant association between Pan-tissue and incidence of CVD in women from African, Caucasian, and Hispanic ancestry [25]. This discrepancy could be due to race, gender, age, etc. differences between the studied populations as well as statistical power. Increases in PhenoAge and GrimAge have also been reported in association with increased risk of CVD in non-diabetic subjects (β = 1.10 and HR = 1.07, respectively) [8,9]. However, our study was again under-powered to detect these effects (power = 0.55 and 0.72, respectively). It is noteworthy that in all these studies including the current study, Pan-tissue was calculated in whole blood which is not the target tissue for diabetic complication. Although it has been shown that DNAm profiles are quite similar in different tissues, there are tissue-specific differentially methylated regions which can affect DNAm calculation [26]. However, accessing target tissues are not feasible especially in large scale studies of living subjects.
To our knowledge, association of DNAm age with risk of developing nephropathy, retinopathy, or neuropathy has not been investigated before. We did not observe any significant association between DNAm age and development of nephropathy and retinopathy. However, we found higher GrimAge and PhenoAge to be associated with higher levels of repeated measures of AER, an indicator of decreasing renal function. These associations remained significant even after adjusting for timeweighed HbA1c indicating that at least part of the effects is independent of HbA1c levels. This result is consistent with previous findings where GrimAge has been associated with albuminuria in non-diabetic subjects [9]. Of the four epigenetic ages, higher GrimAge was also associated with development of both DPN and CAN. Larger number of events and higher statistical power along with relatively larger effect sizes could be among the reasons that we detected association of GrimeAge with neuropathy but not the rest of the complications, although CSME has similar number of cases.
Intensive therapy and keeping HbA1c levels close to the normal range has been associated with decreased risk of developing diabetic complications [27]. In our study, both time-weighted HbA1c (but not treatment group) and T1D duration were associated with higher PhenoAge but not with the other epigenetic ages.
Males had higher Pan-tissue and GrimAge compared to females consistent with previous findings [9,21,25] whereas PhenoAge was higher in women and sex was not associated with Skin & Blood. All Pan-tissue and PhenoAge CpGs are on autosomal chromosome. Only The majority of known risk factors for diabetic complications were not associated with the four epigenetic ages. Two studies have investigated the association of CVD risk factors with Pan-tissue in non-diabetic subjects, but they reported conflicting results [21,25]. In our study, individuals with physically strenuous jobs had significantly lower Pan-tissue compared to individuals with sedentary jobs, and those who achieved one to two times the recommended level of physical activity during leisure time had lower Skin & Blood compared to those who did not achieve the recommended level. Consistent with these findings, a twin study found that although only a small amount of variance in Pan-tissue is explained by non-shared environmental factors in younger individuals, leisure time physical activity can affect Pantissue during adult years [28]. Higher BMI was associated with higher Skin & Blood consistent with what has been observed in non-diabetic subjects before [2]. Current smoking had a large effect on GrimAge; GrimAge increased on average 7 years in current smokers compared to non-smokers. This finding was expected since pack-years were one of the surrogate biomarkers used to generate GrimAge [9]. Triglyceride, HDL, BMI, and physical exercise have been correlated with PhenoAge and/or GrimAge [8,9]. In our study, although some of these factors were associated with Phe-noAge and/or GrimAge in the univariable analysis, these associations did not remain significant in the multivariable analysis [22]. Cell counts (B, CD4T, CD8T, natural killer, eosinophil, and monocyte) and batch (as a categorical variable) were also included in the multivariable analysis. All factors were obtained at DNAm measurement except for stimulated C-peptide which is measured at DCCT eligibility MET metabolic equivalent of task *Time-weighted HbA1c from DCCT baseline to DNAm measurement † Level of activity on the job, at school, or in home making: sedentary such as office work with occasional inter-office walking; moderate activity requires considerable but not constant lifting, walking, bending, pulling, etc. such as homemaker with family and without domestic assistance; and strenuous activity requires almost constant lifting, bending, pulling, scrubbing, etc. such as furniture mover ‡According to the international classification by Ainsworth used by American College of Sports Medicine (ACSM), light, moderate, hard, and very hard activity was allocated 3, 4, 6, and 9 METs, respectively. For each participant, these allocated MET values were multiplied by the time (min) spent in that activity to obtain the MET for that level of activity. The sum of METs from all activities was recorded as the total leisure time activity for each participant. Subjects then were categorized into three groups based on the ACSM recommendation for METs min/week [17,18] We also investigated epigenetic age in a smaller subset (N = 63) using Illumina 450K array in whole blood and 16-17 years later in monocytes which gave us the opportunity to investigate the difference between the two arrays and the change in epigenetic age over time. In 450K whole blood data, on average, Pan-tissue had the highest estimated value whereas in EPIC data GrimAge was estimated higher than the other epigenetic ages. In addition, Skin & Blood, which was on average lower than chronological age in EPIC data, was on average higher than chronological age in 450K data. Therefore, it appears that the type of array can affect the estimated epigenetic age. However, this could also be due to subjects being highly selected (two extremes of HbA1c and complications risk) in 450K data. Epigenetic ages in whole blood and monocyte measured 16-17 years apart were always correlated indicating that DNAm age is consistent over time and in multiple cell types. This is consistent with previous findings showing that Pan-tissue is consistent across life span [23], and a substantial amount (over 70%) of its changes are shared between different tissue/ cell types [29] and also PhenoAge being correlated in different tissues and cell types [8]. In monocyte 450K data, we found a new association which we did not observe in EPIC data: higher stimulated serum C-peptide at DCCT eligibility was associated with lower monocyte Pan-tissue measured decades later at EDIC year 16-17. The association was in the expected direction since preserved beta cell function as measured by stimulated C-peptide is associated with better clinical outcomes (i.e., better glycemic control [30] and lower risk for diabetic complications [31][32][33][34]). The fact that this association was not observed in EPIC data with a much larger sample size could partly be due to the fact that the two sample populations are different: the 450K sample was selected from two extremes of HbA1c and complications risk, whereas the EPIC sample was selected randomly from each cohort/treatment group. In addition, probes may perform slightly differently between the two arrays (450K and EPIC), and as a result measured levels of methylation can be different [35]. In addition, whole blood is a mixture of different cells with different halflives which could dilute the association(s) of individual cell types especially if they make up only a small proportion (e.g., monocytes, 2-8%).
We investigated the physical distance between all CpGs that are included in Pan-tissue, Skin & Blood, and PhenoAge epigenetic age calculations and T1D GWAS loci [36,37]: all CpGs are > 25 Mb away from them. Therefore, it is unlikely that estimated DNAm ages were confounded by methylation levels of the CpGs associated with T1D. However, the CpGs included in GrimAge are not publicly available. There have been two epigenomewide association studies of T1D [38,39]; however, associated CpGs (N = 132) are available for one of them [38]. Of these, only two CpGs are common with Pantissue (cg02047577 (Chr20, 62,587,702 (HG19)) and cg16494477 (Chr5, 170,847,251)): the later CpG and Cell counts (B, CD4T, CD8T, natural killer, eosinophil, and monocyte) were also included in the multivariable whole blood analysis *Stimulated C-peptide at DCCT eligibility †Time-weighted HbA1c from DCCT baseline to monocyte DNAm measurement another site are in common with PhenoAge (cg11903057 (Chr4, 40,198,776) and cg16494477 (Chr5, 170,847,251)). Therefore, it appears unlikely that they can have major effect on DNAm age.

Conclusions
We found that although all four epigenetic ages are correlated with each other and with chronological age, however, there are significant differences between them in subjects with T1D. We also found that DNAm age is consistent over time, but type of array (450K vs. EPIC) and cell type can affect the estimated epigenetic age. None of the epigenetic ages were associated with CVD or retinopathy, but PhenoAge and GrimAge were both associated with decreasing renal function as measured by AER. GrimAge was also associated with both DPN and CAN. Only PhenoAge was positively associated with HbA1c levels and T1D duration, two major risk factors for diabetic complications. Some of the other risk factors of diabetic complications were associated with individual epigenetic ages. Therefore, it seems that the investigated epigenetic ages all work sub-optimally in detecting subjects with T1D who are at higher risk to develop complications. However, PhenoAge and specifically GrimAge performed better that Pan-tissue and Skin & Blood suggesting that including biomarkers associated with aging-related mortality improves the accuracy of DNAm measurement. Nevertheless, only some of the risk factors of diabetic complication which are also among the main factors associated with aging-related mortality in the general population were considered in their development (serum creatinine and glucose in Phe-noAge and smoking pack-years in GrimAge), and major factors such as hypertension, lipid levels, BMI, and HbA1c which is a better indicator of glycemic levels in long-term compared to serum glucose were not considered [8,9]. Including these factors could potentially improve epigenetic age estimation in the general population and specifically in subjects with T1D.

Subjects
The study subjects were from the DCCT/EDIC study. Subjects with T1D aged 13-39 years were recruited into DCCT in 1983-1989 in two cohorts. The primary prevention cohort included participants with 1-5 years of diabetes and no pre-existing retinopathy or nephropathy. The secondary cohort included participants with 1-15 years of diabetes and pre-existing mild retinopathy. Subjects were randomly assigned to receive intensive or conventional diabetes therapy [40]. The DCCT ended in 1993, and subjects subsequently have been followed annually through the EDIC study (Supplementary Figure 11).
Genome-wide DNAm measurement, QC, and normalization  [16][17]. These included 32 subjects from the conventional treatment group with mean DCCT HbA1c > 9.1% (76 mmol/mol) and significant progression of retinopathy and/or nephropathy from the DCCT closeout to EDIC year 10, and 31 subjects from the intensive treatment group with mean DCCT HbA1c < 7.3% (56 mmol/mol) and no development of retinopathy and nephropathy until EDIC year 10 [41]. Three subjects had missing methylation data for monocytes including two subjects from the conventional and one subject from the intensive treatment group. The R package meffil (https://github.com/perishky/ meffil; accessed on December 2018) [42] was used to perform QC and normalization. Samples were removed if their predicted sex based on DNAm did not match their recorded sex or had > 10% probes with detection p value > 0.01 or > 10% probes with bead number < 3. Samples were also excluded if their SNP genotypes did not match with those from GWAS-array (concordance threshold = 0.8). Illumina HumanCoreExome BeadArray (Illumina, San Diego, CA, USA) data imputed to 1000 Genomes (phase 3, v5) was used for this comparison [43]. One sample from EPIC and one sample from 450K monocyte data did not pass QC. Functional normalization ("noob" for dye-bias and background correction followed by "quantile" normalization implemented in meffil) was then performed to account for technical variation in the data [42]. Blood cell proportions were estimated using the method [44] implemented in meffil with "blood gse35069 complete" as reference [42].

T1D complications
CVD was described as any CVD from DNA collection date to EDIC follow-up year 20 [45].

Statistical analysis
Spearman correlation was used to test for correlation between chronological age and DNAm age and between the four epigenetic ages. Paired sample T tests were used to determine if epigenetic ages were significantly different from chronological age and if the four epigenetic ages differed significantly.
EAA was calculated by subtracting chronological age from epigenetic age (EAA = epigenetic age − chronological age). EPIC array was performed in seven batches. Therefore, batch was included in all multivariable analyses of EPIC data. Since whole blood is a combination of 7 different cell types (neutrophil, B cell, CD4T, CD8T, natural killer cell, eosinophil, and monocyte) and their proportions affect DNAm and epigenetic ages, six predicted cell proportions were included as covariates in all multivariable analyses of whole blood. Neutrophils were excluded as the seven cell proportions sum to one.
Cox proportional hazard models were used to test association of DNAm age with development of CVD, nephropathy, and retinopathy using EPIC data. Logistic regression was used to test association of DNAm age with development of neuropathy during EDIC. Subjects who developed complications during DCCT (before DNAm measurement) were excluded from their corresponding analysis. We also investigated association of DNAm age with repeated measures of renal function, eGFR (annual), and AER (biannual, natural log transformed) from DCCT closeout to EDIC follow-up year 18 using LMMs. Two models were fit for both Cox and LMMs. Model 1 included sex, age, and T1D duration at the time of DNAm measurement along with batch and cell proportions. Model 2 included all covariates in model 1 plus repeated cross-sectional measures of HbA1c (Supplementary Table 12).
In univariable analysis, linear regression was used to test the association of different factors with EAA one at a time. In multivariable analysis, all factors plus chronological age were included in the model, and their associations were tested with DNAm age using linear regression. These factors included sex, cohort, treatment group, and stimulated C-peptide at DCCT eligibility as well as T1D duration, time-weighted HbA1c, BMI, systolic and diastolic blood pressure, HDL, LDL, triglyceride, total cholesterol, pulse rate, current smoking, drinking status, physical activity during work, and leisure time at the time of DNAm measurement (Supplementary Table 3). Since HbA1c and treatment group are highly associated, only HbA1c was included in the multivariable analysis of EPIC data. Due to small sample size and being highly selected on multiple traits, associations of only a subset of factors (sex, age, T1D duration, stimulated C-peptide, HbA1c, cohort, and treatment group) were tested in 450K data.
All the statistical analyses were performed using SAS 9.4 (Cary, NC).