Variations in Innate Immune Cell Subtypes Correlate with Epigenetic Clocks, Inflammaging and Health Outcomes

Abstract

Epigenetic clocks in blood have shown promise as tools to quantify biological age, displaying robust associations with morbidity and all-cause mortality. Whilst the effect of cell-type heterogeneity on epigenetic clock estimates has been explored, such studies have been limited to studying heterogeneity within the adaptive immune system. Much less is known about whether heterogeneity within the innate immune system can impact epigenetic clock estimates and their associations with health outcomes. Here, we apply a high-resolution DNAm reference panel of 19 immune cell-types, including young and adult monocyte, natural killer, and neutrophil subsets, demonstrating how shifts within these innate subtypes display associations with epigenetic clock acceleration, inflammaging, and all-cause mortality. The associations of monocyte heterogeneity with inflammation are further validated using transcriptomic and metabolomic data. Additionally, a non-negligible fraction of nucleated red blood cell-like cells in circulation is found to associate with inflammaging, markers of dysfunctional erythropoiesis, and is a major risk factor for all-cause mortality. These results extend findings obtained within the adaptive immune system to innate immune and erythrocyte-like cells, demonstrating how heterogeneity within these other blood cell compartments is also associated with inflammaging, epigenetic clocks, and health outcomes.

Keywords: DNA methylation; aging Biomarkers; epigenetic clocks; health outcomes; inflammaging; innate immune system.

Figure 1

A common axis of DNAm variation reflects maturity and aging across adaptive and innate immune cell‐types. a) Sorting and DNAm profiling of 7 immune cell‐types found in cord‐blood (birth) and 12 immune cell‐types found in adult blood (prefixed with “a”) could be used to identify shifts of “young” and “old” cellular phenotypes during adulthood, within both adaptive and innate immune cell‐types. b) Conceptualization of our hypothesis: individuals aging faster than normal (according to some epigenetic clock) could be driven in part by a subtle shift in the balance between “young” and “old” immune cell subtypes. In the case of an adaptive immune cell‐type, it is well‐known that the adult naïve CD8+ T‐cell fraction decreases with age at the expense of an increase in the adult memory CD8+ T‐cell fraction, and that further tilting of this balance at any age indicates DNAm age acceleration or deceleration, as shown. In the case of an innate immune cell‐type, e.g. monocytes, we hypothesize that the balance between “young” and “old” subtypes could also associate with DNAm age acceleration, for instance, because the younger and older monocyte subtypes reflect non‐inflammatory, respectively, inflammatory subtypes. c) PCA plot of a merged Illumina 450k DNAm cohort of over 220 sorted samples, encompassing the 19 immune cell‐types depicted in a). It is observed how PC2 describes an axis of DNAm variation that correlates simultaneously with maturity and age, and that is common to all adaptive and innate immune cell‐types. PC1 describes an axis of DNAm variation associated with lymphoid versus myeloid status. d) Using the merged DNAm cohort, a 19 immune cell‐type UniLIFE DNAm reference panel was built. Details of its construction and validation are elaborated in Guo et al. The UniLIFE reference panel was applied to a collection of 15 adult whole blood cohorts to estimate corresponding fractions for the 19 immune cell‐types, as shown in the boxplots. Parts of panel a)+b) were generated with Biorender.com.

Figure 2

Shifts within innate immune cell‐types correlate with age‐acceleration. a) For each of the 19 immune cell‐types in the UniLIFE DNAm reference panel, we give the effect sizes of association of the cell‐type fraction with the relative age‐acceleration of Horvath, PhenoAge, and GrimAge2 clocks, as well as with DunedinPACE's pace of aging, for a random effect meta‐analysis model, as assessed over 15 adult whole blood cohorts. Error bars represent 95% confidence intervals. b) For GrimAge2, detailed forest plots for 4 of the cell‐type fractions across the 15 adult blood cohorts: adult naïve CD4+ and CD8+ T‐cells (aCD4Tnv, aCD8Tnv), young (cord‐blood) monocytes (Mono), and adult monocytes (aMono). P‐values for the random effect model are given. Heterogeneity I² indices and P‐values are given above each panel. c) As b) but for DunedinPACE. *Significant under Bonferroni correction.

Figure 3

Elevated adult to young monocyte fractions correlate with inflammaging. a) Flowchart depicting the calculation of the DNAm‐based inflammaging score (Infl.Score) for a sample s. A large EWAS meta‐analysis of measured C‐reactive protein (CRP) levels in serum identified a total of 1765 CRP‐associated CpGs. For an independent DNAm dataset, we first z‐score the CpG DNAm profiles across all samples in the dataset, so that the mean is zero and the standard deviation of each CpG is 1. The InflScore of sample s is then defined as the Pearson correlation of z‐score normalized values with the bidirectional vector of the 1765 CpGs, with +1/‐1 indicating positive/negative correlation with CRP. b) Validation forest plot of the InflScore against chronological age, as assessed over 15 adult whole blood cohorts. P‐value of the random effect model is given. Above the panel, we list the I² heterogeneity index and associated P‐value. c) Boxplots of the InflScore in an EWAS of inflammatory bowel disease, comparing ulcerative colitis and Crohn's Disease cases to healthy controls. P‐values are from a one‐tailed Wilcoxon rank sum test. d) Forest plots of association of young (NK) and adult (aNK) natural killer fractions, nRBC, and adult memory B‐cell fractions with the DNAm‐based CRP Infl.Score. P‐values as in b). e) As d), but for young (Mono) and adult (aMono) monocyte fractions. f) As e), but the association with the DNAm‐based IL6 score. g) As e+f, but for estimated classical (C‐Mono) and non‐classical (NC‐Mono) fractions. *Significant under Bonferroni correction.

Figure 4

Age‐related monocyte subtyping correlates with inflammation in transcriptomic and metabolomic data. a) Venn diagram between genes differentially expressed between cord and adult blood classical monocytes, with genes differentially expressed in classical monocytes between old healthy and old frail individuals. Odds Ratio and P‐value of overlap from a one‐tailed Fisher test are given. b) Gene Set Enrichment Analysis (GSEA) of the 299 overlapping genes from a), displaying the enrichment of inflammation‐related KEGG pathways. A total of 48 inflammation‐related enriched genes were found. c) Scatterplot of AUC‐0.5 (derived from Wilcoxon rank sum test) for these 48 inflammation‐related genes between the adult versus cord blood classical monocytes (x‐axis) and between the frail and healthy old individuals (y‐axis). A Pearson correlation coefficient and P‐value are given. d) As c) but with the y‐axis now labeling the age‐associated t‐statistic as derived from the classical monocytes (n>200000) of the large AIDA scRNA‐Seq dataset. These statistics were derived from a linear mixed model including age and sex as fixed effects and donorID as a random effect. e) Balloon plot of Spearman Rank Correlation Coefficients (SCC) between inflammatory metabolites positively associated with IL‐6 levels with monocyte subtype fractions as estimated in the MGBB‐ABC cohort (n=1653). Cord and adult monocyte fractions were estimated using the 19 cell‐type UniLIFE panel, whereas the classical and non‐classical monocyte fractions were estimated using a 13 cell‐type panel.

Figure 5

Associations of immune cell fractions with disease incidence and mortality in Generation Scotland. a) Heatmap of Cox‐regression Hazard Ratios (HR) and associated P‐values, between selected medical conditions and 19 immune cell‐type fractions in the Generation Scotland (GenS) cohort. b) Further subselection of medical conditions and immune cell‐fractions, displaying the precise HR values and their 95% confidence intervals. c) Forest plots of Cox‐regression HRs, 95% confidence intervals, and P‐values for the 19 immune cell‐type fractions in association with all‐cause mortality in GenS, with the left panel adjusted only for chronological age, and the right panel adjusted for age, other demographic and lifestyle risk factors. The number of samples (nS) and death events (nE) are given below the x‐axis. d) Balloon correlation plot displaying the correlations of the nRBC fraction with clinical complete blood cell (CBC) count measurements in the MGBB cohort. Correlation coefficients are indicated for both Pearson (PCC) and Spearman (SCC). e) Forest plots depicting the association of Horvath, PhenoAge, and GrimAge2 clock extrinsic age acceleration (EAA) and DunedinPACE with all‐cause mortality in the GenS cohort, using a 10‐year follow‐up window. The x‐axis labels the Hazard Ratio (HR), and 95% confidence interval is given for each estimate. The left panel is after adjustment for all demographic and lifestyle risk factors, the middle panel is when also adjusting for the 19 immune‐cell fractions, right panel is when additionally adjusting for baseline co‐morbidity (cancer, CVD, T2D, COPD, depression). The number of samples (nS) and death events (nE) are given below the x‐axis. f) As e) but for cancer and cardiovascular disease (CVD) specific mortality, adjusting for demographics, lifestyle, and 19 IC‐fractions (left panels) and additionally adjusting for baseline cancer/CVD status (right panels). *In all panels, associations are significant after Bonferroni adjustment for multiple testing.