TIME-seq reduces time and cost of DNA methylation measurement for epigenetic clock construction

Abstract

Epigenetic 'clocks' based on DNA methylation have emerged as the most robust and widely used aging biomarkers, but conventional methods for applying them are expensive and laborious. Here we develop tagmentation-based indexing for methylation sequencing (TIME-seq), a highly multiplexed and scalable method for low-cost epigenetic clocks. Using TIME-seq, we applied multi-tissue and tissue-specific epigenetic clocks in over 1,800 mouse DNA samples from eight tissue and cell types. We show that TIME-seq clocks are accurate and robust, enriched for polycomb repressive complex 2-regulated loci, and benchmark favorably against conventional methods despite being up to 100-fold less expensive. Using dietary treatments and gene therapy, we find that TIME-seq clocks reflect diverse interventions in multiple tissues. Finally, we develop an economical human blood clock (R > 0.96, median error = 3.39 years) in 1,056 demographically representative individuals. These methods will enable more efficient epigenetic clock measurement in larger-scale human and animal studies.

Extended Data Fig. 1 |. Biotinylated-RNA bait production and initial hybridization enrichment testing.

a, Schematic of steps involved in production of biotinylated-RNA baits from single-stranded oligo pools for target enrichment in TIME-Seq libraries. The percent of reads overlapping target RRBS mouse rDNA clock CpGs (b) and an IGV browser screenshot of mapped-read pileups (c) using version 1 rDNA baits for enrichment of a TIME-Seq pool. Reads on-target (d) and mouse RRBS blood clock (Petkovich et al., 2017) CpG coverage (e) using mouse-blood specific baits in a pilot experiment targeting non-repetitive clock loci. Dotted line represents coverage cut-off of 10. Pools in both rDNA and blood clock pilot enrichments were sequenced with approximately 1 million paired end (PE) reads each in pool of 16 samples. (f) Adaptor design schematic for comparison of TIME-Seq adaptors with longer barcoded adaptors. Comparison of on-target reads in short TIME-Seq and long cytosine-depleted adaptor designs for both mouse blood clock (g) and (h) rDNA (version 1) baits enrichments.

Extended Data Fig. 2 |. TIME-Seq library and sequencing schematic.

Schematic representation of final library structure (top) and Illumina sequencing (bottom) steps required to sequence TIME-Seq libraries. Index read 1 and read 2 primers are custom primers.

Extended Data Fig. 3 |. Small pilot-experiment sample metrics, correlation of rDNA CpG methylation, and age predictions using a reported RRBS-based rDNA clock.

a, TIME-Seq pilot experimental design using mouse blood DNA from 4 age groups and preparing 2 replicates of each sample with rDNA baits (version 1) as well as RRBS libraries to be sequenced as a fraction of an Illumina MiSeq sequencing run. b, Demultiplexed reads from TIME-Seq pools. c, Mean CpG methylation from reads mapped to the mouse ribosomal DNA meta-locus. Unmethylated lambda phage DNA control is represented as a diamond. d, Percent methylation from reads mapped to ribosomal DNA meta-locus in replicate 1 and replicate 2 in CpGs with coverage of at least 125 reads. e, Replicate correlation from different coverage cutoffs in the rDNA. f, Pileup tracks for samples from a TIME-Seq pool (replicate 1) as well as mapped reads from one sample (mouse ID 3, aged 24 months). Reads are colored by mismatch: blue for T (unmethylated) and red for C (methylated). RRBS rDNA clock coordinates are illustrated on the bottom by black rectangles. g, Percent of reads directly overlapping clock CpGs from TIME-Seq libraries (N = 12; mean from 2 replicates) and shallow-sequenced RRBS libraries (N = 10). h, RRBS rDNA clock predictions using TIME-Seq data enriched for clock loci (N = 12, both replicates) i, Coverage of each clock locus in the original RRBS rDNA clock. CpGs shown in red have a mean coverage of less than 50. Boxplot lengths (panels b, c, g, h) represent the interquartile range (IQR) with the middle line representing median values and the whiskers 1.5 times the IQR.

Extended Data Fig. 4 |. Additional data related to mouse multi-tissue and tissue-specific clock training and testing.

a, Baits overlapping target loci used for mouse clock CpG enrichment. b, Age predictions from the TIME-Seq Mouse Multi-tissue Clock applied to the 157 mouse muscle samples. Pearson correlation between predicted and actual age is shown. c, TIME-Seq White Adipose Clock train (N = 107) and testing set (N = 27) predictions. d, TIME-Seq Kidney Clock train (N = 156) and testing set (N = 38) predictions. For panels c and d, Pearson correlation between predicted and actual age is shown for train and test. The median absolute error is shown for the testing set.

Extended Data Fig. 5 |. Additional Data from validation and benchmarking of TIME-Seq.

a-b, TIME-Seq Mouse rDNA Clock predictions with samples colored for (a) validation library preparation (prep) and (b) cohort of the mouse. Pearson correlation is shown for each panel. c, Correlation between age-adjusted prediction residuals in the validation sets from the different prediction approaches. d, Correlation and significance matrix between ∆Age from each approach and ∆Med^age(blood), that is, the difference in median value from similar aged mice for each blood measurement. The color and size of each circle represent the correlation and p-value significance, respectively. WBC = white blood cell count, NE (%) = percent of neutrophils, LY (%) = percent of lymphocytes, MO (%) = percent of monocytes, EO (%) = percent of eosinophils, BA (%) = percent of basophils, RBC = red blood cell count, Hb = hemoglobin, HCT = hematocrit, MCV = mean corpuscular volume, MCH = mean corpuscular hemoglobin, MCHC = mean corpuscular hemoglobin concentration, RDW = red blood cell distribution width, PLT = platelets, MPV = mean platelet volume. e, Frailty indexes for each of the assayed mice along with Pearson correlation with age. f, Comparison of ∆Age and ∆Med^age(FI) for mice in the validation cohort. Pearson correlation is shown without adjusting for multiple comparisons. g, Comparison between TIME-Seq CpG methylation and RRBS methylation in the same sample and CpG. Pearson correlation between CpG methylation levels is shown.

Extended Data Fig. 6 |. Additional Data for TIME-Seq clocks applied to intervention mice and an in vitro time course.

a, Comparison TIME-Seq Multi-Tissue Clock predictions of high-fat diet mouse liver (N = 12) with standard diet controls (N = 5). b, Comparison of TIME-Seq Liver Clock predictions in OSK-expressing, (+) OSK (N = 5), and control, (−) OSK (N = 9), mice. For panels a-b, statistical comparison between groups was performed using a two-sided Student’s t-test after assessing normality with Shapiro-Wilk’s test. c, Predictions of cell culture samples using the TIME-Seq Mouse Skin Clock. The slope from the linear models fit to data points from each cell line is shown.

Extended Data Fig. 7 |. Comparison of methylation levels from TIME-Seq and BeadChip on the same samples.

Pearson correlation between BeadChip and TIME-Seq DNAme values for each sample and CpG (R = 0.93, p < 2.2e − 16).

Fig. 1 |. TIME-seq enables highly efficient epigenetic age predictions.

a, Schematic of the TIME-seq library preparation for highly multiplexed targeted methylation sequencing to build and measure DNA methylation-based biomarkers. b, Proof-of-concept rDNA clock experiment schematic; 191 mouse blood DNA samples (histogram) were prepared with TIME-seq enriched for rDNA and sequenced for clock training and testing. c, Reads were demultiplexed from each rDNA clock pool. d, Percentage of demultiplexed reads from each sample that mapped to the rDNA meta-locus. e, Mean coverage at the rDNA meta-locus CpGs in rDNA-enriched TIME-seq libraries. f, Mean CpG methylation from each sample in the four pools. In c,d,f, n = 48 for pools 1–3 and n = 47 for pool 4. The boxplot lengths represent the interquartile range (IQR) with the middle line representing median values and the whiskers 1.5 times the IQR. g, Histogram of training (n = 145; red) and testing (n = 37; blue) samples used to develop the TIME-seq mouse rDNA clock. h, TIME-seq mouse rDNA clock age predictions for the training (red) and testing (blue) datasets. Pearson correlation and MedAE on the test are shown in the top left corner. i, Predictions from a TIME-seq rDNA clock developed using only CpGs with at least 50 reads coverage in the RRBS data used to develop the original mouse rDNA clock. Calorie-restricted (CR) mice are represented as red triangles.

Fig. 2 |. Low-cost TIME-seq multi-tissue and tissue-specific clocks applied to 1,137 mouse tissue samples.

a, Circular genome plot illustrating the position and mean coverage of the 6,370 high-coverage CpGs from the TIME-seq libraries in 1,137 mouse tissue samples. b, Principal component analysis (PCA) for the TIME-seq data colored according to their tissue of origin: liver (green), blood (red), skin (blue), muscle (turquoise), kidney (yellow) and WAT (pink). c, TIME-seq mouse multi-tissue clock training (left) and testing (right) predictions plotted against chronological age. d, Linear models fitted to DNA methylation levels at clock CpGs changing with age in the data used to train the TIME-seq mouse multi-tissue clock. Clock CpGs were split according to their clock coefficient sign, which is represented in the transparency of each line. e, Enrichment of protein binding (based on ENCODE data) at genes associated with the TIME-seq mouse multi-tissue clock. f–h, Tissue-specific TIME-seq clock training (left) and testing (right) datasets for the TIME-seq mouse blood clock (f), liver clock (g) and skin clock (h). For all clocks, the Pearson correlation between predicted and actual age is shown. MedAE is shown for the testing set predictions. i, Age-adjusted residuals for liver and skin clock predictions from the same mice predicted in the testing sets. The shading around the regression line represents the 95% confidence interval (CI). j, Age-adjusted residuals for either skin (blue) or liver (green) clock predictions plotted against the predictions from the multi-tissue clock in the same sample.

Fig. 3 |. TIME-seq is a robust and scalable alternative to conventional clock approaches.

a, Experimental schematic for validation of the TIME-seq age prediction methods in an independent cohort of mice with longitudinal time points, paired frailty index and blood composition data. b,c, TIME-seq age predictions in two independent validation library preparations using the TIME-seq mouse blood clock (n = 75) (b) and the TIME-seq mouse multi-tissue clock (n = 74) (c). The lines connect the same mouse at two different ages. Pearson correlations are shown in the top left corner. d,e, Sequencing saturation simulation to estimate clock accuracy from different read numbers in the validation samples for the TIME-seq mouse blood (d) and multi-tissue (e) clocks. The shading around the locally estimated scatterplot smoothing regression line represents the 95% CI. f, Schematic of the benchmarking experiment to compare TIME-seq to Illumina BeadChip and RRBS. g, Comparison of CpG methylation percentage in TIME-seq and BeadChip. Each dot represents the same CpG from the same sample measured by each technology. h, Comparison of ∆Ages for each method and the associated clocks. The boxplot lengths represent the IQR, with the middle line representing the median values and the whiskers 1.5 times the IQR. i, Comparison of cost per sample (top) and total cost (bottom) for the TIME-seq clocks (black), BeadChip (red) and RRBS (blue) across a range of sample scales. The half-filled circles denote points that are overlapping between TIME-seq clocks. m, million; k, thousand.

Fig. 4 |. TIME-seq clocks reflect interventions that slow, accelerate and reverse aging and can be used for in vitro studies.

a, Schematic of dietary restriction treatments started in late life. Groups of mice were treated with a 40% CR diet, MetR diet and ad libitum diet for 6 months starting at 24 months of age. b, Comparison of ∆Ages from TIME-seq blood clock predictions from blood samples from MetR (n = 19), CR (n = 8) and ad libitum (n = 13) fed groups of mice. Group comparison was performed with a pairwise two-sided Wilcoxon test with false discovery rate (FDR) correction for multiple comparisons. c, Schematic of the experiment to assess the effect of a CR diet on mouse liver epigenetic age. Mice were treated with 30% CR or control ad libitum diets starting at 4 months of age for 10 months; then livers were collected. d, Comparison of ∆Age from mice fed ad libitum (n = 4) and CR (n = 5) mice using the TIME-seq mouse liver clock. e, Schematic of HFD treatment and liver collection. Mice were treated for 13 months starting at 4 months. f, Comparison of predicted ages for standard (n = 5) and HFD mice (n = 12) using the TIME-seq liver clock. f, Statistical comparison between groups was performed using a two-sided Student’s t-test after assessing normality with a Shapiro–Wilk test. g, Schematic of AAV treatments with OSK-expressing or control (GFP) cassettes. Livers were collected 1 month after AAV injection from two sets of mice (first set: control, n = 4 treatment, n = 2 aged 15–17 months; second set: control, n = 9; treatment, n = 5 aged 24 months). h, Comparison of ∆Ages in the livers of mice with (OSK⁺) or without (OSK⁻) OSK expression for 1 month. Statistical comparison was performed with a two-sided Wilcoxon test. i, Schematic of the experiment to assess MEFs and mouse adult ear fibroblasts in a cell culture time course lasting 1 month, with collection every 2 weeks. j, TIME-seq multi-tissue clock predictions from cell culture samples collected across the time course. The slope represents the change in predicted age per day in culture based on linear models fitted to data from each cell type.

Fig. 5 |. Highly accurate epigenetic age predictions in 1,056 human blood samples using TIME-seq.

a, Schematic of the experimental design to train, test and validate TIME-seq in 1,056 human blood DNA samples. Wh+H indicates White (Hispanic) ethnicity as denoted in the biobank metadata. b, Coverage of 9,379 CpGs from across the human genome. The colored dots represent CpGs described in the Illumina BeadChip clocks, whereas the smaller gray dots are the other enriched CpGs. c, PCA of the methylation matrix for the training and testing samples, colored according to age from youngest (blue) to oldest (red). d, Predicted ages from the TIME-seq human blood epigenetic clock. Pearson correlation between predicted and actual age is shown for the training (right; R = 0.98, P < 2.2 × 10⁻¹⁶) and testing (left; R = 0.96, P < 2.2 × 10⁻¹⁶) datasets. The MedAE is shown for the testing dataset. e, Annotation of the 405 clock CpGs with the coefficient on the y axis. The x axis (not shown) is the genomic space in the same style as b from left (chromosome 1) to right (chromosome 22). CpGs and gene names are colored with the same color key as in b. Feature annotation is coded according to shape as follows: 5′ UTR (open circle), exon (filled diamond), intergenic (open triangle), intron (filled triangle), noncoding (inverted filled triangle), promoter (filled square) and transcription termination site (inverted open triangle). f,g, Gene Ontology (GO) analysis for the enrichment of biological processes (f) or transcription factor binding sites (g) in genes associated with clock CpGs. h, TIME-seq human blood clock predictions in 260 independently prepared human blood DNA samples. Pearson correlation between predicted and actual age is shown (R = 0.96, P < 2.2 × 10⁻¹⁶). i, Predicted ages from each BeadChip-based clock in a subset of the validation set sample (n = 24). Pearson correlation between predicted and chronological age is shown, as well as MedAE in the subset. For Horvath1, Horvath2 and Hannum, P < 2.2 × 10⁻¹⁶. For PhenoAge, P = 4.7 × 10⁻¹⁵. j, ∆Age values (difference in predicted and chronological age in units of years) for each age prediction method. k, Comparison of Pearson correlations between age-adjusted prediction residuals for each TIME-seq-based prediction and BeadChip-based clocks.