Mouse aging cell atlas analysis reveals global and cell type-specific aging signatures

Abstract
Introduction
Results
Discussion
Materials and methods
Data availability
References
Article and author information
Metrics

Abstract

Aging is associated with complex molecular and cellular processes that are poorly understood. Here we leveraged the Tabula Muris Senis single-cell RNA-seq data set to systematically characterize gene expression changes during aging across diverse cell types in the mouse. We identified aging-dependent genes in 76 tissue-cell types from 23 tissues and characterized both shared and tissue-cell-specific aging behaviors. We found that the aging-related genes shared by multiple tissue-cell types also change their expression congruently in the same direction during aging in most tissue-cell types, suggesting a coordinated global aging behavior at the organismal level. Scoring cells based on these shared aging genes allowed us to contrast the aging status of different tissues and cell types from a transcriptomic perspective. In addition, we identified genes that exhibit age-related expression changes specific to each functional category of tissue-cell types. Altogether, our analyses provide one of the most comprehensive and systematic characterizations of the molecular signatures of aging across diverse tissue-cell types in a mammalian system.

Introduction

Aging leads to the functional decline of major organs across the organism and is the main risk factor for many diseases, including cancer, cardiovascular disease, and neurodegeneration (Niccoli and Partridge, 2012; López-Otín et al., 2013). Past studies have highlighted different hallmarks of the aging process, including genomic instability, telomere attrition, epigenetic alterations, loss of proteostasis, deregulated nutrient sensing, mitochondrial dysfunction, cellular senescence, stem cell exhaustion, and altered intercellular communication (López-Otín et al., 2013; Campisi, 2013; Vijg and Suh, 2013; Nikolich-Žugich, 2018). However, the primary root of aging remains unclear, and the underlying molecular mechanisms are yet to be fully understood.

To gain a better insight into the mammalian aging process at the organismal level, the Tabula Muris Consortium, which we are members of, created the single-cell transcriptomic data set Tabula Muris Senis (TMS) (Tabula Muris Consortium, 2020). TMS is one of the largest expert-curated single-cell RNA sequencing (scRNA-seq) data sets, containing over 300,000 annotated cells from 23 tissues and organs of male and female mice (Mus musculus). The cells were collected from mice of diverse ages, making this data a tremendous opportunity to study the genetic basis of aging across different tissues and cell types. The TMS data is organized into scRNA-seq expression of different tissue-cell type combinations (e.g., B cells in spleen) via expert annotation and clustering.

The original TMS paper explored primarily the cell-centric effects of aging, aiming to characterize changes in cell type composition within different tissues. Here we provide a systematic gene-centric study of gene expression changes occurring during aging across different cell types. The cell-centric and gene-centric perspectives are complementary, as the gene expression can change within the same cell type during aging, even if the cell type composition in the tissue does not vary over time.

Our analysis focused on the TMS FACS data (acquired by cell sorting in microtiter well plates followed by Smart-seq2 library preparation Picelli et al., 2014) because it has more comprehensive coverage of tissues and cell types (Supplementary file 1) and is more sensitive at quantifying gene expression levels as compared to the TMS droplet data. As shown in Figure 1A, the FACS data was collected from 16 C57BL/6JN mice (10 males, 6 females) with ages ranging from 3 months (20-year-old human equivalent) to 24 months (70-year-old human equivalent). It contains 120 cell types from 23 tissues, totaling 110,096 cells. We also used the TMS droplet data (derived from microfluidic droplets) for those tissues where the data is available, to further validate our findings on an additional data set generated by a different method.

Figure 1 with 6 supplements see all

Download asset Open asset

Analysis overview.

(A) Sample description. The TMS FACS data was collected from 16 C57BL/6JN mice (10 males, 6 females) with ages ranging from 3 months (20-year-old human equivalent) to 24 months (70-year-old human equivalent). (B) Significantly aging-dependent genes in all 76 tissue-cell types in the FACS data. The left panels show the number of aging-related genes (discoveries) for each tissue-cell type, broken down into the number of upregulated genes (orange) and the number of downregulated genes (blue), with the numbers on the left showing the ratio (up/down). Tissue-cell types with significantly more up/downregulated genes (ratio >1.5) are highlighted in solid color. Most tissue-cell types have significantly more downregulated aging genes. The right panel shows the number of cells sequenced for each tissue-cell type.

Figure 1—source data 1 Source data for Figure 1.: https://cdn.elifesciences.org/articles/62293/elife-62293-fig1-data1-v2.xlsx
Download elife-62293-fig1-data1-v2.xlsx

We investigated the comprehensive expression signatures of aging across tissues and cell types in the mouse. We performed systematic differential gene expression (DGE) analysis to identify aging-related genes in 76 tissue-cell type combinations across 23 tissues (Figure 1B, Supplementary file 1). Furthermore, we characterized both shared and tissue-cell-specific aging signatures. Our study identified global aging genes (GAGs), namely genes whose expression varies substantially with age in most (>50%) of the tissue-cell types. Interestingly, the changes in expression of these genes are highly concordant across tissue-cell types and exhibit strong bimodality. That is, these genes tend to be either downregulated during aging in most of the tissue-cell types or upregulated across the board. We leveraged this coordinated dynamic to construct an aging score based on the GAGs. We found that the aging score is significantly positively correlated with the chronological age, both in the FACS data and in multiple independent data sets. Moreover, the aging score contrasts the aging status of tissue-cell types with different functionalities and turnover rates, shedding light on the heterogeneous aging process across the 76 tissue-cell types. The score distinguished itself by its single-cell resolution and large data scale, as previous works either studied the biological age at an individual-level (Harries et al., 2011; Holly et al., 2013; Peters et al., 2015; Fleischer et al., 2018; Horvath, 2013; Jylhävä et al., 2017; Petkovich et al., 2017) or focused on a small number of organs (Ori et al., 2015; Arrojo E Drigo et al., 2019). Overall, our analysis highlights the power of scRNA-seq in studying aging and provides a comprehensive catalog of aging-related gene signatures across diverse tissue-cell types.

Results

Identification of aging-related genes

We considered 76 tissue-cell types in the TMS FACS data, 26 tissue-cell types in the TMS droplet data, and 17 tissues in an accompanying bulk RNA-Seq mouse aging study (Schaum et al., 2020) (referred to as the bulk data) with sufficient sample size. We performed DGE analysis for each tissue-cell type separately, treating all cells from the tissue-cell type as samples. We tested if the expression of each gene is significantly related to aging using a linear model treating age as a numerical variable while controlling for sex. We applied an FDR threshold of 0.01 (the number of comparisons corresponds to the number of genes in the tissue-cell type) and an age coefficient threshold of 0.005 (in the unit of log fold change per month, corresponding to 10% fold change from 3 m to 24 m). For details, please refer to the DGE analysis subsection in Materials and methods.

As shown in Figure 1B, the number of significantly age-dependent genes per tissue-cell type ranges from hundreds to thousands. Interestingly, most tissue-cell types have more downregulated aging-related genes than upregulated aging-related genes, suggesting a general decrease in gene expression over aging. This downregulation pattern is unlikely to be confounded by technical factors such as sequencing depth because the 18 m/24 m mice were sequenced at a higher depth (Figure 1—figure supplement 1C–E). By doing separate DGE analyses using mice from one sex or mice from a subset of age groups (3 m/18 m), we further found that such a downregulation pattern was mostly driven by 24 m mice and was not specific to one sex (Figure 1—figure supplement 2). We also observed a similar pattern in the droplet data (Figure 1—figure supplement 3).

In addition, we found that most aging-related genes identified in the analysis have monotonic aging trajectories, meaning that their expressions either increased or decreased monotonically during aging (Figure 1—figure supplements 4,5). However, a subset of genes in the FACS data (13%), while being upregulated during aging overall, increased from 3 m to 18 m and slightly decreased from 18 m to 24 m; those genes are enriched in brown adipose tissue (BAT) B cells, large intestine epithelial cells, and mesenchymal adipose tissue (MAT) mesenchymal stem cells of adipose (Figure 1—figure supplement 4C).

Bimodal effects of aging and global aging genes

We found that most genes are significantly related to aging in at least one tissue-cell type (13,376 in the TMS FACS data and 6233 in the TMS droplet data), consistent with the intuition that aging is a highly complex trait involving many biological processes. The aging-related genes discovered in the FACS data significantly overlap with other important gene sets, including both known human and mouse aging markers as recorded in the GenAge database (Tacutu et al., 2018), senescence genes (Campisi, 2013), transcription factors, eukaryotic initiation factors, and ribosomal protein genes (Figure 2—figure supplement 1). Some of the top overlapping genes, significantly related to aging in most tissue-cell types, include known mouse aging markers Jund, Apoe, and Gpx4 and known human aging markers Jund, Apoe, Fos, and Cdc42 from the GenAge database (Tacutu et al., 2018), and senescence genes Jund, Junb, Ctnnb1, App, and Mapk1. In addition, we found that each tissue-cell type has around 5% aging-related genes that are shared by the GenAge human aging markers. However, we did not find any tissue-cell types that are specifically enriched with these known human aging markers, suggesting that the conservation between mouse aging and human aging is relatively uniform across tissue-cell types (Figure 2—figure supplement 2).

We visualized all aging-related genes (significant in ≥1 tissue-cell type) in Figure 2A, where the color indicates the number of genes. The x-axis shows the weighted proportion of tissue-cell types (out of 76 tissue-cell types) where the gene is significantly related to aging, while the y-axis shows the weighted proportion of tissue-cell types where the gene is upregulated. The tissue-cell type weights used here are inversely proportional to the number of cell types in the tissue, in order to ensure equal representation of the tissues. The visualization makes it clear that there are more downregulated aging genes than upregulated ones, consistent with Figure 1B. Perhaps more strikingly, a bimodal pattern is apparent, in the sense that the aging-related genes tend to have a consistent direction of change during aging across most tissue-cell types. Interestingly, it was also recently reported in other studies that many shared aging-related genes exhibit consistent direction of change during aging across mouse tissues and cell types, including the brain (Ximerakis et al., 2019), kidney, lung, and spleen (Kimmel et al., 2019).

Figure 2 with 4 supplements see all

Download asset Open asset

Tissue-cell level global aging genes (GAGs).

(A) Tissue-cell level aging-related genes with color indicating the number of the genes. The x-axis shows the weighted proportion of tissue-cell types (out of all 76 tissue-cell types) where the gene is significantly related to aging, while the y-axis shows the weighted proportion of tissue-cell types where the gene is upregulated. (**B and C**) Heatmap of age coefficients of the GAGs (panel B) and the top 20 GAGs (panel C). The age coefficients are in the unit of log fold change per month and blue/red represent down/upregulation. (D) Top GO biological pathways for the GAGs.

Figure 2—source data 1 Source data for Figure 2.: https://cdn.elifesciences.org/articles/62293/elife-62293-fig2-data1-v2.xlsx
Download elife-62293-fig2-data1-v2.xlsx

Genes that exhibit age-related expression changes across many cell types have particularly strong bimodal behavior. This motivated us to define global aging genes (GAGs), namely the genes that are significantly related to aging in more than 50% of weighted tissue-cell types (i.e., cell types after normalizing by cell type frequencies using the tissue-cell type weights as described above). We identified 330 GAGs in total, among which 93 are consistently upregulated and 190 consistently downregulated (>80% of weighted tissue-cell types); only 47 have an inconsistent directionality (upregulated in 20–80% of weighted tissue-cell types). We found that the GAGs significantly overlap with genes known in aging-related diseases, including strong overlap with genes related to Alzheimer’s disease (p=2.4e-11), neuroblastoma (p=1.4e-7), fibrosarcoma (p=3.3e-5), and osteoporosis (p=1.5e-4), and relatively weaker overlap with genes related to Huntington’s disease (p=2.5e-3), skin carcinoma (p=3.6e-3), kidney cancer (p=1.2e-3), acute promyelocytic leukemia (p=3.2e-3), acute myeloid leukemia (p=3e-3), endometrial cancer (p=1.6e-3), and hypertension (p=1.9e-3) (please see the Global aging genes subsection in Materials and methods for details). Our results are not sensitive to the specific choice of the 50% threshold for selecting GAGs; using different thresholds produced similar gene ontology (GO) enrichment analysis results (Figure 2—figure supplement 3C) or GAG scores (Figure 3—figure supplement 1A–D) as detailed below.

We visualized the age coefficients for 10 top up/downregulated GAGs (a consistent direction in >80% of weighted tissue-cell types) that are related to aging in the most number of tissue-cell types, as shown in Figure 2C. Many of these genes have been previously shown to be highly relevant to aging. For example, the downregulation of Lars2 has been shown to result in decreased mitochondrial activity and increase the lifespan for C. elegans (Lee et al., 2003). On the other hand, Jund is a proto-oncogene known to protect cells against oxidative stress and its knockout may cause a shortened lifespan in mice (Laurent et al., 2008). Moreover, Rpl13a was observed to be upregulated in almost all tissue-cell types. As a negative regulator of inflammatory proteins, Rpl13a contributes to the resolution phase of the inflammatory response, ensuring that the inflamed tissues are completely restored back to normal tissues. It also contributes to preventing cancerous growth of the injured cells caused by prolonged expression of the inflammatory genes (Zhou et al., 2015; Mazumder et al., 2003). Therefore, it is interesting to observe the upregulation of Rpl13a given that most old mice have severe inflammatory symptoms.

As shown in Figure 2D, Gene Ontology (GO) biological pathway enrichment analysis revealed that the 330 GAGs are associated with apoptosis, translation, biosynthesis, metabolism, and cellular organization. These biological processes are highly relevant to aging (Tower, 2015; Anisimova et al., 2018; Barzilai et al., 2012) and are shared across most cell types, consistent with the intuition that GAGs represent the global aging process across tissue-cell types. In addition, the KEGG pathways associated with the GAGs are consistent with the GO terms and additionally highlighted immune-related pathways and multiple aging-related diseases (Figure 2—figure supplement 4A). Moreover, the findings were supported by similar analyses on the set of 59 GAGs discovered in the droplet data (Figure 2—figure supplement 3A B). We also performed pathway enrichment analysis using the Ingenuity Pathway Analysis (IPA) software (Krämer et al., 2014), which confirmed our findings for the biological processes associated with the GAGs (Figure 2—figure supplement 4B). Of note is the finding that the mTOR pathway, a known aging-associated pathway, is predicted to be inhibited given the expression of the GAGs (Weichhart, 2018; Papadopoli et al., 2019; Johnson et al., 2013; Figure 2—figure supplement 4C). Interestingly, mTOR downregulation has been shown to promote longevity (Stallone et al., 2019; Laplante and Sabatini, 2012), a further indication that the GAGs are related to the aging process.

GAG score contrasts the heterogeneous aging status of tissue-cell types

Following the analysis of GAGs, we next leveraged these marker genes to characterize the holistic aging status of different tissue-cell types. We aggregated the expression of global aging markers into a single score for each cell, referred to as the GAG score (please see the GAG score subsection in Materials and methods for details). We used the FACS data to identify the GAGs because it has more comprehensive coverage of different tissue-cell types. Intuitively, the GAG score tags the global aging process, and it reflects both the chronological age of the organism and the tissue-cell type-specific aging effects. To formally dissect these components, we defined a fixed-effect model with the GAG score being the response variable and various other factors, including the chronological age, sex, and binary-coded tissue-cell types, being explanatory variables. As a sanity check, the chronological age effect on the GAG score is significantly positive (p<1e-100). The model explains 60.2% of the GAG score variance in the TMS FACS data while adding extra interaction terms between age and each tissue-cell type only slightly increased the model fit (explained variance 62.6%), indicating a reasonable fit for the current model (Figure 3—figure supplement 2C). Interestingly, while most young cells have smaller GAG scores and old cells have larger GAG scores, we also found four tissue-cell types whose GAG scores are similar between age groups and eight tissue-cell types that have a subpopulation of cells whose GAG scores are more similar to that of cells from a different age group (Figure 3—figure supplement 3), highlighting the heterogeneity of the aging process across tissue-cell types.

We next considered the tissue-cell type effects on the GAG score. Intuitively, a larger GAG score effect suggests that the corresponding cell type could be molecularly more sensitive to aging compared to other cells in the same animal. As shown in Figure 3A, immune cells and stem cells have higher GAG score effects, while most parenchymal cell types have lower GAG score effects; such a contrast is also statistically significant, as shown in Figure 3B. Indeed, immune cells and stem cells are known to undergo the most substantial changes with aging. Specifically, the aging of the immune system is commonly linked to the impaired capacity of elderly individuals to respond to new infections (Montecino-Rodriguez et al., 2013). Also, adult stem cells are critical for tissue maintenance and regeneration, and the increased incidence of aging-related diseases has been associated with a decline in the stem cell function (Ermolaeva et al., 2018). On the other hand, parenchymal cells like pancreatic cells, neurons, heart myocytes, and hepatocytes have lower aging scores. This could be an indication that these tissue-specialized cell types are more resilient to aging and are able to maintain their functions despite the changes in the animal.

Figure 3 with 8 supplements see all

Download asset Open asset

GAG score.

(A) Tissue-cell GAG score effects with 95% confidence intervals. The color represents the functional category of the tissue-cell type. (**B and C**) Effects of cell functional categories (panel B) and binary cell lifespan (panel C) on the GAG score, meta-analyzed over all tissue-cell types within the group. A positive y-value means that the cells in the group (functional category group for panel B and binary cell lifespan group for panel C) have higher GAG score values than other cells of the same age and sex. 95% confidence intervals and nominal p-values are provided to quantify the differences between categories. In panel C, the average lifespan annotation is also provided for a subset of tissue-cell types where such information is available. (D) Comparison between the tissue-cell GAG score effects estimated from the FACS data (x-axis) and the droplet data (y-axis). Each dot corresponds to a tissue-cell type, and a linear fit is provided showing that the estimates are consistent.

Figure 3—source data 1 Source data for Figure 3.: https://cdn.elifesciences.org/articles/62293/elife-62293-fig3-data1-v2.xlsx
Download elife-62293-fig3-data1-v2.xlsx

We also found that the tissue-cell GAG score effects are in general positively correlated with the cell turnover rate. For example, short-lived cells like skin epidermal cells, monocytes, and T cells (Koster, 2009; Guilliams et al., 2018; Westera et al., 2013) have higher GAG score effects while long-lived cells like neurons, oligodendrocytes, pancreatic β-cells, liver hepatocytes, and heart atrial myocytes (Arrojo E Drigo et al., 2019; Teta et al., 2005; Magami et al., 2002; Bergmann et al., 2015) have very low GAG score effects. To quantify this observation, we assigned a binary cell lifespan label to a subset of cell types where such data is available from literature (Milo et al., 2010; Hobson and Denekamp, 1984; Stewart et al., 1980; Lawson et al., 1992; Arrojo E Drigo et al., 2019; Darwich et al., 2014; Lowry and Zehring, 2017; Geering et al., 2013; Cheshier et al., 1999; Fulcher and Basten, 1997; Scollay et al., 1980) (long for >180 days and short for <90 days, Supplementary file 1); using binary labels instead of the actual values allows us to incorporate more cell types whose exact lifespan information is not available but are known to be long-/short-lived. We found that short-lived cells have significantly higher GAG score effects than long-lived cells (p=1.68e-3, Figure 3C). One possible explanation is that the GAG score is associated with the biological processes related to cell proliferation, development, and death, which are more active in cell types that have a higher turnover rate. This striking difference is also consistent with the intuition that cells that have undergone more divisions (also with higher turnover rates) are ‘‘older’ and could have molecular memories.

Validating the GAG score on external data

We performed several analyses to validate the robustness of the GAG score. First, the GAG score is not sensitive to perturbations of the current scoring method, including using different criterion to select GAGs or not performing cell-wise background correction (Figure 3—figure supplement 1A–D). We also found that estimating tissue-cell GAG score effects using only old cells gave an almost identical result (Figure 3—figure supplement 1H).

Next, we performed a parallel analysis to identify GAGs on the TMS droplet data and found 59 such genes (due to smaller sample size and detection power). We found that 34 genes were shared between the droplet GAGs and the 330 FACS GAGs as described above (p=9e-48). Similar to the FACS GAGs, we also found that the droplet GAGs significantly overlap with genes known in many aging-related diseases, including Alzheimer’s disease (p=2.8e-4), neuroblastoma (p=1.2e-3), and fibrosarcoma (p=1.4e-3). In addition, we considered a set of 261 shared aging genes reported in Kimmel et al., 2019. This scRNA-seq study contains cells from the kidney, lung, and spleen in both young and old mice, and the 261 shared aging genes were defined as genes significantly related to aging in more than five cell types in the paper (Kimmel et al., 2019). We found that 90 genes were shared between Kimmel et al. genes and the FACS GAGs (p=2e-105) and 42 genes were shared between Kimmel et al. genes and the droplet GAGs (p=7e-69). All p-values reported here were computed via Fisher’s exact tests.

Using the GAGs identified from the TMS FACS data, we computed the GAG score and further estimated the GAG score effects for cells in the TMS droplet data, the bulk data (treating each mouse sample as a ‘‘cell’’), the Kimmel et al. data, and the data set from Kowalczyk et al., 2015. This last data set has only three subtypes of hematopoietic stem cells (HSCs) and was therefore omitted in other analyses. We found that the chronological age effect on the GAG score is significantly positive in all four validation data sets (p=7e-10 for the bulk data due to smaller sample size and detection power, and p<1e-100 for the other three data sets). Since the GAGs were selected based on the FACS data, the GAG score is agnostic of the age labels in the four validation data sets, confirming that the GAG score is truly indicative of the aging process.

The tissue-cell GAG score effects estimated from the other four data sets are also in line with those estimated from the FACS data (Figure 3—figure supplements 1,4E F). Specifically, in the droplet data and the Kimmel et al. data, immune cell types have higher GAG score effects while epithelial, endothelial, and parenchymal cell types have lower GAG score effects. In particular for the droplet data, short-lived cell types also have higher but non-significant GAG score effects, due to a smaller number of annotated cell types. While looking at the bulk data, we found that immune-related tissues and organs such as whole blood, spleen, and marrow have the highest GAG score effects. It is interesting to observe that in the Kowalczyk et al. data, the GAG score effects of MPPs (multipotent progenitors) is less than that of ST-HSCs (short-term HSCs) which is less than that of LT-HSCs (long-term HSCs). This is exactly aligned with the differentiation potentials of these three cell types, consistent with the hypothesis that more stem-like cells have higher GAG score effects as observed in the FACS data.

Finally, we found that the tissue-cell GAG score effects are highly consistent between data sets (correlation 0.76 with p=2e-19 between the FACS data and the droplet data in Figure 3D, correlation 0.75 with p=1e-3 between the FACS data and the Kimmel et al. data in Figure 3—figure supplement 1G). In summary, we showed that the GAG score is capable of describing the chronological age as well as the transcriptional changes during the aging process. Furthermore, we could use the tissue-cell GAG score effects to contrast the aging status of cell types with different biological properties, including functional categories and turnover rates. This provides a comprehensive analysis demonstrating how the GAG score captures the heterogeneous molecular effects of aging. Additionally, we repeated both the DGE analysis and the GAG score analysis at the tissue-level by combining all cell types from the same tissue. We observed qualitatively similar findings, supporting the robustness of the analyses. Please see Figure 3—figure supplements 5–8 for more details.

Category-specific aging genes

We next consider genes specific to a subset of tissue-cell types, including functional-category-specific genes, cell type-specific genes, tissue-specific-genes, and tissue-cell type-specific genes. Given a set of tissue-cell types, in order to have an overall meta age coefficient for cells in this tissue-cell type set, we first combined the age coefficients of all tissue-cel types within the set by meta-analysis; similarly, we also computed the outside-set meta age coefficient by meta-analyzing all outside-set tissue-cell types. Then we selected the genes that have significantly different within-set and outside-set meta age coefficients as the set-specific genes (please see the Category-specific aging genes subsection in Materials and methods for more details). Of note, almost no genes identified here are shared by GAGs.

In the original TMS paper (Tabula Muris Consortium, 2020), each tissue-cell type was assigned one of the six functional category labels, namely endothelial, epithelial, immune, stem/progenitor, stromal, and parenchymal cells. When examining the data by functional category, we found that the endothelial, immune, stem, and stromal cells exhibit highly category-specific aging behavior. Indeed, we found a higher number of specific aging genes for these categories (Figure 4B). Moreover, their age coefficients are specific to the respective functional category, as we can see from the clear block structure across tissue-cell types (Figure 4A). In addition, when performing GO biological pathway enrichment analysis on these six sets of genes separately, we only found significant pathways for these four categories (Figure 4C). Among them, endothelial-specific genes were associated with various processes related to angiogenesis and negative regulation of cell migration; the latter suggests decreased endothelial cell functionality during aging because endothelial cell migration is essential to angiogenesis (Lamalice et al., 2007). Also, immune-specific genes were associated with activation of various immune responses, in line with a strong link between the aging process and the immune system (Tabula Muris Consortium, 2020; Nikolich-Žugich, 2018). In addition, stem-specific genes were associated with ossification and diverse angiogenesis processes, and both stem-specific and stromal-specific genes were associated with extracellular matrix and structure organization.

Figure 4 with 3 supplements see all

Download asset Open asset

Functional-category-specific genes.

(A) Age coefficients, in the unit of log fold change per month, for functional-category-specific genes. Both genes in the x-axis and tissue-cell types in the y-axis are ordered by functional categories. For example, the upper-left block corresponds to endothelial-specific genes. (B) Number of functional-category-specific genes for each category. (C) GO biological pathways for functional-category-specific genes, with color representing the negative log10 FDR.

Figure 4—source data 1 Source data for Figure 4.: https://cdn.elifesciences.org/articles/62293/elife-62293-fig4-data1-v2.xlsx
Download elife-62293-fig4-data1-v2.xlsx

Such an analysis also facilitates the discovery of interesting specific genes related to aging. For example, C2cd4b (Figure 4—figure supplements 1 and 2A), a parenchymal-specific gene, has large age coefficients in several pancreatic, mammary gland, large intestine cell types, and almost zero age coefficients in other cell types. The increased expression of C2cd4b has been associated with an increased risk of type 2 diabetes (Kycia et al., 2018) and increased expression of C2cd4b in old mice pancreatic cell types may suggest the increased risk of type 2 diabetes for these mice. In addition, C2cd4b has been shown to lead to sexually dimorphic changes in body weight and glucose homeostasis (Mousavy Gharavy et al., 2021), in line with the fact that the mammary gland is a well-known sexually dimorphic tissue. A second example is Gsn (Figure 4—figure supplements 1 and 2B), which is downregulated in stromal and stem cell types and upregulated in other cell types during aging. This can be explained by its function for making gelsolin, an important protein for cell movement. Not only is cell movement important to immune and endothelial cells, but Gsn has also been shown to be a potential biomarker to aging-related neurodegeneration (Manavalan et al., 2013).

Beyond functional-category-specific genes, we also identified genes specific to several cell types, including B cells, basal cells of epidermis, endothelial cells, macrophages, mesenchymal stem cells, mesenchymal stem cells of adipose, myeloid cells, and skeletal muscle satellite cells, corroborated by their association with related biological processes (Figure 4—figure supplement 3). The method also allowed us to identify genes specific to each tissue. However, we did not find any genes specific to a single tissue-cell type. All the gene sets are available in Supplementary file 3.

Discussion

This study provides a systematic and comprehensive analysis of aging-related transcriptomic signatures by analyzing 76 tissue-cell types in the TMS FACS data. Together with the analysis in the first publication of Tabula Muris Senis (Tabula Muris Consortium, 2020), it forms one of the largest analysis to date of the mammalian aging process at the single-cell resolution. Of particular interest are the 330 global aging genes (GAGs) identified in the study. These genes exhibit aging-dependent expressions in a majority of tissue-cell types in the mouse. The GAGs are enriched with many interesting genes, including known human and mouse aging markers, aging-related disease genes, senescence genes, transcription factors, eukaryotic initiation factors, and ribosomal protein genes. Interestingly, most of the GAGs are strongly bimodal as their expressions either decrease or increase during aging across almost all tissue-cell types, suggesting that these genes have a uniform response to aging, which is robust to the specific tissue or cellular context. Moreover, we find a systematic decrease in expression for most genes as well as a decrease in the number of actively expressed genes, suggesting a turning off of transcription activity as the animal ages. A recent study has observed that the number of expressed genes decreases during cellular differentiation in mouse (Gulati et al., 2020). It is interesting that we quantify a similar phenomenon for aging, despite the substantial longer time-scale of aging compared to differentiation.

While we focused on detecting genes whose changes with age are linear and non-sex-specific, it is also interesting to study aging-related changes that are non-linear or sex-specific. For example, as shown in Figure 1—figure supplement 1, the number of expressed genes changes with age in a non-linear and sex-specific manner, suggesting the existence of such genes. We have validated our findings using the TMS droplet data, the bulk RNA-seq data, and external data sets, and it is important to have further validations in future studies. In particular, the bimodal expression pattern is less apparent in the TMS droplet data, perhaps due to its limited tissue-cell type coverage and relatively shallower sequencing depth.

The remarkable bimodal consistency of the GAGs makes them useful as biomarkers to characterize the aging status of individual cells. We proposed a new aging score, namely the GAG score, based on the GAGs. The tissue-cell type-specific GAG score effects quantify how sensitive each tissue-cell type is to aging and are positively correlated with the cell division rate. For example, immune cells tend to have higher GAG scores than other cells of the same age and sex, which reflects the phenomenon that they undergo many cycles of cell division and also change substantially during the animal’s lifespan. One hypothesis is that the GAG score captures some aspects of the true biological age of the cells, which could be different from the birth age of the animal. An interesting direction of future work is to further investigate this model with functional experiments. In line with this, it would be important to study how some of the transcriptomic changes we quantify here, for example the downregulation of mTOR, point toward healthy aging or how can they inform experiments that can uncover the mechanism to ameliorate the aging effects. In addition, while the GAG score was proposed to capture the global aging status, there also exist cell type-specific aging programs; a combination of both would better characterize the overall aging status of a cell. The construction of cell type-specific aging scores would require a more comprehensive longitudinal catalog of diverse subtypes and states of cells within different cell types. While this is beyond the scope of the current work, it is an exciting direction to pursue when such data becomes available.

The GAG score is also related to the transcriptome age predictors developed in previous works (Harries et al., 2011; Holly et al., 2013; Peters et al., 2015; Fleischer et al., 2018), in the sense that they all use the gene expression information and are predictive of the animals’ chronological age. The commonly used approach in previous works is to train a model (e.g., linear/logistic regression model) to predict the individuals’ chronological age from their gene expression. Instead of a model-fitting algorithm, our GAG score uses the GAGs that were selected from a broad range of tissue-cell types in an unbiased manner, by meta-analyzing the DGE results of 76 tissue-cell types and putting each tissue-cell type on the same footing. This ensures that the genes used by the GAG score capture the shared aging process and are not biased toward certain tissue-cell types. Indeed, the GAGs were shown to be associated with biological processes that are highly relevant to aging (Figure 2), providing better interpretability of the score. In comparison, previous studies focused on only one specific tissue, such as the blood (Harries et al., 2011; Holly et al., 2013; Peters et al., 2015) or dermal fibroblasts (Fleischer et al., 2018), and hence may have selected genes that were biased toward that particular tissue.

Overall, our study provides a comprehensive characterization of aging genes across a wide range of tissue-cell types in mice. In addition to the biological insights, it also serves as a comprehensive reference for researchers working on related topics.

Materials and methods

Data preprocessing

Request a detailed protocol

We considered five data sets, namely the TMS FACS data, the TMS droplet data, the data in Schaum et al., 2020 (referred to as the bulk data), the data in Kimmel et al., 2019, and the data in Kowalczyk et al., 2015. For the TMS FACS data and the TMS droplet data, we filtered out genes expressed in fewer than 3 cells, filtered out cells expressing fewer than 250 genes, and discarded cells with a total number of counts fewer than 5000 for the FACS data and a total number of unique molecular identifiers (UMIs) fewer than 2500 for the droplet data. For the bulk data, we filtered out genes expressed in fewer than five samples, and filtered out samples expressing fewer than 500 genes. We did not filter cells for the other two data sets. For all five data sets, we normalized each sample to have 10,000 reads/UMIs per sample, followed by a log transformation (log(x + 1) where x is the read count). We note that such a procedure is the same as that in the original paper (Tabula Muris Consortium, 2020). We did not correct for batch effects as no substantial batch effects were identified in the original TMS paper (Tabula Muris Consortium, 2020).

DGE analysis

View detailed protocol

As shown in Supplementary file 1, we considered 76 tissue-cell types in 23 tissues with more than 100 cells in both young (3 m) and old (18 m, 24 m) age groups for the TMS FACS data; 26 tissue-cell types in 11 tissues with more than 500 cells in both young (1 m, 3 m) and old (18 m, 21 m, 24 m, 30 m) age groups for the TMS droplet data; and all 17 tissues for the bulk data (Schaum et al., 2020). We required more cells for the TMS droplet data than the TMS FACS data because the droplet data has a much lower sequencing depth (6000 UMIs per cell, as compared to 0.85 million reads per cell for the FACS data). Also, we did not focus on the TMS droplet data in the main results due to its limited tissue and cell type coverage.

We performed a DGE analysis for cells in each tissue-cell type separately. In the DGE analysis for a tissue-cell type, all cells in the tissue-cell type were treated as samples, and a separate test was performed for each gene with the observations being the expressions of the gene across the cells. We identified genes significantly related to aging using a linear model treating age as a continuous variable while controlling for sex, namely,

g e n e e x p r e s s i o n \sim a g e + s e x .

Since the tests were performed at a cell level for each tissue-cell type, the cell numbers will only affect the detection power but will not bias the result. Therefore, it is not a confounding factor and was not controlled for. We used the MAST package (Finak et al., 2015) (version 1.12.0) in R to perform the DGE analysis. The zero counts in the scRNA-seq data were handled by the MAST package and we do not observe other types of missing data. We did not control for the cellular detection rate (Finak et al., 2015) (CDR, corresponding to the number of expressed genes in a cell) because we found that in our data, CDR is positively correlated with age and negatively correlated with technical covariates such as sequencing depth and the number of detected ERCC spike-ins, both when considering all the cells or when stratified by sex (Figure 1—figure supplement 1). As a result, controlling for CDR may remove genuine aging effects. We note that CDR is defined as the number of expressed genes in a cell, which is a fundamental quantity of the data set and is not specific to the MAST package. Therefore, such an observation does not imply a potential issue of MAST. Nonetheless, we found that the age coefficients, estimated with and without CDR correction, were highly correlated (0.89 for the FACS data and 0.93 for the droplet data, Figure 1—figure supplement 6), ruling out the possibility that CDR correction would significantly alter the result.

We used the Benjamini–Hochberg (FDR) procedure (Benjamini and Hochberg, 1995) to control for multiple comparisons, where the number of comparisons corresponds to the number of genes in the tissue-cell type. We applied an FDR threshold of 0.01 and an age coefficient threshold of 0.005 (in the unit of log fold change per month, corresponding to around 10% fold change from 3 m to 24 m) for detecting genes significantly related to aging.

Global aging genes

Request a detailed protocol

We selected a gene as a GAG if it is significantly related to aging in more than 50% of weighted tissue-cell types. Here, the tissue-cell type weights are inversely proportional to the number of cell types in the tissue, in order to ensure equal representation of tissues.

For the overlap between GAGs and the genes known in aging-related diseases, we considered the top 25 aging-related diseases and obtained their related genes from the Human Disease Database (MalaCards) (Rappaport et al., 2013; Rappaport et al., 2014; Rappaport et al., 2017). We then converted the human genes to the corresponding mouse orthologs using g:Profiler (Raudvere et al., 2019) (version 1.2.2). The p-values quantifying the significance of the overlap were computed using Fisher’s exact tests.

Pathway enrichment analysis

Request a detailed protocol

We used g:Profiler (Raudvere et al., 2019) to perform GO biological pathway enrichment analysis. We considered biological pathways with FDR smaller than 0.01. We used Gene Set Enrichment Analysis (GSEA MGSig Database) (Subramanian et al., 2005; Liberzon et al., 2011) to perform the KEGG pathway analysis. We filtered for mouse genes and considered biological pathways with FDR smaller than 0.05. We also used the IPA software (Krämer et al., 2014) to perform canonical pathway analysis (Figure 2—figure supplement 4). For Figure 2—figure supplement 4A, we used an FDR threshold of 1e-5 and a z-score threshold of 0.5.

GAG score

Request a detailed protocol

Given a set of GAGs (e.g., the FACS GAGs), the cell-wise GAG score for cell $i$ is computed as:

Compute the raw GAG score as the average expression of the upregulated GAGs (significantly related to aging in >50% of weighted tissue-cell types and upregulated in >80% of weighted tissue-cell types) minus the average expression of the downregulated GAGs (significantly related to aging in >50% of weighted tissue-cell types and downregulated in >80% of weighted tissue-cell types), That is,
${R a w G A G s c o r e}_{i} = m e a n (u p r e g u l a t e d G A G s)_{i} - m e a n (d o w n r e g u l a t e d G A G s)_{i}$
Following the recipe of DeTomaso et al., 2019, z-normalize the raw GAG score using the expected mean and variance of a random set of genes with the same number of up/downregulated genes:
${G A G s c o r e}_{i} = {R a w G A G s c o r e}_{i} / [{s t d}_{i} * \sqrt{1 / n_{u p} + 1 / n_{d o w n}}],$

where ${s t d}_{i}$ is the standard deviation of the gene expression of cell $i$ , $n_{up}$ is the number of upregulated GAGs, and $n_{down}$ is the number of downregulated GAGs.

For robustness consideration, we checked the expression levels of the genes used for computing the GAG score and found that there are no extreme values. Therefore, the GAG score is unlikely to be dominated by a few highly expressed genes. We also considered another version of GAG score where each gene is weighted by its expression range. We found the two versions produce highly correlated results. Please see Figure 3—figure supplement 2A B for more details.

For estimating the GAG score effects, we model the GAG score of a cell $i$ as being linearly dependent of the chronological age, sex, and tissue-cell type of the cell, namely,

{G A G s c o r e}_{i} \sim {a g e}_{i} + {s e x}_{i} + \sum_{j} {tissue-cell type}_{i j}

Here, ${a g e}_{i}$ is the age of the animal (in months) that cell $i$ comes from, ${s e x}_{i}$ is one if cell $i$ comes from a male and zero otherwise, and ${tissue-cell type}_{i j}$ is one if cell $i$ belongs to tissue-cell type $j$ and zero otherwise. We do not include the intercept term because all binary coded tissue-cell types sum up to one. Finally, we further center both response and explanatory variables and perform an ordinary least square regression to estimate the GAG score effects for the age, sex, and each tissue-cell type.

Meta-analyses in Figure 3B and C and Figure 3—figure supplement 1E F were performed assuming a random effect model (Riley et al., 2011). For comparisons of the tissue-cell GAG score effects between data sets, namely those in Figure 3D and Figure 3—figure supplement 1G, we used all tissue-cell types instead of restricting to the 76 TMS FACS tissue-cell types and 26 TMS droplet tissue-cell types. This increased the number of overlapping tissue-cell types between data sets.

Category-specific aging genes

Request a detailed protocol

We considered identifying functional-category-specific genes, cell type-specific genes, tissue-specific genes, and tissue-cell type-specific genes. For a set of tissue-cell types, for example the set of all immune tissue-cell types, the genes specific to the set (or set-specific genes) are selected as follows. For each gene, we first estimate its within-set meta age coefficient by meta-analyzing the age coefficients of all tissue-cell types in the set assuming a random effect model (Riley et al., 2011). Specifically, it is done by assuming that there is a meta age coefficient for the set of tissue-types, and the age coefficient for each tissue-cell type in the set is a random variable whose mean is equal to the meta age coefficient of the set. Similarly, we estimate the outside-set meta age coefficient by meta-analyzing all tissue-cell types outside the set. Then, we define the set-specific genes to be genes whose:

within-set meta age coefficient is significantly different from its outside-set meta age coefficient (FDR < 0.01);
within-set meta age coefficient is large enough (absolute value > 0.005);
outside-set meta age coefficient is not significantly different from 0 (FDR > 0.01).

The p-values are computed based on the mean and the standard error assuming a normal distribution, and FDR is computed with respect to all genes.

Code availability

Request a detailed protocol

The code for reproducing all results is at https://github.com/czbiohub/tabula-muris-senis/tree/master/2_aging_signature (Zhang and Pisco, 2021; copy archived at swh:1:rev:0fd2ee501f4f3bd0e691b6071aee2c9286f1cf92).

Data availability

All data can be downloaded at https://figshare.com/articles/dataset/tms_gene_data_rv1/12827615.

The following previously published data sets were used

1. The Tabula Muris Consortium
(2020) figshare
Tabula Muris.

https://doi.org/10.6084/m9.figshare.12827615.v1

References

(2018) Protein synthesis and quality control in aging
Aging 10:4269–4288.

https://doi.org/10.18632/aging.101721
- PubMed
- Google Scholar
1. Arrojo E Drigo R
2. Lev-Ram V
3. Tyagi S
4. Ramachandra R
5. Deerinck T
6. Bushong E
7. Phan S
8. Orphan V
9. Lechene C
10. Ellisman MH
11. Hetzer MW
(2019) Age mosaicism across multiple scales in adult tissues
Cell Metabolism 30:343–351.

https://doi.org/10.1016/j.cmet.2019.05.010
- PubMed
- Google Scholar
(2012) The critical role of metabolic pathways in aging
Diabetes 61:1315–1322.

https://doi.org/10.2337/db11-1300
- PubMed
- Google Scholar
1. Benjamini Y
2. Hochberg Y
(1995) Controlling the false discovery rate: a practical and powerful approach to multiple testing
Journal of the Royal Statistical Society: Series B 57:289–300.

https://doi.org/10.1111/j.2517-6161.1995.tb02031.x
- Google Scholar
1. Bergmann O
2. Zdunek S
3. Felker A
4. Salehpour M
5. Alkass K
6. Bernard S
7. Sjostrom SL
8. Szewczykowska M
9. Jackowska T
10. Dos Remedios C
11. Malm T
12. Andrä M
13. Jashari R
14. Nyengaard JR
15. Possnert G
16. Jovinge S
17. Druid H
18. Frisén J
(2015) Dynamics of cell generation and turnover in the human heart
Cell 161:1566–1575.

https://doi.org/10.1016/j.cell.2015.05.026
- PubMed
- Google Scholar
1. Campisi J
(2013) Aging, cellular senescence, and Cancer
Annual Review of Physiology 75:685–705.

https://doi.org/10.1146/annurev-physiol-030212-183653
- PubMed
- Google Scholar
(1999) Vivo proliferation and cell cycle kinetics of long-term self-renewing hematopoietic stem cells
PNAS 96:3120–3125.

https://doi.org/10.1073/pnas.96.6.3120
- PubMed
- Google Scholar
(2014) Meta-analysis of the turnover of intestinal epithelia in preclinical animal species and humans
Drug Metabolism and Disposition 42:2016–2022.

https://doi.org/10.1124/dmd.114.058404
- PubMed
- Google Scholar
1. DeTomaso D
2. Jones MG
3. Subramaniam M
4. Ashuach T
5. Ye CJ
6. Yosef N
(2019) Functional interpretation of single cell similarity maps
Nature Communications 10:1–11.

https://doi.org/10.1038/s41467-019-12235-0
- PubMed
- Google Scholar
1. Ermolaeva M
2. Neri F
3. Ori A
4. Rudolph KL
(2018) Cellular and epigenetic drivers of stem cell ageing
Nature Reviews Molecular Cell Biology 19:594–610.

https://doi.org/10.1038/s41580-018-0020-3
- PubMed
- Google Scholar
1. Finak G
2. McDavid A
3. Yajima M
4. Deng J
5. Gersuk V
6. Shalek AK
7. Slichter CK
8. Miller HW
9. McElrath MJ
10. Prlic M
11. Linsley PS
12. Gottardo R
(2015) MAST: a flexible statistical framework for assessing transcriptional changes and characterizing heterogeneity in single-cell RNA sequencing data
Genome Biology 16:278.

https://doi.org/10.1186/s13059-015-0844-5
- PubMed
- Google Scholar
1. Fleischer JG
2. Schulte R
3. Tsai HH
4. Tyagi S
5. Ibarra A
6. Shokhirev MN
7. Huang L
8. Hetzer MW
9. Navlakha S
(2018) Predicting age from the transcriptome of human dermal fibroblasts
Genome Biology 19:221.

https://doi.org/10.1186/s13059-018-1599-6
- PubMed
- Google Scholar
1. Fulcher DA
2. Basten A
(1997) B cell life span: a review
Immunology and Cell Biology 75:446–455.

https://doi.org/10.1038/icb.1997.69
- PubMed
- Google Scholar
1. Geering B
2. Stoeckle C
3. Conus S
4. Simon HU
(2013) Living and dying for inflammation: neutrophils, eosinophils, basophils
Trends in Immunology 34:398–409.

https://doi.org/10.1016/j.it.2013.04.002
- PubMed
- Google Scholar
(2018) Developmental and functional heterogeneity of monocytes
Immunity 49:595–613.

https://doi.org/10.1016/j.immuni.2018.10.005
- PubMed
- Google Scholar
1. Gulati GS
2. Sikandar SS
3. Wesche DJ
4. Manjunath A
5. Bharadwaj A
6. Berger MJ
7. Ilagan F
8. Kuo AH
9. Hsieh RW
10. Cai S
11. Zabala M
12. Scheeren FA
13. Lobo NA
14. Qian D
15. Yu FB
16. Dirbas FM
17. Clarke MF
18. Newman AM
(2020) Single-cell transcriptional diversity is a hallmark of developmental potential
Science 367:405–411.

https://doi.org/10.1126/science.aax0249
- PubMed
- Google Scholar
1. Harries LW
2. Hernandez D
3. Henley W
4. Wood AR
5. Holly AC
6. Bradley-Smith RM
7. Yaghootkar H
8. Dutta A
9. Murray A
10. Frayling TM
11. Guralnik JM
12. Bandinelli S
13. Singleton A
14. Ferrucci L
15. Melzer D
(2011) Human aging is characterized by focused changes in gene expression and deregulation of alternative splicing
Aging Cell 10:868–878.

https://doi.org/10.1111/j.1474-9726.2011.00726.x
- PubMed
- Google Scholar
1. Hobson B
2. Denekamp J
(1984) Endothelial proliferation in tumours and normal tissues: continuous labelling studies
British Journal of Cancer 49:405–413.

https://doi.org/10.1038/bjc.1984.66
- PubMed
- Google Scholar
1. Holly AC
2. Melzer D
3. Pilling LC
4. Henley W
5. Hernandez DG
6. Singleton AB
7. Bandinelli S
8. Guralnik JM
9. Ferrucci L
10. Harries LW
(2013) Towards a gene expression biomarker set for human biological age
Aging Cell 12:324–326.

https://doi.org/10.1111/acel.12044
- PubMed
- Google Scholar
1. Horvath S
(2013) DNA methylation age of human tissues and cell types
Genome Biology 14:R115.

https://doi.org/10.1186/gb-2013-14-10-r115
- Google Scholar
(2013) mTOR is a key modulator of ageing and age-related disease
Nature 493:338–345.

https://doi.org/10.1038/nature11861
- PubMed
- Google Scholar
(2017) Biological age predictors
EBioMedicine 21:29–36.

https://doi.org/10.1016/j.ebiom.2017.03.046
- PubMed
- Google Scholar
(2019) Murine single-cell RNA-seq reveals cell-identity- and tissue-specific trajectories of aging
Genome Research 29:2088–2103.

https://doi.org/10.1101/gr.253880.119
- PubMed
- Google Scholar
1. Koster MI
(2009) Making an epidermis
Annals of the New York Academy of Sciences 1170:7–10.

https://doi.org/10.1111/j.1749-6632.2009.04363.x
- PubMed
- Google Scholar
1. Kowalczyk MS
2. Tirosh I
3. Heckl D
4. Rao TN
5. Dixit A
6. Haas BJ
7. Schneider RK
8. Wagers AJ
9. Ebert BL
10. Regev A
(2015) Single-cell RNA-seq reveals changes in cell cycle and differentiation programs upon aging of hematopoietic stem cells
Genome Research 25:1860–1872.

https://doi.org/10.1101/gr.192237.115
- PubMed
- Google Scholar
(2014) Causal analysis approaches in ingenuity pathway analysis
Bioinformatics 30:523–530.

https://doi.org/10.1093/bioinformatics/btt703
- PubMed
- Google Scholar
1. Kycia I
2. Wolford BN
3. Huyghe JR
4. Fuchsberger C
5. Vadlamudi S
6. Kursawe R
7. Welch RP
8. Albanus RD
9. Uyar A
10. Khetan S
11. Lawlor N
12. Bolisetty M
13. Mathur A
14. Kuusisto J
15. Laakso M
16. Ucar D
17. Mohlke KL
18. Boehnke M
19. Collins FS
20. Parker SCJ
21. Stitzel ML
(2018) A common type 2 diabetes risk variant potentiates activity of an evolutionarily conserved islet stretch enhancer and increases C2CD4A and C2CD4B expression
The American Journal of Human Genetics 102:620–635.

https://doi.org/10.1016/j.ajhg.2018.02.020
- PubMed
- Google Scholar
(2007) Endothelial cell migration during angiogenesis
Circulation Research 100:782–794.

https://doi.org/10.1161/01.RES.0000259593.07661.1e
- PubMed
- Google Scholar
1. Laplante M
2. Sabatini DM
(2012) mTOR signaling in growth control and disease
Cell 149:274–293.

https://doi.org/10.1016/j.cell.2012.03.017
- PubMed
- Google Scholar
1. Laurent G
2. Solari F
3. Mateescu B
4. Karaca M
5. Castel J
6. Bourachot B
7. Magnan C
8. Billaud M
9. Mechta-Grigoriou F
(2008) Oxidative stress contributes to aging by enhancing pancreatic angiogenesis and insulin signaling
Cell Metabolism 7:113–124.

https://doi.org/10.1016/j.cmet.2007.12.010
- PubMed
- Google Scholar
(1992) Turnover of resident microglia in the normal adult mouse brain
Neuroscience 48:405–415.

https://doi.org/10.1016/0306-4522(92)90500-2
- PubMed
- Google Scholar
1. Lee SS
2. Lee RY
3. Fraser AG
4. Kamath RS
5. Ahringer J
6. Ruvkun G
(2003) A systematic RNAi screen identifies a critical role for mitochondria in C. elegans longevity
Nature Genetics 33:40–48.

https://doi.org/10.1038/ng1056
- PubMed
- Google Scholar
(2011) Molecular signatures database (MSigDB) 3.0
Bioinformatics 27:1739–1740.

https://doi.org/10.1093/bioinformatics/btr260
- PubMed
- Google Scholar
(2013) The hallmarks of aging
Cell 153:1194–1217.

https://doi.org/10.1016/j.cell.2013.05.039
- PubMed
- Google Scholar
1. Lowry LE
2. Zehring WA
(2017) Potentiation of natural killer cells for cancer immunotherapy: a review of literature
Frontiers in Immunology 8:1061.

https://doi.org/10.3389/fimmu.2017.01061
- PubMed
- Google Scholar
1. Magami Y
2. Azuma T
3. Inokuchi H
4. Kokuno S
5. Moriyasu F
6. Kawai K
7. Hattori T
(2002) Cell proliferation and renewal of normal hepatocytes and bile duct cells in adult mouse liver
Liver 22:419–425.

https://doi.org/10.1034/j.1600-0676.2002.01702.x
- PubMed
- Google Scholar
1. Manavalan A
2. Mishra M
3. Feng L
4. Sze SK
5. Akatsu H
6. Heese K
(2013) Brain site-specific proteome changes in aging-related dementia
Experimental & Molecular Medicine 45:e39.

https://doi.org/10.1038/emm.2013.76
- PubMed
- Google Scholar
1. Mazumder B
2. Sampath P
3. Seshadri V
4. Maitra RK
5. DiCorleto PE
6. Fox PL
(2003) Regulated release of L13a from the 60S ribosomal subunit as a mechanism of transcript-specific translational control
Cell 115:187–198.

https://doi.org/10.1016/S0092-8674(03)00773-6
- PubMed
- Google Scholar
1. Milo R
2. Jorgensen P
3. Moran U
4. Weber G
5. Springer M
(2010) BioNumbers--the database of key numbers in molecular and cell biology
Nucleic Acids Research 38:D750–D753.

https://doi.org/10.1093/nar/gkp889
- PubMed
- Google Scholar
(2013) Causes, consequences, and reversal of immune system aging
Journal of Clinical Investigation 123:958–965.

https://doi.org/10.1172/JCI64096
- PubMed
- Google Scholar
1. Mousavy Gharavy SN
2. Owen BM
3. Millership SJ
4. Chabosseau P
5. Pizza G
6. Martinez-Sanchez A
7. Tasoez E
8. Georgiadou E
9. Hu M
10. Fine NHF
11. Jacobson DA
12. Dickerson MT
13. Idevall-Hagren O
14. Montoya A
15. Kramer H
16. Mehta Z
17. Withers DJ
18. Ninov N
19. Gadue PJ
20. Cardenas-Diaz FL
21. Cruciani-Guglielmacci C
22. Magnan C
23. Ibberson M
24. Leclerc I
25. Voz M
26. Rutter GA
(2021) Sexually dimorphic roles for the type 2 diabetes-associated C2cd4b gene in murine glucose homeostasis
Diabetologia 64:850–864.

https://doi.org/10.1007/s00125-020-05350-x
- PubMed
- Google Scholar
1. Niccoli T
2. Partridge L
(2012) Ageing as a risk factor for disease
Current Biology 22:R741–R752.

https://doi.org/10.1016/j.cub.2012.07.024
- PubMed
- Google Scholar
1. Nikolich-Žugich J
(2018) The twilight of immunity: emerging concepts in aging of the immune system
Nature Immunology 19:10–19.

https://doi.org/10.1038/s41590-017-0006-x
- PubMed
- Google Scholar
1. Ori A
2. Toyama BH
3. Harris MS
4. Bock T
5. Iskar M
6. Bork P
7. Ingolia NT
8. Hetzer MW
9. Beck M
(2015) Integrated transcriptome and proteome analyses reveal Organ-Specific proteome deterioration in old rats
Cell Systems 1:224–237.

https://doi.org/10.1016/j.cels.2015.08.012
- PubMed
- Google Scholar
1. Papadopoli D
2. Boulay K
3. Kazak L
4. Pollak M
5. Mallette F
6. Topisirovic I
7. Hulea L
(2019) mTOR as a central regulator of lifespan and aging
F1000Research 8:998.

https://doi.org/10.12688/f1000research.17196.1
- PubMed
- Google Scholar
1. Peters MJ
2. Joehanes R
3. Pilling LC
4. Schurmann C
5. Conneely KN
6. Powell J
7. Reinmaa E
8. Sutphin GL
9. Zhernakova A
10. Schramm K
11. Wilson YA
12. Kobes S
13. Tukiainen T
14. Ramos YF
15. Göring HH
16. Fornage M
17. Liu Y
18. Gharib SA
19. Stranger BE
20. De Jager PL
21. Aviv A
22. Levy D
23. Murabito JM
24. Munson PJ
25. Huan T
26. Hofman A
27. Uitterlinden AG
28. Rivadeneira F
29. van Rooij J
30. Stolk L
31. Broer L
32. Verbiest MM
33. Jhamai M
34. Arp P
35. Metspalu A
36. Tserel L
37. Milani L
38. Samani NJ
39. Peterson P
40. Kasela S
41. Codd V
42. Peters A
43. Ward-Caviness CK
44. Herder C
45. Waldenberger M
46. Roden M
47. Singmann P
48. Zeilinger S
49. Illig T
50. Homuth G
51. Grabe HJ
52. Völzke H
53. Steil L
54. Kocher T
55. Murray A
56. Melzer D
57. Yaghootkar H
58. Bandinelli S
59. Moses EK
60. Kent JW
61. Curran JE
62. Johnson MP
63. Williams-Blangero S
64. Westra HJ
65. McRae AF
66. Smith JA
67. Kardia SL
68. Hovatta I
69. Perola M
70. Ripatti S
71. Salomaa V
72. Henders AK
73. Martin NG
74. Smith AK
75. Mehta D
76. Binder EB
77. Nylocks KM
78. Kennedy EM
79. Klengel T
80. Ding J
81. Suchy-Dicey AM
82. Enquobahrie DA
83. Brody J
84. Rotter JI
85. Chen YD
86. Houwing-Duistermaat J
87. Kloppenburg M
88. Slagboom PE
89. Helmer Q
90. den Hollander W
91. Bean S
92. Raj T
93. Bakhshi N
94. Wang QP
95. Oyston LJ
96. Psaty BM
97. Tracy RP
98. Montgomery GW
99. Turner ST
100. Blangero J
101. Meulenbelt I
102. Ressler KJ
103. Yang J
104. Franke L
105. Kettunen J
106. Visscher PM
107. Neely GG
108. Korstanje R
109. Hanson RL
110. Prokisch H
111. Ferrucci L
112. Esko T
113. Teumer A
114. van Meurs JB
115. Johnson AD
(2015) The transcriptional landscape of age in human peripheral blood
Nature Communications 6:8570.

https://doi.org/10.1038/ncomms9570
- PubMed
- Google Scholar
(2017) Using DNA methylation profiling to evaluate biological age and longevity interventions
Cell Metabolism 25:954–960.

https://doi.org/10.1016/j.cmet.2017.03.016
- PubMed
- Google Scholar
(2014) Full-length RNA-seq from single cells using Smart-seq2
Nature Protocols 9:171–181.

https://doi.org/10.1038/nprot.2014.006
- PubMed
- Google Scholar
1. Rappaport N
2. Nativ N
3. Stelzer G
4. Twik M
5. Guan-Golan Y
6. Stein TI
7. Bahir I
8. Belinky F
9. Morrey CP
10. Safran M
11. Lancet D
(2013) MalaCards: an integrated compendium for diseases and their annotation
Database 2013:bat018.

https://doi.org/10.1093/database/bat018
- PubMed
- Google Scholar
1. Rappaport N
2. Twik M
3. Nativ N
4. Stelzer G
5. Bahir I
6. Stein TI
7. Safran M
8. Lancet D
(2014) MalaCards: a comprehensive Automatically-Mined database of human diseases
Current Protocols in Bioinformatics 47:1–24.

https://doi.org/10.1002/0471250953.bi0124s47
- PubMed
- Google Scholar
1. Rappaport N
2. Twik M
3. Plaschkes I
4. Nudel R
5. Iny Stein T
6. Levitt J
7. Gershoni M
8. Morrey CP
9. Safran M
10. Lancet D
(2017) MalaCards: an amalgamated human disease compendium with diverse clinical and genetic annotation and structured search
Nucleic Acids Research 45:D877–D887.

https://doi.org/10.1093/nar/gkw1012
- PubMed
- Google Scholar
1. Raudvere U
2. Kolberg L
3. Kuzmin I
4. Arak T
5. Adler P
6. Peterson H
7. Vilo J
(2019) G:profiler: a web server for functional enrichment analysis and conversions of gene lists (2019 update)
Nucleic Acids Research 47:W191–W198.

https://doi.org/10.1093/nar/gkz369
- PubMed
- Google Scholar
(2011) Interpretation of random effects meta-analyses
BMJ 342:d549.

https://doi.org/10.1136/bmj.d549
- PubMed
- Google Scholar
1. Schaum N
2. Lehallier B
3. Hahn O
4. Pálovics R
5. Hosseinzadeh S
6. Lee SE
7. Sit R
8. Lee DP
9. Losada PM
10. Zardeneta ME
11. Fehlmann T
12. Webber JT
13. McGeever A
14. Calcuttawala K
15. Zhang H
16. Berdnik D
17. Mathur V
18. Tan W
19. Zee A
20. Tan M
21. Pisco AO
22. Karkanias J
23. Neff NF
24. Keller A
25. Darmanis S
26. Quake SR
27. Wyss-Coray T
(2020) Ageing hallmarks exhibit organ-specific temporal signatures
Nature 583:596–602.

https://doi.org/10.1038/s41586-020-2499-y
- PubMed
- Google Scholar
(1980) Thymus cell migration quantitative aspects of cellular traffic from the Thymus to the periphery in mice
European Journal of Immunology 10:210–218.

https://doi.org/10.1002/eji.1830100310
- PubMed
- Google Scholar
(2019) mTOR and aging: an old fashioned dress
International Journal of Molecular Sciences 20:2774.

https://doi.org/10.3390/ijms20112774
- Google Scholar
(1980) Proliferation kinetics of the mouse bladder after irradiation
Cell Proliferation 13:75–89.

https://doi.org/10.1111/j.1365-2184.1980.tb00451.x
- Google Scholar
1. Subramanian A
2. Tamayo P
3. Mootha VK
4. Mukherjee S
5. Ebert BL
6. Gillette MA
7. Paulovich A
8. Pomeroy SL
9. Golub TR
10. Lander ES
11. Mesirov JP
(2005) Gene set enrichment analysis: a knowledge-based approach for interpreting genome-wide expression profiles
PNAS 102:15545–15550.

https://doi.org/10.1073/pnas.0506580102
- PubMed
- Google Scholar
1. Tabula Muris Consortium
(2020) A single-cell transcriptomic atlas characterizes ageing tissues in the mouse
Nature 583:590–595.

https://doi.org/10.1038/s41586-020-2496-1
- PubMed
- Google Scholar
1. Tacutu R
2. Thornton D
3. Johnson E
4. Budovsky A
5. Barardo D
6. Craig T
7. Diana E
8. Lehmann G
9. Toren D
10. Wang J
11. Fraifeld VE
12. de Magalhães JP
(2018) Human ageing genomic resources: new and updated databases
Nucleic Acids Research 46:D1083–D1090.

https://doi.org/10.1093/nar/gkx1042
- Google Scholar
1. Teta M
2. Long SY
3. Wartschow LM
4. Rankin MM
5. Kushner JA
(2005) Very slow turnover of beta-cells in aged adult mice
Diabetes 54:2557–2567.

https://doi.org/10.2337/diabetes.54.9.2557
- PubMed
- Google Scholar
1. Tower J
(2015) Programmed cell death in aging
Ageing Research Reviews 23:90–100.

https://doi.org/10.1016/j.arr.2015.04.002
- PubMed
- Google Scholar
1. Vijg J
2. Suh Y
(2013) Genome instability and aging
Annual Review of Physiology 75:645–668.

https://doi.org/10.1146/annurev-physiol-030212-183715
- PubMed
- Google Scholar
1. Weichhart T
(2018) mTOR as regulator of lifespan, aging, and cellular senescence: a Mini-Review
Gerontology 64:127–134.

https://doi.org/10.1159/000484629
- PubMed
- Google Scholar
(2013) Closing the gap between T-cell life span estimates from stable isotope-labeling studies in mice and humans
Blood 122:2205–2212.

https://doi.org/10.1182/blood-2013-03-488411
- PubMed
- Google Scholar
1. Ximerakis M
2. Lipnick SL
3. Innes BT
4. Simmons SK
5. Adiconis X
6. Dionne D
7. Mayweather BA
8. Nguyen L
9. Niziolek Z
10. Ozek C
11. Butty VL
12. Isserlin R
13. Buchanan SM
14. Levine SS
15. Regev A
16. Bader GD
17. Levin JZ
18. Rubin LL
(2019) Single-cell transcriptomic profiling of the aging mouse brain
Nature Neuroscience 22:1696–1708.

https://doi.org/10.1038/s41593-019-0491-3
- PubMed
- Google Scholar
Software
1. Zhang MJ
2. Pisco AO
(2021) tabula-muris-senis, version swh:1:rev:0fd2ee501f4f3bd0e691b6071aee2c9286f1cf92
Software Heritage.

https://archive.softwareheritage.org/swh:1:dir:7dc782970300a97e9bca9038ba34728c857a0638;origin=https://github.com/czbiohub/tabula-muris-senis;visit=swh:1:snp:40213ae8e5dd05dbdbdeb6bb45166c8436b328f1;anchor=swh:1:rev:0fd2ee501f4f3bd0e691b6071aee2c9286f1cf92/
1. Zhou X
2. Liao WJ
3. Liao JM
4. Liao P
5. Lu H
(2015) Ribosomal proteins: functions beyond the ribosome
Journal of Molecular Cell Biology 7:92–104.

https://doi.org/10.1093/jmcb/mjv014
- PubMed
- Google Scholar

Article and author information

Author details

Martin Jinye Zhang
1. Department of Electrical Engineering, Stanford University, Palo Alto, United States
2. Department of Epidemiology, Harvard T.H. Chan School of Public Health, Boston, United States
3. Program in Medical and Population Genetics, Broad Institute of MIT and Harvard, Cambridge, United States
Contribution
Conceptualization, Formal analysis, Investigation, Visualization, Methodology, Writing - original draft, Writing - review and editing

For correspondence
jinyezhang@hsph.harvard.edu

Competing interests
No competing interests declared

"This ORCID iD identifies the author of this article:" 0000-0003-0006-2466
Angela Oliveira Pisco

Chan-Zuckerberg Biohub, San Francisco, United States

Contribution
Conceptualization, Data curation, Supervision, Investigation, Writing - original draft, Writing - review and editing

For correspondence
angela.pisco@czbiohub.org

Competing interests
No competing interests declared

"This ORCID iD identifies the author of this article:" 0000-0003-0142-2355
Spyros Darmanis

Chan-Zuckerberg Biohub, San Francisco, United States

Contribution
Data curation, Writing - review and editing

Competing interests
No competing interests declared
James Zou
1. Department of Electrical Engineering, Stanford University, Palo Alto, United States
2. Chan-Zuckerberg Biohub, San Francisco, United States
3. Department of Biomedical Data Science, Stanford University, Palo Alto, United States
Contribution
Supervision, Funding acquisition, Investigation, Methodology, Writing - original draft, Writing - review and editing

For correspondence
jamesyzou@gmail.com

Competing interests
No competing interests declared

"This ORCID iD identifies the author of this article:" 0000-0001-8880-4764

Funding

Chan-Zuckberg Biohu

James Zou

National Science Foundation (CCF 1763191)

James Zou

National Institutes of Health (R21 MD012867-01)

James Zou

National Institutes of Health (P30AG059307)

James Zou

Silicon Valley Community Foundation

James Zou

National Institutes of Health (R01 MH115676)

Martin Jinye Zhang

The funders had no role in study design, data collection and interpretation, or the decision to submit the work for publication.

Acknowledgements

We would like to thank S Quake, R Sinha, R Sit, J Cool, B van de Geijn, H Shi, and X Xu for feedback. MJZ and JZ are supported by NSF CCF 1763191, NIH R21 MD012867-01, NIH P30AG059307, and grants from the Silicon Valley Foundation and the Chan-Zuckerberg Initiative. MJZ is also supported by NIH R01 MH115676.

Copyright

This article is distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use and redistribution provided that the original author and source are credited.