The mediating role of DNA methylation clocks in associations of race, ethnicity, education, income, and occupation with mortality: findings from NHANES 1999-2002

Descriptive characteristics and survival analysis

The unweighted sample statistics for this study sample are presented in Table 1. Among the 2,402 participants in the analysis, 947 deaths (39%) were observed up to Dec 31, 2019. Of these, 164 were attributed to cardiovascular-related causes, while 197 were due to cancer. The racial and ethnic composition of the sample included 39% White Americans, 22% Black Americans, 36% Hispanic, and 3.5% individuals identifying as another race or mixed race. In terms of socioeconomic characteristics, 46% of participants had not completed high school, 38% had household income between two and five times of the poverty thresholds, and 42% reported their occupation as semi-routine blue-collar. The sample was 49% female and 25% foreign-born. The age of the participants ranged from 50 to 84 years, with the mean age as 65 years (SD=9).

The correlations between all the mediators we examined are shown in Figure 1. Most clocks had strong correlations with chronological age (r>0.7), with somewhat more moderate correlations for Weidner and Vidal-Bralo (r=0.5 and r=0.6, respectively), and Horvath predicted telomere showed moderate negative association with chronological age as the correlation coefficient was -0.6. Only DunedinPoAm and Yang epigenetic age showed weaker correlations with chronological age, r=0.04 and r=0.2, respectively.

Figure 1. Raw correlation matrix of potential mediators and age in NHANES 1999–2002 (n = 2,402). Mediators include chronological age, physiological functioning clocks, pace-of-aging clock, clocks trained on other aging biomarkers, as well as behavioral, anthropometric, and clinical variables.

The estimated coefficients from the survival analysis models (without any mediators) are presented in Supplementary Table 1. Compared to White participants, Black and Hispanic participants exhibited significantly higher risks of all-cause mortality, with hazard ratios (HR) of 1.51 (95%CI: 1.25-1.82) and 1.40 (95%CI: 1.12-1.75), respectively, after adjusting for gender, age, age², and nativity. All categories of lower educational attainment and lower household income were associated significantly with higher risks on all-cause mortality compared to the most advantaged groups, with HRs ranging from 1.63 (95%CI: 1.24-2.14) to 3.00 (95%CI: 2.18-4.12) with the most disadvantaged groups showing the biggest ratios for mortality. Individuals in high-skilled blue-collar and semi-routine blue-collar occupations showed significantly higher all-cause mortality risk compared to white-collar and professionals (HR=1.64, 95%CI: 1.19-2.27 and HR=1.55, 95%CI: 1.13-2.12, respectively).

The survival analysis of DNA methylation measures and clinical mediators and mortality are presented in Supplementary Table 2. Among the chronological age clocks, most of the clocks showed significant associations with all-cause mortality, and GrimAge, GrimAge2, and DunedinPoAm showed the strongest effect sizes with 1.11 (95%CI: 1.09-1.13), 1.10 (95%CI: 1.09-1.12), and 43.7 (95%CI:16.9-113), respectively. Weidner and SkinBlood epigenetic age showed less strong associations, and Yang epigenetic age and Horvath predicted telomere showed null effects. Only GrimAge, GrimAge2, and DunedinPoAm showed significant associations with both cardiovascular-specific and cancer-specific mortality due to smaller sample size. In addition, Hannum, Horvath, Vidal-Bralo, and PhenoAge epigenetic age showed significant associations with cancer-specific mortality with a HR range of 1.04 (95%CI: 1.01-1.07) to 1.05 (95%CI: 1.02-1.07).

The mediating roles of DNA methylation measures and clinical factors in the association between race and ethnicity and all-cause mortality

The mediation analysis of DNA methylation clocks and clinical biomarkers for the association between race and ethnicity and all-cause mortality are shown in Figure 2 and Supplementary Table 3. Multiple mediation signals were observed when comparing Black and White participants, as well as for the comparison between Hispanic and White participants. Physiologically trained clocks GrimAge, GrimAge2, and DunedinPoAm showed moderate mediation. Among Black vs. White participants comparisons, GrimAge, GrimAge2, and DunedinPoAm mediated a proportion of 34% (95%CI: 19%-67%), 46% (95%CI: 27%-89%), and 24% (95%CI: 10%-49%), respectively. Unexpectedly, all chronological age clocks except Lin epigenetic age exhibited significant negative mediation ranging from -3% (95%CI: -8%- -0.3%) to -19% (95%CI: -38%- -10%), and Horvath predicted telomere exhibited -44% (95%CI: -79%- -27%) mediation. These results indicated that although compared to White participants, Black participants had overall higher risk of mortality, they also had on average higher Horvath predicted telomere values. These higher values of DNA methylation predicted telomere length are associated with reduced risks for mortality and protective for Black participants to show even greater mortality disparities. Comparing Hispanic to White participants, GrimAge2 and DunedinPoAm pace of aging exhibited significant positive mediation of 39% (95%CI: 18%-101%) and 28% (95%CI: 13%-70%), respectively. Most chronological age clocks did not exhibit significant mediation effects, except for Vidal-Bralo, which showed a significant negative mediation effect of -11% (95%CI: -44%- -2.3%). No DNA methylation or clinical biomarkers demonstrated significant mediation effects in the comparison between other races and White participants due to small sample size.

Figure 2. Average mediated proportions for DNA methylation and clinical mediators on the relationship between race and ethnicity and all-cause mortality, NHANES 1999-2002 (n=2,402). Coefficients marked with one asterisk indicate statistical significance (p < 0.05), and those with three asterisks indicate significance after Bonferroni correction (p < 0.0021). The x-axis is ordered as chronological age clocks, physiological functioning clocks, pace-of-aging clock, clocks trained on other aging biomarkers, and behavioral, anthropometric, and clinical mediators.

The mediating roles of DNA methylation measures and clinical factors in the association between education and all-cause mortality

The mediation analysis for educational mortality disparities is shown in Figure 3 and Supplementary Table 4. The patterns of mediation were consistent across all comparisons of educational attainment levels. Physiologically trained clocks showed stronger effects compared to other types of clocks. GrimAge, GrimAge2, and DunedinPoAm epigenetic age exhibited significant mediation effects across all comparisons, ranging from 11% (95%CI: 2.8%-28%) to 32% (95%CI: 17%-64%), which were higher than the smoking pack-years, the leading clinical biomarkers. Horvath predicted telomere mediated 5.3% (95%CI: 1.7%-10%) and 5.7% (95%CI: 1.9%-13%) of the mortality disparities for participants with less than a high school education and those with a high school degree or educational development program (GED), respectively. In addition, PhenoAge epigenetic age mediated 6.5% (95%CI: 2.4%-14%) of the difference between participants with less than a high school education and those with a college degree, and Hannum epigenetic age showed a marginally significant positive mediation of 4.4% (95%CI: 0.3%-13.4%) between participants with a high school degree or GED and those with a college degree or higher.

Figure 3. Average mediated proportions for DNA methylation and clinical mediators on the relationship between education attainment and all-cause mortality, NHANES 1999-2002 (n=2,402). Coefficients marked with one asterisk indicate statistical significance (p < 0.05), and those with three asterisks indicate significance after Bonferroni correction (p < 0.0021). The x-axis is ordered as chronological age clocks, physiological functioning clocks, pace-of-aging clock, clocks trained on other aging biomarkers, and behavioral, anthropometric, and clinical mediators.

The mediating roles of DNA methylation measures and clinical factors in the association between household income and all-cause mortality

The mediation analysis for income-related disparities in all-cause mortality is presented in Figure 4 and Supplementary Table 5. The results followed patterns similar to those observed for educational attainment. PhenoAge, GrimAge, GrimAge2, Horvath predicted telomere, and DunedinPoAm showed significant positive proportion of mediation effects with smaller effects compared to the results from education attainment, but still greater effects compared to all the clinical mediators. GrimAge, GrimAge2, and DunedinPoAm epigenetic age significantly mediated effects across all group comparisons ranging from 14% (95%CI: 6.6%-29%) to 37% (95%CI: 22%-58%). In addition, PhenoAge and Horvath predicted telomere significantly mediated part of the group comparisons, ranging from 3.7% (95%CI: 0.3%-11%) to 8.3% (95%CI: 2.5%-17%).

Figure 4. Average mediated proportions for DNA methylation and clinical mediators on the relationship between household income and all-cause mortality, NHANES 1999-2002 (n=2,402). Coefficients marked with one asterisk indicate statistical significance (p < 0.05), and those with three asterisks indicate significance after Bonferroni correction (p < 0.0021). The x-axis is ordered as chronological age clocks, physiological functioning clocks, pace-of-aging clock, clocks trained on other aging biomarkers, and behavioral, anthropometric, and clinical mediators.

The mediating roles of DNA methylation measures and clinical factors in the association between occupation and all-cause mortality

The mediation of DNA methylation clocks on occupation-related disparities is shown in Figure 5 and Supplementary Table 6. Compared to the white-collar and professional participants, individuals in semi-routine white-collar and those who had never worked did not show significant mediations of any DNA methylation clocks and clinical biomarkers. For the comparison between high-skilled and semi-routine blue-collar and white-collar/professional participants, physiologically trained clocks again showed stronger effects compared to other types of clocks and clinical mediators. GrimAge, GrimAge2, and DunedinPoAm epigenetic age mediated 25% (95%CI: 14%-67%) up to 52% (95%CI: 26%-128%). In addition, Hannum, PhenoAge, and Horvath predicted telomere mediated 10% (95%CI: 2.2%-27%), 9.0% (95%CI: 2.1%-35%), and 6.1% (95%CI: 0.9%-16%) of the mortality disparities, respectively, when comparing high-skilled blue-collar workers to white-collar professionals.

Figure 5. Average mediated proportions for DNA methylation and clinical mediators on the relationship between occupation status and all-cause mortality, NHANES 1999-2002 (n=2,402). Coefficients marked with one asterisk indicate statistical significance (p < 0.05), and those with three asterisks indicate significance after Bonferroni correction (p < 0.0021). The x-axis is ordered as chronological age clocks, physiological functioning clocks, pace-of-aging clock, clocks trained on other aging biomarkers, and behavioral, anthropometric, and clinical mediators.

Secondary analysis of cause-specific mortality disparities

The exposure survival analysis for cause-specific mortality is also presented in Supplementary Table 1. For cardiovascular-specific mortality, due to the smaller simple size, only the comparisons between Black and White participants and between participants with a poverty income ratio (PIR) of less than 1 and those with a PIR of 5 or above remained significant with a HR of 2.12 (95%CI: 1.25-3.58) and 4.33 (95%CI: 1.63-11.51), respectively. For cancer-specific mortality, only models for educational attainment and occupational status showed significant results. Compared to participants with a college degree or higher, those with less than a high school education (HR=3.53, 95%CI: 1.43-8.70), a high school diploma or GED (HR=3.21, 95%CI: 1.41-7.29), and some college (HR=3.34, 95%CI: 1.52-7.33) had significantly higher risks. Similarly, participants in high-skilled blue-collar (HR=2.63, 95%CI: 1.50-4.63) and semi-routine blue-collar (HR=1.95, 95%CI: 1.13-3.34) occupations had significantly higher cancer-specific mortality risks compared to white-collar professionals. The mediator survival analysis for both cardiovascular- and cancer- specific mortality are shown in Supplementary Table 2. GrimAge, GrimAge2, and DunedinPoAm epigenetic age exhibited significant effects similar to the all-cause mortality results. For cancer-specific mortality only, more clocks showed significant associations, with HRs ranging from 1.04 (95%CI: 1.01-1.07) to 1.05 (95%CI: 1.02-1.07), including Hannum, Horvath, Vidal-Bralo, and PhenoAge.

The mediation analysis for cause-specific mortality had lower statistical power than for all-cause mortality (Supplementary Tables 7-14). For racial disparities, comparing Black and White participants, GrimAge, GrimAge2, and DunedinPoAm epigenetic age significantly mediated 17% (95%CI: 6.6%-56%), 23% (95%CI: 8.7%-63%), and 10% (95%CI: 2.0%-37%) of the cardiovascular-specific mortality disparity, respectively, though these effects were smaller than that of sedentary behavior (mediated proportion=26%; 95%CI: 11%-92%) (Supplementary Table 7). However, for cardiovascular disparities between those with a PIR <1 and those with a PIR ≥5, GrimAge and GrimAge2 epigenetic age mediated 23% (95%CI: 11%-50%) and 26% (95%CI: 11%-62%) which were higher than other clinical mediators (Supplementary Table 9). In cancer-specific mortality disparities, GrimAge, GrimAge2, and DunedinPoAm epigenetic age also demonstrated significant mediation for education- and occupation-related disparities, with mediated proportion ranging from 6.5% (95%CI: -3.4%-21%) to 30% (95%CI: 12%-66%). These three DNA methylation epigenetic age measures mediated even larger proportions for income-related cancer-specific mortality disparities with 76% (95%CI: 30%-326%), 80% (95%CI: 28%-416%), and 54% (95%CI: 12%-241%), respectively (Supplementary Table 13).

Sensitivity analysis of all-cause mortality stratified by subgroups

The mediation analysis stratified by gender is shown in Supplementary Tables 15–22. Although no significant difference was found, it is interesting to notice that the mediation signals for all-cause mortality disparities between Hispanic and White participants and people from different levels of household incomes were primarily driven by women, in comparison to all-cause mortality disparities between Black and White participants, people from different levels of education attainments, and individuals from different levels of occupation status were driven by men. The mediation analysis stratified by race and ethnicity is presented in Supplementary Tables 23–31. Again, although no significant difference shown across different races and ethnicities, DNA methylation biomarkers showed the most signals among White participants and then followed by Black participants. No significant signals were observed for Hispanic participants.

For more sensitivity analysis, we provided the results for all ten components of GrimAge2 in Supplementary Tables 32–35. These results showed that overall mediation effects of GrimAge2 were larger than that of any individual component, with the pack-years component showing the strongest effect among them. We also provided sensitivity analysis using completed cases only in Supplementary Tables 36–39. The results showed similar patterns compared to our main results using the imputed dataset. In addition, we provided testing results for exposure-mediator interaction and non-linear square term effects of each mediator on all-cause mortality in Supplementary Tables 40, 41. As the results demonstrated there was no consistent evidence across groups showing significant interactions between exposures and mediators in Supplementary Table 40. Supplementary Table 41 also suggested there was no evidence showing any non-linear effects of DNA methylation clocks on the all-cause mortality outcomes, although among clinical mediators, BMI, total cholesterol, and C-reactive protein presented significant square terms with betas that ranged from -0.037 to 0.066.

The NHANES participants

We analyzed the participants from the U.S. National Health and Nutrition Examination Survey (NHANES) 1999-2000 and 2001-2002 cycles who had publicly available DNA methylation clock measurements. Eligible participants were non-Hispanic White, non-Hispanic Black, Hispanic, and other races (including multi-racial participants), aged 50 years or older who had a blood sample collected at baseline and consented to future genetic research. Further details on the methodology of DNA methylation data are available on the NHANES website [45]. Participants were included in this study if they had non-zero DNA methylation specific weights (N=2,532). We excluded individuals who were 85 years and older (N=130) since NHANES does not report the exact ages for this group. After exclusion, the final study sample consisted of 2,402 participants.

Socioeconomic position measures

Socioeconomic position measures were derived from NHANES baseline questionnaires, with all variables being self-reported and publicly available. Educational attainment was assessed as the highest self-reported level of education completed by participants. This variable was classified as less than high school, high school graduate or general educational development program (GED), some college or an associate’s degree, and a college degree or higher.

Household income was measured using the poverty income ratio (PIR), defined as the ratio of a participant’s household income to the federal poverty levels from Department of Health and Human Services’ (HHS) poverty guidelines. To facilitate comparisons, household income was categorized into four groups: below the poverty line (PIR < 1), PIR 1-2, PIR 2-5, and PIR ≥ 5 for comparison with other work [46, 47].

The longest held occupation of the participants was used to define their occupational status. When data on the longest held occupation was missing, the reported current occupation was used to fill the missing values. Occupational status was categorized into five groups: white-collar and professional, semi-routine white-collar, high-skilled blue-collar, semi-routine blue-collar, and never worked, according to the National Statistics Socio-Economic Classification [48].

Epigenetic aging and clinical mediators

Thirteen DNA methylation clocks available in NHANES were included in this analysis. Clocks trained on chronological age included Hannum [49], Horvath [50], Weidner [51], Lin [52], VidalBralo [53], SkinBlood [54], and Zhang [55]. Clocks trained to predict morbidity or mortality included PhenoAge [56], GrimAge [57], and GrimAge2 [58]. In addition, the pace of aging clock DunedinPoAm [59], which was trained on aging rates measured across multiple organ systems, was included. Other clocks trained on aging biomarkers such as Yang [60] which was trained on cancer biomarkers, and Telomere [61] which was trained on leukocyte telomere length were also included.

To assess the mediation effects for traditional clinical measures from established biological pathways, additional clinical measures were included in tests of proportion mediated. For behavioral variables, we included physical activity (sedentary vs. active), the 2015 health eating index (HEI) [62], smoking pack-years, and drinking behaviors (abstainer vs. alcohol drinker). Anthropometric measures included waist-to-thigh ratio and body mass index (BMI). Clinical biomarkers included total cholesterol, high-density lipoprotein (HDL), low-density lipoprotein (LDL), glucose, and C-reactive protein (CRP). Physical activity status was defined as sedentary if participants did not report any vigorous, moderate, or muscle-strengthening activities over the past 30 days [63].

Mortality outcomes and covariates

Mortality data were obtained from the National Center for Health Statistics publicly available mortality files, which track deaths through linkage with the National Death Index (NDI). Participants were classified as deceased if their records were linked to the NDI by December 31, 2019, and were assumed to be alive otherwise. The underlying cause of death as recorded in the NDI has been validated with an estimated classification error of approximately 5% [64]. The maximum follow-up period was 17 years. The age at death or loss of follow-up was used as the survival outcome. All-cause mortality was examined as the primary event with loss to follow-up considered as censored, while cause-specific mortality from either cardiovascular disease or cancer were evaluated as the secondary and exploratory event with all other deaths and loss to follow-up both considered as censored. As shown in the directed acyclic graph (DAG) presented in Supplementary Figure 1, covariates controlled in the model included gender (women vs. men), baseline age and baseline age² (to account for potential non-linear effects of age on predictor and outcome variables), and nativity (U.S.-born vs. foreign-born).

Statistical modeling

Descriptive analyses were performed using means and standard deviations for continuous variables, and frequencies with percentages for categorical variables. A nonparametric random forest multiple imputation method was used to impute the missingness in exposures, covariates, and clinical mediators (missForest package in R [65, 66]). All the participants from 1999-2002 cycle (N = 21,004) were used to impute missing values. This package is designed to handle data that are Missing at Random (MAR). Correlation analyses to show the clusters across all DNA methylation biomarkers were conducted using R functions rcorr. Survival analyses were conducted using a Weibull model implemented in the survey package in R, incorporating NHANES DNA methylation specific survey weights to account for the complex survey design. Hazard ratios (HR) and 95% confidence intervals (CI) were reported. Only linear terms of mediators were included in the models, and Wald-test was used to test the non-linear square term effects of mediators in the associations between mediator and outcome.

Mediation analyses were conducted using the mediation package in R [67], estimating average mediated proportions along with their 95% CIs. This package outputs mediation proportions ranging from negative infinity to positive infinity. A negative mediated proportion indicates that the exposure and mediator have effects on the outcome in opposite directions, whereas a mediated proportion greater than 100% indicates that the mediator captures an effect in the same direction that is larger than the total effect of the exposure on the outcome. Exposure may influence an outcome through multiple pathways, producing both beneficial and detrimental effects that operate in different directions. As a result, the overall effect estimate may be attenuated relative to the true causal effect in that direction. This phenomenon is particularly relevant for social exposures, as some adverse social conditions may also be associated with social cohesion and resilience [5, 14]. A more detailed discussion is included in Appendix 1.

This mediation model uses a potential outcome framework to estimate mediation effects. We used a parametric approach within this package, which means the estimates will be calculated from the regression parameters and the standard errors will be calculated from bootstrap method. The assumptions of this approach include 1) no unobserved pre-treatment confounding between social positions and mortality outcomes besides the covariates already controlled; 2) no unobserved pre-treatment and post-treatment confounding for the associations of DNA methylation clocks and mortality outcomes; 3) sufficient overlap in the distributions of social position and DNA methylation clocks across covariates. A more detailed discussion of how well our study met the assumption is included in Appendix 2. The package built-in function test.TMint was used to test if there were interactions between exposures and mediators. Results with a p-value of <0.05 were considered as marginally significant, and a p-value of <0.0021 (0.05/24, as there was total 24 mediators included in the study) were considered as significant results that passed multiple testing. Sensitivity analyses conducted among different gender and racial subgroups were also presented. All statistical analyses were performed using R version 4.5.0. The R code and data used for these analyses are provided in a GitHub repository (https://github.com/shywatson/SES_DNAm_NHANES).