Multi-omic rejuvenation of naturally aged tissues by a single cycle of transient reprogramming.

The expression of the pluripotency factors OCT4, SOX2, KLF4, and MYC (OSKM) can convert somatic differentiated cells into pluripotent stem cells in a process known as reprogramming. Notably, partial and reversible reprogramming does not change cell identity but can reverse markers of aging in cells, improve the capacity of aged mice to repair tissue injuries, and extend longevity in progeroid mice. However, little is known about the mechanisms involved. Here, we have studied changes in the DNA methylome, transcriptome, and metabolome in naturally aged mice subject to a single period of transient OSKM expression. We found that this is sufficient to reverse DNA methylation changes that occur upon aging in the pancreas, liver, spleen, and blood. Similarly, we observed reversion of transcriptional changes, especially regarding biological processes known to change during aging. Finally, some serum metabolites and biomarkers altered with aging were also restored to young levels upon transient reprogramming. These observations indicate that a single period of OSKM expression can drive epigenetic, transcriptomic, and metabolomic changes toward a younger configuration in multiple tissues and in the serum.

Oct4/Sox2/Klf4/Myc

intraperitoneal injection

reduced representation bisulfite sequencing

differentially methylated

principal component analyses

differentially methylated regions

differentially expressed genes

minichromosome maintenance

NFE2 like BZIP transcription factor 2

aspartate aminotransferase

alanine aminotransferase

heat shock transcription factor 4

neutrophil‐derived extracellular traps

INTRODUCTION

The simultaneous expression of four specific factors, OCT4, SOX2, KLF4, and MYC (OSKM), also known as “Yamanaka factors,” in adult differentiated cells is able to shut off their transcriptional programs for cell identity, activate the transcription of pluripotency genes, and establish a new identity equivalent to embryonic stem cells (Takahashi & Yamanaka, 2006). This process, generally known as reprogramming, erases molecular and cellular traits of aging acquired by somatic cells throughout their lifespan (Lapasset et al., 2011; Liu et al., 2011; Mahmoudi & Brunet, 2012; Petkovich et al., 2017; Rando & Chang, 2012). Moreover, reprogrammed cells can subsequently differentiate into somatic cells that are now rejuvenated relative to their parental cells (Lapasset et al., 2011; Liu et al., 2011). Among the various molecular changes associated with aging, DNA methylation at specific CpG sites has turned out to be tightly linked to aging and, particularly, to biological aging rather than chronological aging (Field et al., 2018; Horvath & Raj, 2018). Interestingly, examination of aging‐associated DNA methylation has revealed that rejuvenation of this aging trait occurs progressively during reprogramming, being initiated at the early stages of the process and continuing until full reprogramming (Gill et al., 2021; Olova et al., 2019).

The expression of OSKM in mice recapitulates the process of reprogramming (Abad et al., 2013; Mosteiro et al., 2016; Ohnishi et al., 2014). Upon OSKM expression in vivo, a fraction of cells within tissues shut off their cell identity markers and progressively activate the pluripotency program (Abad et al., 2013). The completion of reprogramming in vivo manifests by the formation of teratomas, a tumor overgrowth formed by pluripotent cells differentiating into multiple cell lineages (Abad et al., 2013; Ohnishi et al., 2014). Of note, interruption of the process of reprogramming at its early stages is fully reversible and does not result in a detectable risk of teratoma. Remarkably, cycles of short OSKM expression followed by recovery result in rejuvenation, both in vivo and in vitro (Ocampo et al., 2016; Sarkar et al., 2020). Rejuvenation by multiple cycles of OSKM has been demonstrated at various levels. In cells, there is a reduction in aging‐associated DNA damage and epigenetic changes, including DNA methylation at specific CpG sites (Gill et al., 2021; Ocampo et al., 2016; Sarkar et al., 2020). In mice, cycles of OSKM in adult individuals improve their capacity to respond to tissue injury (Chen et al., 2021; Ocampo et al., 2016). Also, viral transduction of OSK in the retina of old mice has been shown to reduce aging‐associated DNA methylation and to improve vision (Lu et al., 2020). Finally, cycles of OSKM in progeric mice with constitutive DNA damage extend significantly their lifespan (Ocampo et al., 2016).

Here, we perform a multi‐omic and multi‐tissue analysis of naturally aged mice exposed to a single cycle of transient OSKM expression. By studying the effects of a single cycle of OSKM, we aim to capture the more direct effects of transient reprogramming. We also examine, not only DNA methylation changes associated with aging, but also transcriptomic and serum metabolites. Finally, we observe in vivo rejuvenation in the pancreas where OSKM is highly expressed, but, interestingly, also in liver, spleen, and peripheral blood where OSKM is weakly expressed.

RESULTS

OSKM promotes epigenetic rejuvenation in pancreas

To induce transient in vivo reprogramming, we used a previously reported strain of mice in which the expression of OCT4, SOX2, KLF4, and MYC (OSKM) can be temporarily induced by doxycycline supplementation in the drinking water (Abad et al., 2013; Mosteiro et al., 2016). In these mice, the pancreas is the most susceptible tissue to undergo reprogramming, followed by the intestine and stomach (Abad et al., 2013). Reprogrammed tissues present focal areas of dysplasia in which the tissue alters its normal architecture and cells lose differentiation markers (Abad et al., 2013). At a more advanced stage, tissues present compact areas of undifferentiated cells expressing pluripotency markers, such as NANOG (Abad et al., 2013). To evaluate the effects of OSKM activation in old mice, we treated 55 weeks old reprogrammable mice for one week with a low dose of doxycycline (0.2 mg/ml). The age of 55 weeks was chosen in an effort to detect age‐related differences caused by functional decline, rather than to alterations in cell composition that often occur at older ages; and the specific conditions for the treatment with doxycycline were chosen to avoid the formation of teratomas according to our previous experience (Abad et al., 2013). Upon OSKM induction for 1 week, we observed clear histological changes in the pancreas (Figure S1a). This treatment, however, was not sufficient to achieve pluripotency, as judged by minimal or undetectable expression of pluripotency markers Nanog, Oct4 (endogenous), and Tfe3 in pancreas (Figure S1b,c). Of note, the detection of pluripotency markers in pancreas required two weeks of OSKM expression (Figure S1b). Importantly, the histological changes observed in the pancreas after one week of OSKM induction were reversed two weeks after removal of doxycycline (Figure S1a). Therefore, one week of OSKM expression allows for transient and reversible changes in the histology of the pancreas without achieving pluripotency. It is important to mention that while pancreas, intestine, and stomach manifest histological changes during OSKM expression, other tissues, like liver or spleen, do not present observable histological changes except for signs of extramedullary hematopoiesis in the spleen (Figure S1a). It is also worth mentioning that the OSKM cassette is not homogeneously expressed in all cell types. For example, acinar cells of the pancreas rapidly and broadly express SOX2, 24h after i.p. injection of doxycycline, as detected by immunohistochemistry; however, this was not the case of other cell types of the pancreas, including most endocrine cells (Figure S1d).

Changes in the epigenome, and particularly DNA methylation, are robustly linked to aging (Field et al., 2018; Horvath, 2013; Horvath & Raj, 2018; Stubbs et al., 2017). To address the impact of transient OSKM expression on epigenetic aging in vivo, we used MspI‐based Reduced Representation Bisulfite Sequencing (RRBS) (Meissner et al., 2005). We performed this analysis on genomic DNA from pancreata of reprogrammable mice at the age of 55 weeks that had been treated as described above, referred to as “old‐OSKM” group (n = 5). As control groups (n = 5 per group), we used reprogrammable mice of 55 weeks of age (“old” group) and of 13 weeks (“young” group) without treatment with doxycycline (Figure S1e). Of all the methylation sites probed by RRBS, we focused on those located at regulatory elements, particularly promoters and enhancers. In total, we identified 11.272 promoters, which represent about 35% of the active promoters in pancreas (Shibata et al., 2018; Yue et al., 2014), defined as H3K27ac‐rich regions around transcription start sites; similarly, we identified 5.737 enhancers, defined as non‐promoter H3K27ac‐rich regions (Shibata et al., 2018; Yue et al., 2014), which represent about 42% of the active enhancers in pancreas. By comparing old and young groups, we identified a list of differentially methylated (DM) promoters (Figure S1f). Focusing on these aging‐DM promoters, we performed principal component analyses (PCA) that clearly segregated young from old mice. Interestingly, the methylation profile of old‐OSKM mice was placed by PCA between the young and old groups, revealing partial epigenetic rejuvenation of aging‐sensitive promoters (Figure 1a). Furthermore, we divided these aging‐DM promoters into two groups depending on their age‐associated gain or loss of methylation, respectively. This yielded a subset of promoters that are hypermethylated with aging and demethylated by OSKM (19 out of 51 in total, 37%), and a subset of promoters that are hypomethylated with aging and remethylated by OSKM (17 out of 42 in total, 40%) (Figure 1b, Figure S1g‐h, Table S1). As a notable example, the Hnf1a promoter, a key transcription factor for the development of the pancreas (Lau et al., 2018), was among the promoters hypermethylated with aging and OSKM induced a trend of partial demethylation (Figure 1b).

FIGURE 1

Transient OSKM reprogramming partially rejuvenates the methylation profile of old pancreas. (a) Principal component analysis (PCA) of aging‐related differentially methylated (DM) promoters of young (13 weeks, n = 5), old (55 weeks, n = 5), and old‐OSKM (55 weeks, n = 5) pancreas. (b) DM promoters were divided into hyper‐ and hypomethylated during aging, and shown is the number of these promoters that alters their methylation profile due to transient OSKM activation. Hnf1a promoter is a representative example of an age‐associated hypermethylated promoter that becomes demethylated in old‐OSKM pancreas. (c) PCA of aging‐associated DM enhancers of young, old, and old‐OSKM pancreas. (d) DM enhancers were classified into hyper‐ and hypomethylated during aging, and shown is the subset of these enhancers that alters their methylation profile due to transient OSKM activation. Arid5b enhancer is a representative example of an age‐related hypermethylated enhancer that becomes demethylated in old‐OSKM pancreas. (e) The methylation status of aging‐hypermethylated or hypomethylated promoters and (f) enhancers that were found to be OSKM‐demethylated or remethylated, respectively, was evaluated immediately after OSKM cessation (day 7) or after 14 days of recovery (day 21). (g) Methylation levels measured by bisulfite pyrosequencing of four CpGs, located in regions hypermethylated with aging in the pancreas (n = 4 to 6). Bars in b, d, and g represent the standard deviation (SD) of the data. Statistical significance was evaluated using one‐way ANOVA with Tukey's multiple comparison method, and comparisons are indicated as *p < 0.05 and **p < 0.01

Similarly, we generated a list of differentially methylated enhancers by comparing the methylation levels of young versus old pancreata (Figure S1i). As before, PCA analysis showed that aging‐DM enhancers in OSKM mice were distinct from old mice and displaced toward the young group (Figure 1c). We then separated those enhancers that are hypermethylated with aging and demethylated by OSKM (17 out of 43 in total), and those enhancers hypomethylated with aging and remethylated by OSKM (8 out of 32 in total) (Figure 1d, Figure S1j‐k, Table S2). For example, DM enhancer at Chr10: 68231400–68232600 showed a trend to be demethylated after OSKM activation (Figure 1d). Among the potential target genes in proximity to this enhancer, it is Arid5b that encodes a pro‐inflammatory RNA‐binding protein induced by NF‐κB (Nyati et al., 2020).

We wondered to what extent the observed changes in methylation occurred during the period of OSKM expression (1 week) or during the post‐recovery period (2 weeks). For this, we analyzed the aging‐DM regions in a group of 55 weeks old reprogrammable mice at the end of OSKM expression (d7) and compared it with the same regions after recovery (d21). We observed that the majority of OSKM‐induced demethylation in the pancreas happened after turning off OSKM expression, that is, during the 2 weeks of recovery period. In contrast, remethylation events were already present at d7 and half of them were preserved during recovery, while the other half were lost (Figure 1e,f).

To confirm the effects of OSKM on methylation with a different technique, individual CpG sites were selected among the above‐identified DM regions for validation by bisulfite pyrosequencing. For this validation, we used very old mice (around 100 weeks): wild‐type mice treated with doxycycline (“very old,” n = 6), and reprogrammable mice treated with doxycycline for one week followed by two weeks of recovery (“very old‐OSKM+2w,” n = 4) or four weeks of recovery (“very old‐OSKM+4w,” n = 4) (Figure S1e). We selected those individual CpGs within the RRBS with the highest methylation changes with aging and optimal sequence context for pyrosequencing (Table S3). We tested a total of eleven CpGs. Among them, nine CpGs were expected to gain methylation with aging (Figure 1g, Figure S1l,m). Out of them, eight showed the expected increase in methylation with aging, specifically, those near genes Tpt1 (2 close CpGs), Stambpl1, Dnmbp (2 consecutive CpGs), Akt1, Gm12339, and Gm17678. Interestingly, all eight regions presented reduced methylation in very old‐OSKM mice treated with doxycycline, in some cases becoming indistinguishable from young mice (Figure 1g, Figure S1l). Of note, reversion of methylation was consistently more profound after 4 weeks of recovery than after 2 weeks of recovery. This, together with the RRBS data above (Figure 1e,f), is another indication that the recovery period is important for the demethylation of aging‐DMRs. We also tested two aging‐hypomethylated CpGs (1700011F14Rik and Ptprj) and pyrosequencing confirmed that both were hypomethylated with aging. However, we could not detect remethylation in these two positions after transient OSKM activation (Figure S1n).

In summary, RRBS analysis has revealed a total of 93 promoters and 75 enhancers in which methylation changes with aging, and out of them a total of 61 (36%) were reversed toward a younger state by OSKM expression in old mice (Figure 1b,d). Using a separate cohort of old mice and a different technique (pyrosequencing), we confirmed rejuvenation in 8 out of 11 tested regions that showed aging‐dependent methylation changes. We conclude that a single period of OSKM expression is able to rejuvenate a fraction of the methylation changes that occur with aging in promoters and enhancers.

Transcriptional rejuvenation of the pancreas

To analyze the effects of aging and OSKM at the transcriptional level, we performed RNAseq of the pancreata of young (13 weeks, n = 4), old (55 weeks, n = 5), and old‐OSKM (55 weeks, 1 cycle of OSKM followed by 2 weeks of recovery, n = 4) (Figure S1e). By comparing the transcriptome of old versus young samples, we identified a list of significantly differentially expressed genes (aging‐DEGs, n = 217) (Table S4). PCA analysis of these aging‐associated DEGs clearly separated old and young samples as expected. Remarkably, PCA placed old‐OSKM profiles between old and young ones, suggesting that transient OSKM activation partially restores aging‐related transcriptomic changes (Figure 2a). To visualize the behavior of individual genes, we plotted the aging‐DEGs and we colored those genes that were affected by OSKM. Interestingly, the large majority of genes upregulated by OSKM in old mice (pink dots) corresponded to genes with reduced expression upon aging (Figure 2b, Figure S2a). Conversely, those genes downregulated by OSKM in old mice (blue dots) were genes with upregulated expression upon aging (Figure 2b, Figure S2a). Looking then at the whole transcriptome, we first interrogated gene‐sets well‐established to change with aging, such as mTOR (López‐Otin et al., 2013) and DNA replication (Dumit et al., 2014). As expected, the mTOR signaling gene‐set was upregulated in old versus young control mice (Figure 2c), whereas the DNA replication gene‐set was downregulated (Figure 2d). Notably, old‐OSKM samples behaved like young samples when compared to untreated old control mice (Figure 2c‐d, Table S5). To better analyze the behavior of gene‐sets across the three groups of samples, we performed a pattern analysis based on a Normal‐Normal hierarchical model (gaga) (Rossell, 2009) to identify those gene‐sets that (i) change significantly between young and old samples, and (ii) are similarly expressed in young and old‐OSKM samples. This analysis identified a total of 179 gene‐sets rejuvenated by OSKM (see pattern 1 in Table S6). When the same analysis was performed after randomizing the samples (i.e., samples were randomly assigned to the three experimental groups: young, old, and old‐OSKM), only 44 gene‐sets were found (see pattern 1 in Table S7). More importantly, the gene‐sets rejuvenated by OSKM included important aging‐related processes, such as mTOR upregulation, insulin increase, reduction of NADPH and pyrimidine synthesis, and decline of mitochondrial processes, such as fatty acid oxidation, tricarboxylic acid cycle, and oxidative phosphorylation (Brink et al., 2009; Fontana et al., 2010; López‐Otin et al., 2013) (Figure 2e, Figure S2b‐f). All these gene‐set alterations were ameliorated in old‐OSKM mice (Figure 2e, Figure S2b‐f). DNA replication and repair are known to be reduced with aging (Gorbunova et al., 2007; Moskalev et al., 2013). We observed that gene‐sets related to the minichromosome maintenance (MCM) helicase, the DNA replication machinery, mismatch repair, and base excision repair, were all reduced in our old mice and upregulated to young levels in old‐OSKM mice (Figure 2f, Figure S2g‐h). Another important feature of aging is impaired protein homeostasis (Powers et al., 2009). Again, this process was improved in old‐OSKM mice (Figure 2g, Figure S2i). Finally, loss of collagens occurs with age in pancreas (Riopel & Wang, 2014) and, remarkably, old‐OSKM mice increased their levels of collagen‐related genes (Figure 2h, Figure S2j). Overall, transient OSKM activation appears to orchestrate a positive reconfiguration of the transcriptome that opposes key hallmarks of aging.

FIGURE 2

Transient OSKM reprogramming rejuvenates the transcriptome of old pancreas. (a) Principal component analysis (PCA) of aging‐related differentially expressed genes (DEGs: fold change > 1.5 and raw pval < 0.01) including young, old, and old‐OSKM pancreas. (b) Representation of these aging‐DEGs (217 genes in total) colored by their alteration of expression induced by OSKM: OSKM‐upregulated genes are depicted in pink and OSKM‐downregulated genes are depicted in blue, while the names of the top ten genes, either upregulated or downregulated with aging, are also depicted (c) Enrichment analysis based on ROAST (Efron & Tibshirani, 2007; Wu et al., 2010) was performed comparing young versus old, and old versus old‐OSKM pancreas depicting: (c) mTOR signaling pathway (KEGG_04150) and (d) DNA Replication (KEGG_03030). Statistical significance in c and d was evaluated using Kolmogorov–Smirnov test. (e‐h) Z‐score representation of the expression profile with GS adjustment for the indicated gene‐sets among young, old, and old‐OSKM pancreas. Gene‐sets have been selected for following this pattern: Young ≠ Old & Young = Old OSKM after a Normal‐Normal hierarchical model (gaga) (Rossell, 2009) (5% FDR)

Evidence of rejuvenation in tissues with low reprogramming

Given its rejuvenating potential at the epigenetic and transcriptomic levels in the pancreas, we wondered whether other tissues that are modestly affected by OSKM expression would nevertheless present some indication of reversion of the epigenomic changes associated to aging. For this, we performed a similar methylation profiling by RRBS analysis in the liver (5.3% of the promoters and 3.5% of the enhancers were covered) and spleen (25% of the promoters and 22.4% of the enhancers were covered), two tissues that modestly upregulate OSKM expression upon doxycycline treatment (Figure S1c). As before, we identified a list of DM promoters and enhancers during aging. Similar to the pancreas, PCA analysis of the methylation patterns of liver DM promoters and enhancers indicated that old‐OSKM mice partially recovered a younger methylation pattern although the separation between experimental groups was not as clear as in the case of the pancreas (Figure 3a). From the total combined number of aging‐DM promoters and enhancers (n = 108), a substantial fraction (61%) underwent rejuvenation (Figure 3b, Figure S3a‐f, Table S8). This included Foxa3 which serves as a pioneer transcription factor for the maintenance of liver‐specific transcription (Iwafuchi‐Doi et al., 2016) (Figure 3b). Other notable promoters were those of Hoxd10, a tumor suppressor gene whose promoter hypermethylation has been linked to hepatocellular carcinoma (Guo et al., 2017), and Thy1 whose expression stimulates liver regeneration (Ichinohe et al., 2017) (Figure 3b). In the pancreas, we observed a temporal pattern for OSKM‐induced methylation changes (remethylation during OSKM expression and demethylation during the recovery period). However, this temporal pattern was not evident in the liver (Figure S3g).

FIGURE 3

Old livers present rejuvenated features after transient OSKM reprogramming. (a) Principal component analysis (PCA) of aging‐related differentially methylated (DM) promoters including young, old, and old‐OSKM livers. (b) DM promoters are divided into hyper‐ and hypomethylated during aging, and shown is the number of these promoters that alter their methylation profile due to transient OSKM activation. Hox10, Foxa3, and Thy1 promoters are representative examples of age‐associated hypermethylated promoters that become demethylated in old OSKM livers. (c) The expression of global aging genes identified by mouse Aging Cell Atlas (Zhang et al., 2021) was evaluated in very old livers (100 weeks; group 1 consists of 5 wild‐type mice as control, 5 reprogrammable mice activating OSKM for 1 week and recovering for 2 weeks, and 5 reprogrammable mice activating OSKM for 1 week and recovering for 4 weeks) compared to young (13 weeks; n = 5) control livers. (d) The expression of Nrf2 and ApoM, two aging‐associated genes, was measured in very old livers (100 weeks; group 1 as above together with group 2 consisting of other 5 wild‐type mice as control and 7 reprogrammable mice activating OSKM for 1 week and recovering for 4 weeks). (e) The levels of aspartate aminotransferase (AST/GOT) and alanine aminotransferase (ALT/GPT) were measured in the serum of wild‐type (n = 7) and reprogrammable (n = 6) mice before doxycycline treatment, as well as after 1 week of doxycycline and 2 weeks of recovery. Statistical significance was evaluated using one‐way ANOVA with Tukey's multiple comparison method, and comparisons are indicated as *p < 0.05, **p < 0.01 and ***p < 0.001

To obtain further insights into the aging‐associated transcriptome of old‐OSKM livers, we examined the top genes of a recently reported aging signature based on the Mouse Aging Cell Atlas (Tabula Muris Consortium, 2020) that applies to multiple tissues including the liver and spleen but not to the pancreas (Zhang et al., 2021). In a pilot test, we measured by qRT‐PCR a total of 15 mRNAs from this signature in young and old livers. Five of these genes were confirmed to be downregulated in the liver of our old mice (Figure 3c). Interestingly, all these five genes recovered young levels of expression in the livers of old‐OSKM mice (Figure 3c). These findings were further validated in an independent cohort of very old mice (100 weeks) (Figure S3h). Finally, we tested the levels of two other genes associated with aging, namely, Nrf2, a key regulator of cellular redox homeostasis (Schmidlin et al., 2019), and ApoM, a high‐density lipoprotein (HDL) that promotes vascular homeostasis (Ding et al., 2020). We confirmed that both genes are downregulated with aging in the liver of our very old mice, and this was reverted upon one cycle of OSKM (Figure 3d).

To associate the above findings with a physiological readout of liver function, we measured the levels of aspartate aminotransferase (AST/GOT) and alanine aminotransferase (ALT/GPT) in an independent cohort of very old mice (different from the two cohorts mentioned above). Interestingly, the serum levels of transaminases of OSKM mice were significantly lower after 1 week of OSKM activation and 2 weeks of recovery reflecting an improved liver function (Figure 3e).

Senescent cells increase with aging and may account for up to 17% of the hepatocytes in extremely old mice (~3 years of age) (Wang et al., 2009). Cellular senescence is a potent barrier for reprogramming (Haridhasapavalan et al., 2020), and we wondered if our transient reprogramming protocol could reduce cellular senescence in vivo. For this, we examined the senescence marker p16INK4a (Cdkn2a) and the senescence‐associated cytokines Mcp1 and Cxcl2, known to increase with aging in mouse liver (Yousefzadeh et al., 2020). As expected, all of these were increased with aging; however, their levels did not decline in any of the two cohorts of very old‐OSKM mice (Figure 3c, Figure S3h,i). The levels of p21CIP1 expression (Cdkn1a) were not informative because their levels did not change with aging nor with OSKM (Figure S3h,i). We also measured the number of γH2AX‐positive cells, a marker of in vivo senescence, and again the levels in the liver were similar in very old mice with or without 1 cycle of OSKM (Figure S3j). These observations suggest that a single period of transient OSKM in vivo may not be sufficient to rejuvenate the transcriptome of the senescent cells present in the aged liver.

In the spleen, PCA analysis of aging‐DM promoters and enhancers did not allow to infer an epigenetic rejuvenation of this tissue (Figure S4a,e), although up to 163 promoter and enhancer regions showed evidence of rejuvenation according to their average methylation level (Figure S4b‐d,f‐h, Table S9). To further explore the possibility of epigenetic rejuvenation in hematopoietic cells, we focused on the methylation of a specific intragenic region of the Hsf4 gene. This particular region has been identified by independent studies as differentially methylated during aging in various tissues (Beerman et al., 2013; Han et al., 2018; Taiwo et al., 2013), and it is extraordinary due to the large magnitude of change in percentage of methylation from young to old individuals (Han et al., 2018). First, we measured by pyrosequencing the methylation levels of the Hsf4 region in the blood of mice of different ages, confirming a linear change of large magnitude, from ~20% methylation in young mice (~10 weeks) to ~55% methylation in very old mice (~100 weeks) (Figure 4a). Then, we measured Hsf4 methylation in the blood of very old mice, reprogrammable and non‐reprogrammable littermates, before and after a period of 5 weeks (all mice were treated with doxycycline for 1 week followed by 4 weeks of recovery). Interestingly, it was possible to detect a 4% increase in Hsf4 methylation in control non‐reprogrammable mice (Figure 4b). In contrast, during the same period of time, reprogrammable mice reduced their average methylation levels by ~4%, resulting in a total difference of 8% methylation between the reprogrammable and the non‐reprogrammable mice (Figure 4b).

FIGURE 4

Evidences of OSKM‐induced rejuvenation in hematopoietic cells. (a) Correlation of the average methylation of 3 close CpGs located in the intragenic region of Hsf4 gene (Chr8; at the positions 105271000,105271005, and 105271015), as measured in the blood of mice, with their chronological age (in weeks) of the mice. (b) Δ change of the methylation levels in these CpGs of Hsf4 before and after a period of 5 weeks in very old (100 weeks) wild‐type mice (n = 4) and reprogrammable mice (n = 3). Both experimental groups were treated for 1 week with doxycycline and recovered for 4 weeks. Statistical significance was evaluated using Mann–Whitney nonparametric t‐test. (c) The expression of global aging genes identified by mouse Aging Cell Atlas (Zhang et al., 2021) was evaluated in very old spleens (100 weeks; 5 wild‐type mice as control, 7 reprogrammable mice activating OSKM for 1 week and 4 weeks of recovery) compared to young (13 weeks; n = 4) control spleens. (d) The expression of Nrf2, an aging‐associated gene, was measured in the same group of mice. Statistical significance was evaluated using one‐way ANOVA with Tukey's multiple comparison method, and comparisons are indicated as *p < 0.05, **p < 0.01, and ***p < 0.001

We also tested the above‐mentioned aging signature based on the Mouse Aging Cell Atlas in the spleen (Zhang et al., 2021). Interestingly, OSKM transient expression rescued the age‐associated decline of seven of these genes while the senescent marker p16INK4a (Cdkn2a) was unchanged (Figure 4c, Figure S4i). Finally, Nrf2 also declined in very old control spleens and its expression was reset to young levels in very old‐OSKM spleens (Figure 4d).

We conclude that a single cycle of OSKM expression has a detectable rejuvenating effect on tissues such as liver, spleen, and blood that do not express high levels of OSKM and do not manifest obvious signs of histological alterations. Conceivably, some of the observed effects on liver, spleen, and blood could be secondary to the direct rejuvenating actions of OSKM in other tissues, such as the pancreas.

Serum metabolomic profiling reveals systemic benefits

The aforementioned results indicating epigenetic and transcriptional rejuvenation in several tissues prompted us to identify possible systemic signs of rejuvenation in the serum. To address this, we performed mass spectrometry‐based metabolomics on the sera of female mice (to reduce sex‐related variations) from two independent experiments (each experiment analyzed separately by mass spectrometry). Each experiment included a group of young mice and a group of very old (~100 weeks) reprogrammable mice. The sera of very old reprogrammable mice were analyzed longitudinally, that is, before and after a single cycle of reprogramming (1 week of doxycycline and 2 or 4 weeks of recovery). When combining the two independent experiments, a total of 23 metabolites were identified as significantly changed between young and very old mice (Figure S5a,b, Table S10). Out of these aging‐associated metabolites, 4 were reversed after reprogramming (Figure 5a). These metabolites are 4‐hydroxyproline, thymine, trimethyl‐lysine, and indole‐3‐propionic acid, and some of them had been previously linked to aging. Serum 4‐hydroxyproline has been reported to decline with aging (Seo et al., 2016). This modified amino acid amounts to 14% of all the amino acids in collagen and its presence in the serum is considered to reflect total collagen levels (Gabr et al., 2017). Notably, our transcriptomic analyses identified collagen synthesis as a gene‐set downregulated with aging and rescued by OSKM in the pancreas (Figure 2h). Thymine has been found to extend lifespan in C. elegans (Wan et al., 2019). The other two metabolites, trimethyl‐lysine and indole‐3‐propionic acid, have not been studied in the context of aging.

FIGURE 5

Metabolomic analysis in the serum of very old mice reveals systemic beneficial effects upon transient OSKM activation. (a) Metabolomic analyses were performed on the sera of very old (100 weeks) female reprogrammable mice from two independent experiments (each experiment analyzed separately by mass spectrometry); Group 1: n = 6 and Group 2: n = 6. The sera of these mice were analyzed longitudinally, that is, before and after a single cycle of reprogramming. The Δ change of the levels of 4 metabolites that were identified to change with aging, being either upregulated (4‐hydroxyproline, thymine, trimethyl‐lysine) or downregulated (indole‐3‐propionic acid), are depicted after activating OSKM for 1 week and recovering for 2 to 4 weeks. Statistical significance was evaluated using a paired t‐test as values were confirmed to follow a normal distribution using the Shapiro–Wilk test. Comparisons are indicated as *p < 0.05, **p < 0.01, and ***p < 0.001

DISCUSSION

In this work, we report that a single cycle of transient OSKM expression in naturally aged mice can elicit epigenetic, transcriptomic, and serum metabolomic changes that partially restore younger patterns. We have focused on the effects of a single cycle of OSKM expression. The idea behind this choice was to identify and quantify early events induced by a single cycle of OSKM expression. Moreover, we have used a relatively low dose of doxycycline to avoid dramatic changes in cell identity thereby minimizing the risk of teratoma formation. This protocol of induction was sufficient to cause histologically detectable alterations in the pancreas, but not in the spleen or liver. Another relevant aspect of our experimental design is that we allowed for a recovery period of 2 to 4 weeks post‐OSKM expression. This recovery period was sufficient for a full restoration of normal histology in the pancreas. Conceivably, a fraction of the rejuvenation changes induced by OSKM may be unstable and therefore lost during the recovery period. Additionally, as we will argue below, some rejuvenation changes may indeed occur during the recovery period.

In the case of the pancreas, we assessed whether methylation changes in aging‐associated differentially‐methylated (DM) regions were already present after 1 week of OSKM expression. Notably, about half of the changes observed at the time of switching off OSKM disappeared after 2 weeks. This illustrates the existence of unstable changes induced by OSKM and the importance of allowing a recovery period to identify durable changes. We also observed interesting differences in the behavior of aging‐DM regions depending on whether the epigenetic rejuvenation consists on losses or gains of DNA methylation. The majority of gain‐of‐methylation events observed after recovery were already present at the time of extinguishing OSKM expression. In contrast, the majority of loss‐of‐methylation events were absent after OSKM expression and, therefore, were acquired during the recovery period. In line with these observations, pyrosequencing measurement of individual aging‐methylated CpGs consistently showed more demethylation 4 weeks post‐OSKM expression compared to 2 weeks. Together, we speculate that gains of methylation are closely linked to OSKM expression, whereas demethylation events are linked to the re‐establishment of normal chromatin during the recovery period. Conceivably, the reversion phase may be as important as the OSKM expression phase for molecular and cellular rejuvenation.

Analysis of the global transcriptome of the pancreas also allowed us to select for aging‐related differentially expressed genes (aging‐DEGs). Similar to the methylation data, one cycle of OSKM (1 week of induction and 2 weeks of recovery) was sufficient to reset the levels of aging‐DEGs to a younger state. Interestingly, transient OSKM activation reversed well‐known transcriptional changes associated to aging. This is the case of insulin and mTOR signaling, mitochondrial metabolic processes, and DNA repair. We did not observe however evidence of reduction in senescence markers. Senescent cells are generally resistant to be reprogrammed (Haridhasapavalan et al., 2020). Therefore, a single period of transient OSKM may not be sufficient to rejuvenate the transcriptome of the senescent cells present in aged tissues. This does not exclude the possibility that multiple cycles of OSKM expression may achieve, directly or indirectly, the reduction of senescent cells in aged tissues.

We have not been able to establish associations between specific aging‐DM regions and the mRNA levels of the associated genes. This lack of association between aging‐DM regions and transcription is, however, a general finding made by multiple researchers studying aging‐DM regions (see detailed discussion in (Horvath & Raj, 2018)). For example, aging‐DM regions often occur in bivalent, silent, promoters where gain of methylation does not impact the already silent transcriptional state (Horvath & Raj, 2018).

In the case of the spleen and liver, the methylation changes observed in aging‐DM regions were not as pronounced as in the case of the pancreas, as judged by principal component analyses. This is not entirely surprising if we consider the mild expression of the OSKM transgene in these tissues. Nevertheless, spleen and liver did show evidence of rejuvenation of a transcriptional signature of aging derived from the Mouse Aging Cell Atlas (Zhang et al., 2021). It seems, therefore, that transcriptional rejuvenation may be easier to achieve compared to epigenetic rejuvenation.

To evaluate molecular features of aging in the blood, we have focused on three closely located CpGs in an exon of the Hsf4 gene that have called the attention of at least three independent groups (Beerman et al., 2013; Han et al., 2018; Taiwo et al., 2013). A remarkable property of this methylated region is the large magnitude of change in percentage of methylation from young to old individuals (Han et al., 2018). Indeed, we found that the degree of methylation of these sites in peripheral blood is tightly linked to the age of the mice. We have measured the methylation levels of Hsf4 in the blood of very old mice (~100 weeks of age) before and after a period of 5 weeks. In the case of control mice, we observed an increase in methylation during this period of 5 weeks. This was in contrast to mice exposed to one cycle of reprogramming (1 week of OSKM expression followed by 4 weeks of recovery) in which Hsf4 methylation was decreased or maintained.

We also detected lower levels of serum transaminases (ALT/GPT and AST/GOT) after one cycle of reprogramming. Serum transaminases are generally considered an indication of liver injury. However, there are growing evidences suggesting that other tissues also contribute to the total levels of serum transaminases (Kobayashi et al., 2020).

As a readout of systemic rejuvenation, we focused on serum metabolites. We analyzed the serum metabolome in two cohorts of very old mice (~100 weeks of age) in a longitudinal design, that is, serum was analyzed before treatment and after one cycle of OSKM expression (1 week of doxycycline followed by 4 weeks of recovery). We have detected four serum metabolites that changed with aging and were reversed by OSKM in the two cohorts of old mice: 4‐hydroxyproline, thymine, trimethyl‐lysine, and indole‐3‐propionic acid. Serum 4‐hydroxyproline has been previously reported to decline with aging (Seo et al., 2016), and it is considered to reflect total collagen content in the organism (Gabr et al., 2017). Of note, the transcriptional network for collagen synthesis was downregulated with aging and rescued by OSKM in the pancreas. Thymine has been found to extend lifespan in C. elegans (Wan et al., 2019), and this is consistent with the aging‐associated reduction that we have observed for this metabolite and its reversion by OSKM. Little is known about trimethyl‐lysine and indole‐3‐propionic acid in connection with aging. Trimethyl‐lysine is a modification abundant in chromatin and we speculate that its presence in the serum may derive from neutrophil‐derived extracellular traps (NETs), a process by which neutrophils extrude their chromatin into the bloodstream (Hazeldine & Lord, 2018). Interestingly, the production of NETs by stimulated neutrophils declines with aging (Hazeldine & Lord, 2018), which would be consistent with a reduction in trimethyl‐lysine.

Finally, it should be mentioned that OSKM expression may have non‐cell‐autonomous effects that potentially underlie some of our observations. For example, molecular rejuvenation of a subset of cells may influence the rejuvenation of other non‐reprogrammed cells, locally or distantly, through modifications in the microenvironment including the secretion of soluble factors. This may also result in changes in cell composition that have not been addressed in this work.

In this work, we report stable reversion of molecular features of natural aging by a single cycle of OSKM expression. The observed rejuvenation features include DNA methylation, transcription, and serum metabolites. These molecular manifestations of rejuvenation likely reflect a complex mixture of direct and indirect effects on multiple tissues. We hope this serves as the basis for future studies to dissect the mechanisms underlying OSKM‐driven rejuvenation in vivo. Also, it may provide benchmarking to recapitulate OSKM‐like rejuvenation with pharmacological or nutritional interventions.

EXPERIMENTAL PROCEDURES

Animal procedures

Animal experimentation was performed between two different institutes: at the Spanish National Cancer Research Centre CNIO in Madrid and at the Institute of Research in Biomedicine IRB in Barcelona, according to protocols approved by the CNIO‐ISCIII Ethical Committee for Research and Animal Welfare (CEIyBA) in Madrid, and by the Animal Care and Use Ethical Committee of animal experimentation of Barcelona Science Park (CEEAPCB) and the Catalan Government in Barcelona. We used the reprogrammable mice known as i4F‐B which carries a ubiquitous doxycycline‐inducible OSKM transgene, abbreviated as i4F, and inserted into the Pparg gene (Abad et al., 2013). Mice of both sexes were used, and of different ages; young (females, 13 weeks), old (females, 55 weeks), and very old (males and females, 100 weeks). 0.2 mg/ml of Doxycycline hyclate BioChemica (PanReac) was administered in the drinking water supplemented with 7.5% sucrose for a period of 7 days. Mice were sacrificed two or four weeks after doxycycline removal. In the case of the intraperitoneal injection of doxycycline, young mice received 2.5 mg of doxycycline dissolved in 100 μl of saline and sacrificed 24 h later.

DNA isolation

All pancreas, liver, and spleen tissues were snap‐frozen directly after collection. Genomic DNA was extracted from frozen tissues using the DNeasy Blood & Tissue Kit (Qiagen). In the case of blood samples, genomic DNA isolation was conducted according to a standard phenol‐chloroform extraction protocol after red blood cell lysing. All samples were collected, processed for their methylation status, and further analyzed.

Reduced Representation Bisulfite Sequencing (RRBS) library preparation

RRBS libraries were prepared as described previously (Stubbs et al., 2017). 500ng of genomic DNA was digested with MspI (Thermo Scientific). Following digestion, DNA fragments were end‐repaired and T‐tailed with Klenow fragment lacking 5′ → 3′ and 3′ → 5′ exonuclease activity (NEB). In‐house adapters containing 8‐nucleotide unique molecular identifier (UMI) sequences were ligated onto the DNA fragments with T4 DNA ligase (NEB). Excess adapters were removed by AMPure XP beads (Agencourt, 0.8X ratio). Libraries were bisulfite converted with the EZ‐96 DNA Methylation‐Direct MagPrep kit and cleaned up using an automated liquid handling platform (Agilent Bravo). The bisulfite‐converted libraries were amplified for 12 cycles with KAPA HiFi HotStart Uracil+ReadyMix (Roche). The amplified libraries were cleaned up with AMPure beads (0.8X ratio) and quality checked using an Agilent High Sensitivity DNA kit on an Agilent 2100 Bioanalyzer system. Libraries passing quality checks were sequenced on an Illumina HiSeq 2500 sequencer for 75 bp paired‐end sequencing. Reads were trimmed with Trim Galore (version 0.5.0). Trimmed reads were aligned to the mouse genome (GRCm38) with Bismark (version 0.20.0). Aligned reads were deduplicated based on UMI sequence and chromosomal position with Umibam (version 0.0.1). Bismark coverage files were generated from the deduplicated bam files using Bismark Methylation Extractor.

RRBS DNA methylation analysis

Promoters were defined as −2000 bp to +500 bp of the transcription start site. Pancreas, spleen, and liver enhancers were defined using published H3K27ac data (Shibata et al., 2018; Yue et al., 2014). H3K27ac peaks were called using MACS and peaks that did not overlap promoter regions were classified as enhancers. Enhancers were linked to genes based on proximity. Enhancers were linked to the nearest gene within 1 Mb. The methylation analysis was carried out in Seqmonk. CpG sites were only carried forward in the analysis if they were covered by at last 5 reads in all the samples of that tissue type. The methylation level of each promoter and enhancer was calculated using the mean of all the methylation levels of the CpGs within those regions. Differential methylation analysis was carried out using logistic regression. Regions were identified as significant if their P‐value was below 0.05 and if they demonstrated a minimum methylation difference of at least 10% DNA methylation. Rejuvenated regions were those in which the mean methylation in the old‐OSKM group was closer to the mean of the young group than to the mean of the old group.

Bisulfite pyrosequencing

DNA methylation patterns were validated by bisulfite pyrosequencing. Bisulfite modification of DNA was performed with the EZ DNA methylation‐gold kit (Zymo Research) following manufacturer's instructions. The sets of primers for PCR amplification and sequencing were designed using the specific software PyroMark assay design (version 2.0.01.15). Mm10 genome was used for the alignments. The primers used were as follows: 1_133696298‐311 forward primer: 5′‐AGAGTGTTAGAGTTGGAGAGAT‐3′, 1_133696298‐311 reverse primer: 5′‐[Btn]AAAAAAACCTCTAACCTCCATATATC‐3′ (Sequencing: 5′‐TGGAGAGATTTTGAAGTT‐3′), 2_90532378 forward primer: 5′‐GGTTTGGAATTTGGTTTATGTATAGA‐3′, 2_90532378 reverse primer: 5′‐[Btn]CTCAAACCAAAAAACCCTAATCTCC‐3′ (Sequencing: 5′‐TTTTTTTTAAGAGAATAGGTTATA‐3′), 11_75813125 forward primer: 5′‐GGGTTAGGTTGAGTTTTTTAGAATGAAT‐3′, 11_75813125 reverse primer: 5′‐[Btn]ACCCTCTCCATCTATACCTACTCC‐3′ (Sequencing: 5′‐AGGGTTTTGTATTTAATTTTTATT‐3′), 12_112670818 forward primer: 5′‐AGAGGTTGGGATTGGTAAGGATT‐3′, 12_112670818 reverse primer: 5′‐[Btn]TTACCCCAAAACAAACATCTCACCC‐3′ (Sequencing: 5′‐GGTAAGGATTATTTTAGGGTT‐3′), 14_121078071 forward primer: 5′‐ATGAGATGATTTTAGTTAAGGATTTAGTT‐3′, 14_121078071 reverse primer: 5′‐[Btn]ATTACACAAAACTTCCAACTTACT‐3′ (Sequencing: 5′‐ATAGGAATAAAAGTTTTTTGATAA‐3′), 14_75862317‐319–323 forward primer: 5′‐GAAGGAAGGGAATTTTGAGATTTG‐3′, 14_75862317‐319–323 reverse primer: 5′‐[Btn]TCTCCCAAAACTATACCATCACCA‐3′ (Sequencing: 5′‐ TTTGTTTTTGTTTTTTTATTATAAG‐3′), 19_34218524 forward primer: 5′‐TGATTTGGTTAAAGGTAGAAAAGTAAGA‐3′, 19_34218524 reverse primer: 5′‐[Btn]AACCTTTATAAACAATCAATAAATATACCT‐3′ (Sequencing: 5′‐GGTTTTTGTGGATAAATATTA‐3′), 19_43890303‐291 forward primer: 5′‐GGATTTAGGTGGGTTTTATTTAGAAAATG‐3′, 19_43890303‐291 reverse primer: 5′‐[Btn]TCCCAATACCCACAATCCCTTTTT‐3′ (Sequencing: 5′‐GTTTTATTTAGAAAATGGTTTGGA‐3′), 19_58600208 forward primer: 5′‐AGAGGAAATAATTTTATAGTGTAGTAAGA‐3′, 19_58600208 reverse primer: 5′‐[Btn]CAAATTCTCAACCATAAAAATCACTCTA‐3′ (Sequencing: 5′‐TGGGTGGAAATGTGA‐3′), mHsf4 forward primer: 5′‐GTGYGTYGTAAGGTGGGATAAATTGTAGAAAAAATG‐3′, mHsf4 reverse primer: 5′‐[Btn]TCCRTACTCTCCTACACTCCTCTCAAAACTTA‐3′ (Sequencing: 5′‐ATGGTGTTTTTTGTTTGTAG‐3′). After PCR amplification, pyrosequencing and quantification were performed using PyroMark Q24 reagents, equipment, and software (Qiagen).

RNA isolation and analysis of mRNA levels

Total RNA was extracted from pancreas samples using guanidine thiocyanate, followed by acid phenol‐chloroform extraction. For the rest of the tissue samples (liver, spleen), total RNA was isolated with Trizol (Invitrogen), following provider's recommendations. Up to 5 μg of total RNA was reverse transcribed into cDNA using iScriptTM Advanced cDNA Synthesis Kit for RT‐qPCR (Bio‐Rad). Quantitative real‐time PCR was performed using GoTaq® qPCR Master Mix (Promega) in a QuantStudio 6 Flex thermocycler (Applied Biosystem) or 7900HT Fast (Applied Biosystem). For input normalization, we used the housekeeping genes Gapdh, β‐actin, and 18S. The primers used were as follows: Gapdh forward primer: 5′‐TTCACCACCATGGAGAAGGC‐3′, Gapdh reverse primer: 5′‐CCCTTTTGGCTCCACCCT‐3′; 18S forward primer: 5′‐ GTAACCCGTTGAACCCCATT‐3′, 18S reverse primer: 5′‐CCATCCAATCGGTAGTAGCG‐3′; β‐actin forward primer: 5′‐TGTTACCAACTGGGACGACA‐3′, β‐actin reverse primer: 5′‐GGGGTGTTGAAGGTCTCAAA‐3′; Nanog forward primer: 5′‐CAAGGGTCTGCTACTGAGATGCTCTG‐3′, Nanog reverse primer: 5′‐TTTTGTTTGGGACTGGTAGAAGAATCAG‐3′; Oct4 (endogenous) forward primer: 5′‐TCTTTCCACCAGGCCCCCGGCTC‐3′, Oct4 (endogenous) reverse primer: 5′‐TGCGGGCGGACATGGGGAGATCC‐3′; Tfe3 forward primer: 5′‐TGCGTCAGCAGCTTATGAGG‐3′, Tfe3 reverse primer: 5′‐AGACACGCCAATCACAGAGAT‐3′; E2A‐c‐Myc forward primer, 5‐GGCTGGAGATGTTGAGAGCAA‐3, E2A‐c‐Myc reverse primer 5‐AAAGGAAATCCAGTGGCGC‐3; Tm9sf2 forward primer: 5′‐GCAACGAGTGCAAGGCTGATA‐3′, Tm9sf2 reverse primer: 5′‐CCCCGAATAATACCTGACCAAGA‐3′; Cd164 forward primer: 5′‐GTGTTTCCTGTGTTAATGCCAC‐3′, Cd164 reverse primer: 5′‐CACAAGTCAGTGCGGTTCAC‐3′; Ddx5 forward primer: 5′‐CGGGATCGAGGGTTTGGTG‐3′, Ddx5 reverse primer: 5′‐GCAGCTCATCAAGATTCCACTTC‐3′; P4hb forward primer: 5′‐GCCGCAAAACTGAAGGCAG‐3′, P4hb reverse primer: 5′‐GGTAGCCACGGACACCATAC‐3′; Lars2 forward primer: 5′‐CATAGAGAGGAATTTGCACCCTG‐3′, Lars2 reverse primer: 5′‐GCCAGTCCTGCTTCATAGAGTTT‐3′; Serinc3 forward primer: 5′‐ GTCCCGTGCCTCTGTAGTG‐3′, Serinc3 reverse primer: 5′‐CAAGACACAATAGTGCCAAGGAA‐3′; Nptn forward primer: 5′‐CGCTGCTCAGAACGAACCAA‐3′, Nptn reverse primer: 5′‐ GCTGGAAGTGAGGTTACACTG‐3′; Cdkn2a forward primer: 5′‐CGAACTCTTTCGGTCGTACCC‐3′, Cdkn2a reverse primer: 5′‐CGAATCTGAACCGTAGTTGAGC‐3′; Cdkn1a forward primer: 5′‐TCTGAGCGGCCTGAAGATTC‐3′, Cdkn1a reverse primer: 5′‐CTGCGCTTGGAGTGATAGAA‐3′; Cxcl2 forward primer: 5′‐CTCAAGGGCGGTCAAAAAGT‐3′, Cxcl2 reverse primer: 5′‐TTTTTCTTTCTCTTTGGTTCTTCC‐3′; Mcp1 forward primer: 5′‐ATTGGGATCATCTTGCTGGT‐3′, Mcp1 reverse primer: 5′‐CCTGCTGTTCACAGTTGCC‐3′; Nrf2 forward primer: 5′‐CTGAACTCCTGGACGGGACTA‐3′, Nrf2 reverse primer: 5′‐CGGTGGGTCTCCGTAAATG‐3′; ApoM forward primer: 5′‐TAACTCCATGAATCAGTGCCCT‐3′, ApoM reverse primer: 5′‐ CCCGCAATAAAGTACCACAGG‐3′.

Bulk RNA sequencing and analysis

RNA and libraries preparation

Total RNA was extracted from pancreas tissues as described above, and treated with the NEBNext rRNA depletion kit (New England Biolabs) to deplete the ribosomal RNA (rRNAs). For the library preparation, we used the NEBNext Ultra II RNA library prep kit for Illumina (New England Biolabs) (9 cycles of amplification). Sequencing was performed on a HiSeq 2500 Sequencing System of Illumina; the type of sequencing was 50 nt Paired‐End and 55 million reads were obtained per sample.

Pre‐processing

Paired‐end reads were aligned to the mm10 genome UCSC Genome Browser using STAR 2.3.0e (Dobin et al., 2013) with the following parameter values: outFilterMismatchNoverLmax = 0.05; outFilterMatchNmin = 25; the rest of parameters were set to their default values. RPKM estimations per isoforms were computed with the R package Casper (Rossell et al., 2014) which were then aggregated at the gene level (entrez). Reads mapping to five or more locations were excluded (parameter keep.multihits set to FALSE). The resulting RPKM matrix was log2‐transformed, quantile normalized, and a priori corrected by the total number of reads in each sample; for doing so, a linear model was fitted to the expression matrix gene‐wise, in which the group condition was included as explanatory variable.

Functional enrichment analysis

Pathway enrichment analysis for group comparisons of gene expression was performed using a modification of ROAST (Wu et al., 2010), a rotation‐based approach implemented in the R package limma (Ritchie et al., 2015) that is especially suitable for small size experiments and is based on limma differential expression. Such modifications were implemented to accommodate the re‐standardized maxmean statistic in the ROAST algorithm (Efron & Tibshirani, 2007), in order to enable it for competitive testing (Goeman & Bühlmann, 2007). For doing so, genes were annotated according to Gene Ontology (GO)(The Gene Ontology Consortium, 2019), Broad Hallmarks (Liberzon et al., 2015) and Kegg (Kanehisa & Goto, 2000) gene‐set collections. GO and Kegg terms were retrieved from R package org. Mm.eg.db (Carlson, 2019), while Broad Hallmark sets were translated to mouse homologous genes using the R package biomaRt (Durinck et al., 2009). For visualization and interpretation purposes, results were represented with Komolgorov–Smirnov based statistic usually used in Gene Set Enrichment Analysis (GSEA) (Subramanian et al., 2005).

Pattern analysis at gene‐set level

To better analyze the behavior of gene‐sets across the three groups of samples, we used a hierarchical model based on a mixture of Normal distributions (Normal‐Normal) (Yuan & Kendziorski, 2006) as implemented in the gaga (Rossell, 2009) R package. This model allows to classify genes according to their expression pattern across the groups of interest. To get a measure of the pathway activity in the transcriptomic data, we summarized each pathway as a gene expression signature. For doing so, expression values were centered and scaled gene‐wise according to the mean and the standard deviation computed across samples, which were then averaged across all genes included in a given gene‐set. In addition, a global signature was computed using all the genes in the expression matrix and used for a priori correction of gene‐set scores by fitting a linear model, in which the group condition was included as explanatory variable. This strategy has proven to be useful to alleviate systematic biases due to the gene‐correlation structure present in the data, and to adjust by the expectation under gene randomization, that is, the association expected for a signature whose genes have been chosen at random (Efron & Tibshirani, 2007; Mestres et al., 2018). Pathway scores were merged with gene level expression, and a gaga (Rossell, 2009) pattern analysis was performed using a Normal‐Normal model fit, according to their observed distribution. In these analyses, the threshold for statistical significance was set at 5% FDR. All statistical analyses were carried out using R and Bioconductor [M19].

Serum analysis of transaminases

Serum was obtained from WT (n = 7) and reprogrammable (n = 6) very old (100w) mice after 6 h of fasting. Samples were collected before doxycycline treatment, as well as after 1 week with doxycycline and 2 weeks of recovery, and they were analyzed by a Spinlab 100 (Spinreact,9‐9059) machine for measuring the levels of transaminases (ALT/GPT and AST/GOT).

Metabolomics

Sample preparation serum

A volume of 25 µl of serum were mixed with 250 µl of a cold solvent mixture with ISTD (MeOH/Water/Chloroform, 9/1/1, −20°C), in a 1.5 ml microtube, before being vortexed and centrifugated (10 min at 15,000 g, 4°C). The upper phase of the supernatant was split into three parts: 50 µl was used for GC‐MS experiment in injection vial, 30 µl was used for the SCFA (Short Chain Fatty Acids) UHPLC‐MS method, and 50 µl was used for other UHPLC‐MS experiments (Viltard et al., 2019).

Widely targeted analysis of intracellular metabolites gas chromatography (GC) coupled to a triple quadrupole (QQQ) mass spectrometer

GC‐MS/MS method was performed on a 7890B gas chromatograph (Agilent Technologies, Waldbronn, Germany) coupled to a triple quadrupole 7000C (Agilent Technologies, Waldbronn, Germany) equipped with a high sensitivity electronic impact source (EI) operating in positive mode (Viltard et al., 2019).

Targeted analysis of bile acids by ion‐pairing ultra‐high performance liquid chromatography (UHPLC) coupled to a Triple Quadrupole (QQQ) mass spectrometer

Targeted analysis was performed on a RRLC 1260 system (Agilent Technologies, Waldbronn, Germany) coupled to a Triple Quadrupole 6410 (Agilent Technologies) equipped with an electrospray source operating in positive mode. Gas temperature was set to 325°C with a gas flow of 12 L/min. Capillary voltage was set to 4.5 kV (Viltard et al., 2019).

Targeted analysis of polyamines by ion‐pairing ultra‐high performance liquid chromatography (UHPLC) coupled to a Triple Quadrupole (QQQ) mass spectrometer

Targeted analysis was performed on a RRLC 1260 system coupled to a Triple Quadrupole 6410 (Agilent Technologies) equipped with an electrospray source operating in positive mode. The gas temperature was set to 350°C with a gas flow of 12 L/min. The capillary voltage was set to 3.5 kV (Viltard et al., 2019).

Targeted analysis of short‐chain fatty acid by ion‐pairing ultra‐high‐performance liquid chromatography (UHPLC) coupled to a 6500+ QTRAP mass spectrometer

Targeted analysis was performed on a RRLC 1260 system coupled to a 6500+ QTRAP (Sciex, Darmstadt, Germany) equipped with an electrospray ion source(Viltard et al., 2019).

Untargeted analysis of intracellular metabolites by ultra‐high performance liquid chromatography (UHPLC) coupled to a Q‐Exactive mass spectrometer. Reversed‐phase acetonitrile method

The profiling experiment was performed with a Dionex Ultimate 3000 UHPLC system (Thermo Scientific) coupled to a Q‐Exactive (Thermo Scientific) equipped with an electrospray source operating in both positive and negative mode and full scan mode from 100 to 1200 m/z. The Q‐Exactive parameters were as follows: sheath gas flow rate 55 au, auxiliary gas flow rate 15 au, spray voltage 3.3 kV, capillary temperature 300°C, S‐Lens RF level 55 V. The mass spectrometer was calibrated with sodium acetate solution dedicated to low mass calibration (Viltard et al., 2019).

Differential expression in metabolomics data

To assess differences of metabolites levels in serum samples, a linear mixed effect was fitted for each metabolite separately in which the biological specimen was included as a random effect. Wald tests derived from the models were used to assess statistical significance. These analyses were performed with R packages lme4 (Bates et al., 2015), lmerTest (Kuznetsova et al., 2017), and multcomp (Hothorn et al., 2008).

Histological analysis

Tissue samples were fixed overnight in 10% neutral‐buffered formalin (4% formaldehyde in solution; Sigma), paraffin‐embedded, sectioned at a thickness of 3 μm, mounted in Superfrost®plus slides and dried. Slides were deparaffinized in xylene and re‐hydrated through a series of graded ethanol until water. Sections were stained with hematoxylin and eosin (HE) for the visualization of the tissue architecture. For immunohistochemistry, sections were stained with a Rabbit mAb Phospho‐Histone H2A.X (Ser139) (20E3) (Cell Signaling, ref: 9718). Slides were dehydrated, cleared, and mounted with toluene‐free mounting medium for microscopic evaluation. Whole digital slides were acquired with a slide scanner and images captured with the NanoZoomer Digital Pathology software (NDP.view2).

Statistical analysis

Mice were randomly allocated to their experimental groups, except from the cohort of very old mice (100 weeks) which were distributed according to their pre‐determined type (mouse genotype), and therefore, there was no randomization. Quantitative PCR data were obtained from independent biological replicates (n values correspond to the number of mice; technical replicates of PCR were not considered in the n value). Statistical analyses were carried out using GraphPad Prism v8.0 (GraphPad software) or as it is indicated in each specific method and stated at the figure legends.

CONFLICT OF INTEREST

D.E.M.H is a founder and shareholder at Chronomics Limited. G.K. is founder of Samsara Therapeutics and advisor of The Longevity Labs. W.R. is consultant and shareholder of Cambridge Epigenetix. M.S. is founder, shareholder, and advisor of Senolytic Therapeutics, Inc., Iduna Therapeutics, Inc, and Rejuveron Senescence Therapeutics, AG. The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.

AUTHOR CONTRIBUTIONS

D.C. designed most experiments, performed animal experimentation, collected samples, performed RNA analysis, contributed to bioinformatic data analysis, prepared the figures, and co‐wrote the manuscript. D.G. and W.R. designed, performed, analyzed, and interpreted the RRBS methylome data. L.M. designed and performed some reprogramming experiments, collected tissue samples, and extracted genomic DNA. R.G.U. and M.F.F. designed, performed, analyzed, and interpreted the bisulfite pyrosequencing experiments. A.B. and C.S.O‐A. designed and contributed to bioinformatic analysis. M. Abad performed the first rejuvenation studies and provided advise throughout the project. D.E.M.H contributed to the original conception and design of the study. M. Aguilera performed all the histological stainings. N.P. supervised and interpreted the pathological analyses, and blindly scored the histological stainings. S.D., F.A., N.N., and G.K. performed, analyzed, and interpreted the metabolomics data. M.S. designed and supervised the study, analyzed the data, and co‐wrote the manuscript. All authors discussed the results and commented on the manuscript.

Fig S1‐S5

Click here for additional data file.

Table S1

Click here for additional data file.

Table S2

Click here for additional data file.

Table S3

Click here for additional data file.

Table S4

Click here for additional data file.

Table S5

Click here for additional data file.

Table S6

Click here for additional data file.

Table S7

Click here for additional data file.

Table S8

Click here for additional data file.

Table S9

Click here for additional data file.

Table S10

Click here for additional data file.