Epigenetic editing at individual age-associated CpGs affects the genome-wide epigenetic aging landscape

Abstract

Aging is reflected by genome-wide DNA methylation changes, which form the basis of epigenetic clocks, but it is largely unclear how these epigenetic modifications are regulated and whether they directly affect the aging process. In this study, we performed epigenetic editing at age-associated CpG sites to explore the consequences of interfering with epigenetic clocks. CRISPR-guided editing targeted at individual age-related CpGs evoked genome-wide bystander effects, which were highly reproducible and enriched at other age-associated regions. 4C-sequencing at age-associated sites revealed increased interactions with bystander modifications and other age-related CpGs. Subsequently, we multiplexed epigenetic editing in human T cells and mesenchymal stromal cells at five genomic regions that become either hypermethylated or hypomethylated upon aging. While targeted methylation seemed more stable at age-hypermethylated sites, both approaches induced bystander modifications at CpGs with the highest correlations with chronological age. Notably, these effects were simultaneously observed at CpGs that gain and lose methylation with age. Our results demonstrate that epigenetic editing can extensively modulate the epigenetic aging network and interfere with epigenetic clocks.

Fig. 1. Epigenetic editing is stable but not coherent within the target region.

a, Schematic presentation of two CRISPR-guided epigenetic editors used in this study. Deficient CAS9 protein is linked to DNMT3A/3L at the C terminus (dCAS9-DNMT3A/3L, coexpressed with eGFP for selection), and CRISPRoff comprising DNMT3A/3L at the N terminus, and tagBFP and a KRAB domain at the C terminus. Single guide RNAs were designed 50–200 bp distant of the target CpG. b, Manhattan plot illustrating the DNAm changes upon targeting the age-associated hypermethylated region in PDE4C (EPIC Illumina BeadChip data of chromosome 19; mean of three replicas; 14 days after transfection). Significant DNAm changes are highlighted in red: delta mean DNAm > 0.1 and an adj. P < 0.05 (limma P value, Benjamini–Hochberg adjusted). The target CpG in PDE4C is highlighted in black. c, Bisulfite amplicon sequencing of 26 neighboring CpGs at the target region of PDE4C. The bar plot depicts the frequency of methylated CpGs on individual reads, indicating that even on modified DNA strands not all neighboring CpGs become coherently methylated. 61.7% and 42.9% of reads have higher methylation levels than observed in the wild type, for dCAS9-DNMT3A and CRISPRoff, respectively. d, Time-course experiment of DNAm at PDE4C measured by bisulfite barcoded amplicon sequencing. The lines resemble the 26 different CpGs in the amplicons. e, Frequency of DNAm patterns in bisulfite barcoded amplicon sequencing data of PDE4C. Reads are clustered by their DNAm pattern. The binding region of one gRNA is indicated. The second gRNA binds two base pairs next to the CpG #26 and might explain the low methylation gain at #25 and #26.

Fig. 2. Genome-wide bystander effects of epigenome editing.

a,b, Significant DNAm changes upon epigenetic editing at PDE4C with either dCas9-DNMT3A (a) or CRISPRoff (b). Volcano plots depict significantly hypermethylated (red) and hypomethylated (blue) CpGs and their numbers are indicated (n = 3; 14 days after transfection; limma P value, Benjamini–Hochberg adjusted). c, Multivariate density estimate comparing DNAm changes in the dCAS9-DNMT3A and CRISPRoff experiments. Pearson correlation (R² = 0.53) indicates that there is a high reproducibility of the bystander modifications, even with the different epigenetic modification approaches. d, Enrichment of the bystander effects in relation to CpG islands (CGI), and at shore and shelf regions surrounding CGIs. Enrichment was calculated in relation to all CpGs on the array and was highly significant for all categories: chi-squared test, degrees of freedom (d.f.) = 1, not adjusted for multiple testing, North Shelf P = 9.51 × 10⁻¹⁰ (CRISPRoff) and P = 2.80 × 10⁻¹² (DNMT3A), South Shelf P = 4.76 × 10⁻⁹ (CRISPRoff) and P = 8.79 × 10⁻¹³ (DNMT3A), all other P < 10⁻¹⁵. e,f, Relative frequency of nucleotides next to CpGs with epigenetic bystander modifications normalized to CpGs with identical CpG density (EPIC manifest). Guanine and cytosine were overrepresented at the −1 and +1 flank positions, whereas thymine and adenine were enriched at the −2 and +2 positions (separate chi-squared tests of the four proximate genomic positions, d.f. = 1, not adjusted for multiple testing, all four P < 10⁻¹⁵). g,h, Relative frequency of nucleotides next to CpGs that become either hypomethylated with age (g; 4,389 CpGs), or hypermethylated with age (h; 5,328 CpGs) in a large-scale epigenome-wide association study. Frequency was normalized to the entirety of CpGs on the BeadChip. Adenine and thymidine were overrepresented in the −1 and +1 flank positions of hypomethylated CpGs (separate chi-squared tests of the two proximate genomic positions, d.f. = 1, not adjusted for multiple testing, both P < 10⁻¹⁵). i,j, Distribution of epigenetic bystander modifications in the dCAS9-DNMT3A experiments (i) and the CRISPRoff experiments (j) was analyzed for all CpGs on the array, for the 4,389 CpGs with age-associated hypomethylation and 5,328 CpGs with age-associated hypermethylation. Age-associated hypermethylation was enriched at CpGs that gain DNAm upon epigenetic editing (both P < 10⁻¹⁵, two-sided, two-sample Kolmogorov–Smirnov test). k,l, Cumulative distribution of the density functions in i,j to better visualize that bystander effects upon targeting PDE4C are enriched at other age-hypermethylated CpGs.

Fig. 3. Chromatin conformation contributes to epigenetic bystander modifications.

a, ATAC-seq data of HEK293T cells was used to estimate chromatin accessibility at epigenetic bystander modifications upon modification of PDE4C with either dCas9-DNMT3A or CRISPRoff. The higher ATAC-seq signals indicate that bystander modifications are enriched at open chromatin regions (two-tailed, unpaired t-test). b, In analogy, we tested whether CpGs with age-associated DNAm changes are also associated to chromatin accessibility in the ATAC-seq data HEK293T cells. The 4,389 age-hypomethylated CpGs revealed overall lower ATAC signal, whereas the 5,328 age-hypermethylated CpGs were enriched at open chromatin (two-tailed unpaired t-test). c, Intrinsic 4C-sequencing in HEK293T cells for three age-associated sites in PDE4C, FHL2 and MEIS1-AS3 (triplicates per viewpoint). The Manhattan plots depict the number of reads (rpm, reads per million) at the cis-interacting regions on the corresponding chromosomes. The highest signal was observed around the 4C viewpoints. d, Bystander effects of PDE4C become enriched with increasing coverage in 4C-sequencing of PDE4C-interacting sites. The association is best described by an exponential model (adjusted R²: 0.586, model P value < 1.801 × 10⁻⁷). Different symbols correspond to the three replicates. e–g, Overlap of 4,389 age-hypo and 5,328 age-hypermethylated CpGs at interacting genomic regions in 4C-sequencing data (chi-squared tests were performed for the threshold of 100 rpm: PDE4C: age-hypermethylated sites P = 0.0004 and age-hypomethylated sites P = 0.0019 (e); FHL2: age-hypermethylated sites P = 0.0008 and age-hypomethylated sites P = 0.0088 (f); and MEIS1-AS3: age-hypermethylated sites P = 0.0008 and age-hypomethylated sites P = 0.009 (g); all d.f = 1).

Fig. 4. Multiplexed epigenetic editing at age-hyper- and age-hypomethylated genomic regions.

a, Scheme of multiplexed epigenetic editing at five genomic regions that gain DNAm with aging. b, Scatter-plot of DNAm changes across the three replicas. CpGs corresponding to the targeted genes are highlighted in the corresponding color. None of the bystander modifications reached statistical significance. c, Cumulative distribution function comparing the entirety of CpGs from the BeadChip with age-associated CpGs (4,389 age-hypo and 5,328 age-hypermethylated CpGs). d, Correlation of differential gene expression with differential methylation. CpGs and transcripts were matched by gene IDs. Pairs related to one of the five target genes are highlighted in the corresponding color. e, Scheme of multiplexed epigenetic editing at five genomic regions that lose DNAm with aging. f,g, Three days after transfection, DNAm changes were analyzed by a volcano plot (f) and a scatter-plot (g). The volcano plot shows DNAm changes and the limma P value (Benjamini–Hochberg adjusted). CpGs corresponding to the five target genes are highlighted by the corresponding colors (n = 3). h, At 15 days after transfection, the gains in DNAm in the target regions were hardly observed anymore. i, At day 3 after transfection, Gaussian kernel density estimate at age-associated CpGs (4,389 age-hypo and 5,328 age-hypermethylated CpGs) showed that bystander modifications were enriched at genomic regions that gain methylation with age. j, Cumulative distribution of DNAm at age-associated CpGs at day 3 after transfection. Notably, bystander effects are underrepresented at age-hypo and overrepresented at age-hypermethylated CpGs (both P < 10⁻¹⁵, two-sided, two-sample Kolmogorov–Smirnov test). k,l, At day 15 after transfection kernel density estimate and cumulative distribution did not reveal bystander effects anymore.

Fig. 5. CRISPR-epigenetic editing in human primary T cells.

a, Scheme of multiplexed epigenetic editing at five genomic regions that gain DNAm with aging. b, Scatter-plot of Illumina BeadChip data showing clear mean DNAm changes at FHL2 and KLF14 (n = 2). c,d, Gaussian kernel density estimates and cumulative distribution function showing significantly different methylation profiles at hyper- and hypo-CpGs (two-sided, two-sample Kolmogorov–Smirnov test P < 10⁻¹⁵). e, Scheme of multiplexed epigenetic editing at five genomic regions that lose DNAm with aging. f, Scatter-plot of Illumina BeadChip data revealing no mean methylation change after 10 days (n = 3). g,h, Gaussian kernel density estimates and cumulative distribution function. The genome-wide epigenetic landscape is less affected when targeting age-hypomethylated CpGs (10 days of culture). i,j, Estimation of epigenetic age with eight different epigenetic clocks:: Zhang, Vidal-Bralo, Lin, Horvath 1 (multi-tissue), Horvath 2 (skin and blood), Hannum and PhenoAge and updated PhenoAge. Bar plots depict the deviation of epigenetic age predictions upon targeted modification at the age-hypermethylated (i) and the age-hypomethylated CpGs (j), as compared to the scramble guide RNA controls (difference in years). k–n, Density plots of DNAm profiles for different strata of age-relatedness. Pearson correlation of individual CpGs with chronological age was determined in four T cell datasets^–. Upon targeting the hyper-CpGs or hypo-CpGs there is a clear enrichment of bystander modifications at sites with increasing positive (k,m), and negative (l,n) correlation with age (all P < 10⁻¹⁵, Kruskal–Wallis test).

Fig. 6. Epigenetic editing at age-hypermethylated CpGs in mesenchymal stromal cells.

a, Scheme of multiplexed epigenetic editing in MSCs at five genomic regions that gain DNAm with aging. b, Scatter-plot of DNAm changes across the three biological replicas. CpGs associated to targeted genes (10-kbp region) are highlighted in the corresponding color. While there is a maximum gain of DNAm at a CpG in PDE4C, the modifications are lower overall and not equal across all CpGs in the target region. c,d, Cumulative distribution functions of DNAm profiles at CpGs that gain or lose DNAm with aging in T cells. Bystanders increase for strata with higher age-relatedness (positive correlation: P < 10⁻⁸; NS, negative correlation; Kruskal–Wallis test).

Extended Data Fig. 1. DNA methylation in amplicons of the target region.

a) DNA methylation at the most age-associated CpG site in PDE4C was measured by digital droplet PCR in peripheral blood, as described before (n = 80, R² = 0.76; Han et al., BMC Biology 18,71 (2020)). b) Schematic presentation of the genomic region in PDE4C with guide RNA binding sites (blue) and the bisulfite amplicon sequencing region with corresponding methylation changes. c) Pyrosequencing measurement of HEK293T cells after epigenetic modification at PDE4C (three replicates). DNAm levels at the seven CpGs in the pyrosequencing amplicon are depicted at day 14 (***P < 0.001; **** P < 0.00001). d) Pyrosequencing of the same samples over a time course of up to 100 days confirmed that the modified DNAm is stable at this site. e) Correlation of DNAm was analyzed at the 26 neighboring CpGs within the bisulfite barcoded amplicon sequencing reads. The heatmap depicts Pearson correlations between individual CpGs. f) Heatmaps showing the complete time-course experiment upon modification with dCAS9-DNMT3A, corresponding to Fig. 1e.

Extended Data Fig. 2. Further analysis of epigenetic bystander modifications.

a,b) Scatterplots of Illumina BeadChip data demonstrate DNAm changes in HEK293T cells upon targeted methylation at PDE4C with either dCas9-DNMT3A (a; n = 3), or CRISPRoff (b; n = 3). Significant CpGs are highlighted in analogy to Fig. 2a, b. c,d) Due to sequence homology, 123 genomic regions were identified as potential off-target binding sites. Within a window of 5 kb, we identified 300 CpGs in these regions. Predicted off-target binding sites for PDE4C gRNAs did not reveal clear enrichment in bystander modifications of dCas9-DNMT3A (c), or CRISPRoff (d). e) Enrichment of the bystander effects in relation genomic regions: TSS1500: 1500 bp upstream of transcription start site (TSS); TSS200: 200 bp upstream of TSS; UTR: untranslated region. Enrichment was calculated in relation to all CpGs on the array and was highly significant for all categories (Chi² test). f,g) Cumulative distribution function analog to Fig. 2k, l comparing DNAm profiles of all CpGs, age-hypermethylated CpGs and a matched control, which has the same distribution of beta-values as age-hypermethylated sites in HEK293T control cells. Enrichment is also highly significant accounting for beta-values (Fig. f: P = 7.57*10⁻¹⁰; Fig. g: P = 2.59*10⁻¹⁰, two-sided, two-sample Kolmogorov–Smirnov-test).

Extended Data Fig. 3. Epigenetic editing at FHL2.

a) Methylation at the CpG site cg22454769 in FHL2 reveals a strong correlation with chronological age in blood, as previously described (n = 80, R² = 0.72; Han et al., BMC Biology 18,71 (2020)). b) The CpG site cg22454769 in FHL2 (centered at genomic position 0) revealed gain of DNAm upon targeting PDE4C with either CRISPRoff or dCAS9-DNMT3A (Illumina BeadChip data; both n = 3). c) The genomic region in FHL2 was targeted with dCAS9-DNMT3A in HEK293T cells and successful gain in DNAm was demonstrated by pyrosequencing after 14 days (cg22454769 corresponds to CpG 2; n = 2). d) Scatter-plot of EPIC BeadChip data confirming the specificity of epigenetic editing at FHL2 (CpGs within a 10 kbp window around cg22454769 are highlighted; n = 2). e,f) Gaussian kernel density estimate and cumulative distribution function showing enrichment of bystander effects at 5076 age-hypermethylated and 4127 age-hypomethylated sites (P < 10⁻¹⁵ for age-hypermethylated and P = 0.0009 age-hypomethylated CpGs; two-sided, two-sample Kolmogorov–Smirnov-test). g) Scatter-plot of EPIC datasets comparing methylation changes as compared to corresponding controls after targeting either PDE4C or FHL2 with dCAS9-DNMT3A (R = 0.32; P < 10⁻¹⁵, Pearson correlation test). The dotted lines indicate at least 10% gain of DNAm upon targeting either PDE4C or FHL2, and 488 CpGs passed both thresholds.

Extended Data Fig. 4. Further analysis of ATAC-seq and 4C-sequencing data.

a) Violin-plot of ATAC-seq signals comparing chromatin accessibility at all CpGs and hypermethylated bystander positions of FHL2 (delta DNAm >10%; n = 2, P < 10⁻¹⁵, unpaired t-test), in analogy to Fig. 3a. b) Multivariate kernel density estimate of ATAC-seq and DNA methylation changes (dCAS9-DNMT3A). While differentially methylated positions exhibited in tendency higher chromatin accessibility the overall correlation was rather low. In analogy, we analyzed the 4389 age-hypo and 5328 age-hypermethylated CpGs. c) Cis-interactions (on the same chromosome) of 4C-Sequencing with around 20 million reads per sample. Particularly the two viewpoints at positively age-correlated regions in PDE4C and FHL2 exhibited cis-interactions spanning a mean of 271 and 672 kbp, respectively. In contrast, the 4 C reads of age-hypomethylated region in MEIS1-AS3 comprised only 89 kbp. d) Trans-interactions (across different chromosomes) of 4C-sequencing. The viewpoints at positively age-correlated regions like PDE4C and FHL2 exhibited trans-interactions spanning a mean of 6.5 and 13.0 mbp, respectively. e) Bystander effects of FHL2 (abs. difference beta-value > 0.1; n = 2) become enriched with increasing coverage in 4C-sequencing of FHL2 interacting sites (n = 3, different symbols correspond to replica; R² = 0.81, linear model P < 10⁻¹⁵).

Extended Data Fig. 5. Pyrosequencing validation of targeted epigenetic modification.

a) Genomic tracks of differential DNAm at day 7 in a 10k base window centered at the corresponding target site (dashed vertical line) in the five targeted genes. b) Pyrosequencing measurements of demonstrating stable DNAm levels at PDE4C and FHL2 (triplicates). HEK293T cells were either transfected guides for PDE4C, FHL2, ELOVL2, KLF14, and TEAD1 together with dCas9-DNMT3A, or scramble guide RNA and dCas9-DNMT3A. The epigenetic modifications both regions were stable for at least 17 days. DNAm was significantly higher in experimental conditions across all timepoints (P value < 0.05, unpaired t-test). c) Gaussian kernel density estimate comparing the entirety of CpGs from the BeadChip with age-associated CpGs (4389 age-hypo and 5328 age-hypermethylated CpGs). In these experiments hardly any reproducible bystander effects were observed. d) Transcriptomic changes upon multiplexed epigenetic editing at the age-hypermethylated CpGs. 10 genes were significantly up- and 11 genes significantly down-regulated (adjusted P < 0.05; >2 fold change).

Extended Data Fig. 6. Targeted modifications at age-hypomethylated regions.

a) Genomic tracks of differential DNAm at day 3 in a 10k base window centered at the corresponding target site (dashed vertical line) at five CpGs that lose DNAm with age. b) Kernel density estimate revealing no methylation changes at 692 predicted off-target binding sites (sequence homology) compared to all CpGs on BeadChip. c) Multivariate kernel density estimate and linear regression (Pearson correlation) comparing DNAm changes with epigenome editing between PDE4C and hypo-CpGs.

Extended Data Fig. 7. Association of bystander effects in T cells.

a,b) Cumulative distribution functions comparing DNAm profiles at H3K27me3 peaks (GSE176621: H3K27me3 Mint-ChIP-seq of T cell male adult; 69806 CpGs) and all other CpGs (P values < 10⁻¹⁵, two-sided, two-sample Kolmogorov–Smirnov-test). c,d) Cumulative distribution functions comparing DNAm profiles at the top 4764 CpGs associated to EZH2 binding sites (GSE253773: EZH2 CHIP-seq of CD4 T cells). There are hardly any changes in DNAm patterns at EZH2 target sites. e-h) Cumulative distribution functions corresponding to Fig. 5k–n. The bystander modifications upon epigenetic editing at either age-hypermethylated CpGs or age-hypomethylated CpGs were clearly enriched at CpGs with increasing positive correlation (e,g), or decreasing negative correlation (f,h) with chronological age (R).

Extended Data Fig. 8. Enrichment of bystander modifications in mesenchymal stromal cells.

a,b) Transcriptomic changes upon epigenetic editing in MSCs was analyzed for two of the three donors in independent experiments. a) The scatter-plot depicts reads-per million (RPM) in control versus experiment. b) Comparison of DNAm changes (n = 3) and transcriptomic changes (n = 2) revealed no clear association between DNAm changes and gene expression changes. c,d) In our previous work, we have identified 646 CpGs that gain and 2442 CpGs that lose DNAm during culture expansion of MSCs. Since these CpGs may be associated with replicative senescence, they are referred to as senescence-hyper and senescence-hypo, respectively. We have analyzed if these CpGs are enriched in bystander modifications upon epigenetic editing at age-hypermethylated sites in MSCs. a) Gaussian kernel density estimate and b) cumulative distribution function demonstrates that bystander effects are slightly enriched in senescence-hyper (P < 10⁻⁰⁶) and senescence-hypo-CpGs (P < 10⁻⁹, two-sided, two-sample Kolmogorov–Smirnov-test). e,f) Gaussian kernel density estimate and cumulative distribution function comparing the entirety of CpGs from the BeadChip with age-associated CpGs in blood. Bystander effects are significantly enriched at age-hypermethylated sites (P < 10⁻¹⁵, two-sided, two-sample Kolmogorov–Smirnov-test). g,h) Density plots of gaussian kernel density estimate of bystanders at CpGs that gain or lose DNAm with aging in T cells, corresponding to Fig. 6c, d.

Extended Data Fig. 9. Hypothetical model of the epigenetic aging network.

The findings of this study indicate that epigenetic editing at age-hypermethylated CpGs, utilizing dCAS9-DNMT3A or CRISPRoff, leads to epigenetic bystander modifications that are notably enriched at other age-associated CpG sites—specifically, those CpGs that either gain or lose DNA methylation with aging. A similar phenomenon was observed when targeting hypomethylated CpGs. These bystander modifications appear to be concentrated in open chromatin regions within interacting chromatin domains. This suggests the existence of an epigenetic network in which differentially methylated regions can communicate with one another. Although the exact mechanism underlying this interaction remains unclear, such connections may contribute to the stabilization of the epigenetic state and facilitate coherent epigenetic modifications. At the single-cell level, we propose that age-associated DNAm is influenced not only by stochastic changes at individual CpGs but also by a genome-wide crosstalk of epigenetic modifications.