Dynamic regulation of genetic pathways and targets during aging in Caenorhabditis elegans

Identification of significantly aging-dependent pathways

Based on the GSEA approach, we have identified several significantly related genetic pathways during aging in C. elegans. There were significantly 42 up- and 42 down-regulated pathways from the stage of L4 to D6, as well as 26 up- and 56 down-regulated pathways from the stage of D6 to D15 (Figure 1A and Additional file 1). By the comparison of each group, we have also identified some overlapping pathways significantly associated in both of these two periods, which were detailed in Table 1 and Figure 1B. In summary, there were totally 17 overlapping pathways significantly up-regulated from L4 to D6 but down-regulated from D6 to D15, and most of them were related to Signal Transduction (7/17) and Lipid Metabolism (4/17) according to KEGG pathway maps in the KEGG database (http://www.genome.jp/kegg/). Interestingly, 6 different signaling pathways were significantly par-ticipating in this dynamic regulation, including MAPK signaling pathway, mTOR signaling pathway, Wnt signaling pathway, TGF-beta signaling pathway, ErbB signaling pathway and Hedgehog signaling pathway.

Moreover, there were totally 36 overlapping pathways significantly down-regulated both from L4 to D6 and from D6 to D15, and most of them were related to Amino Acid Metabolism (9/36) and Carbohydrate Metabolism (12/36); there were totally 19 overlapping pathways significantly up-regulated both from L4 to D6 and from D6 to D15, and most of them were related to Translation (4/19) and Transcription (3/19) as well as Replication and Repair (5/19). All of these data revealed the dynamic re-gulation in multiple pathways during aging in C. elegans.

Dynamic regulation in multiple signaling pathways

In the study of Youngman MJ et al. [14], a decline in PMK-1 p38 mitogen-activated protein kinase (MAPK) pathway, one of conserved pathways involved in pathogen defense, was observed by analyzing gene expression levels using unpaired t test in synchronized populations of N2 worms during mid-to-late adulthood (at Day 6 of adulthood and Day 15 of adulthood). It was consistent with our result on MAPK signaling pathway, which was identified as significantly down-regulated from D6 to D15 of adulthood (highlighted in red color in Additional file 1). Surprisingly, MAPK signaling pathway was shown to be significantly up-regulated during young-to-mid adulthood (at the late larval L4 stage and Day 6 of adulthood) in our study, which was not reported by Youngman MJ et al. study [14]. Furthermore, MAPK signaling pathway was identified as truly aging-dependent pathways based on the hierarchical clustering of both samples and involved genes (Figure 2A for the data from L4 to D6 and Fig2B from D6 to D15). Obviously, using the expression data of involved genes in MAPK signaling pathway, six samples were clustered into two groups with samples at the same stage (L4-1, L4-2 and L4-3 in group of L4; D6-1, D6-1 and D6-3 in group of D6; D15-1, D15-1 and D15-3 in group of D15). Actually, MAP kinase (MAPK)/Receptor Tyrosine Kinase (RTK)/RasGTPase signaling pathways were used repeatedly during metazoan development, which were reported to be required for larval viability and for many different developmental processes, including induction of vulval, uterine, spicule, P12 and excretory duct cell fates, control of sex myoblast migration and axon guidance, and promotion of germline meiosis [15]. The core components or regulators of C. elegans RTK/Ras/MAPK signaling pathway were also identified as the aging-dependent targets [16-18], such as LET-23, LIN-1 and LIN-45 (Figure 2 and Additional file 2). Among these candidate targets, four interesting genes, including ced-2, lin-45, ZK856.8 and F59D6.7, were significantly up-regulated from L4 to D6 and down-regulated from D6 to D15, showing dynamic expression patterns. Their significances of gene expression changes between two adjacent stages were shown in Table 2. In C. elegans, the ced-2 gene encodes a Src homology (SH) 2 and 3-containing adaptor protein, homologous to human CrkII, which is required for phagocytosis during programmed cell death and for migration of the distal tip cells of the somatic gonad. It was previously reported that ced-2 mutants may exhibit persistant corpses indicating defects in phagocytosis and distal tip cell migration defects in the gonadal arms [19]. The lin-45 gene encodes an ortholog of the vertebrate protein RAF which is required for larval viability, fertility and the induction of vulval cell fates [18]. Based on RNAi studies of three genes constituting the ERK cascade (lin-45/RAF1, mek-2/MEK1/2, andmpk-1/ERK1/2) as well as lip-1 encoding a MAPK phosphatase that inactivates MPK-1, a novel ERK-MAPK-mediated signaling pathway has been identified, which promotes longevity through two candidate transcription factors, SKN-1 and DAF-16 [20]. Additionally, the genes of ZK856.8 and F59D6.7 would be two important targets for further study of aging in MAPK signaling pathway.

The mammalian target of rapamycin (mTOR), also known as mechanistic target of rapamycin has been identified as a key modulator of ageing and age-related disease. The fact that inhibition of mTOR signaling pathway may extend lifespan in model organisms and confer protection against a growing list of age-related pathologies has been clearly approved [21]. In our study, this pathway was also identified as one of truly aging-dependent pathways with dynamic regulation, which was significantly up-regulated from L4 to D6 and down-regulated from D6 to D15 (Table 1). Like MAPK signaling pathway, based on the expression of candidate genes in mTOR signaling pathway, the samples at the same time point can be clustered together (Figure 3A for the data from L4 to D6 and Figure 3B from D6 to D15). There were 24 genes involved in mTOR signaling pathway both in the period of young-to-mid and mid-to- late adulthood, which were both clustered into 2 groups (the group of A and B in Figure 3A & 3B). Among these genes appearing in both the group of B in Figure 3A and the group A in Figure 3B, besides lin-45 mentioned above, the gene of pdk-1 was identified with dynamic expression pattern, which was significantly up-regulated from L4 to D6 (p=4.58E-04) and down-regulated from D6 to D15 (p=1.19E-02) (Table 2).In C. elegans, pdk-1 encodes the 3-phosphoinositide-dependent kinase 1 ortholog (PDK-1), which is a component of the DAF-2/insulin receptor-like signaling pathway and accordingly, functions to regulate such processes as dauer larvae formation as well as longevity [21]. The increased lifespan in daf-2 insulin/IGF-1 receptor mutants is accompanied by up-regulation of the MDL-1 Mad basic helix-loop-helix leucine zipper transcriptionfactor, which is an inhibitor of cell proliferation and growth that slows progression of an age-related pathology in C. elegans [22]. The previous studies have illustrated that mTOR signaling pathway plays the central role, whereas other pathways act on it by activating or antagonizing either at upstream or downstream [23-26].

Wnt signaling pathways have extremely diverse functions in animals, including the control of gene expression, cell behavior, cell adhesion, and cell polarity [19]. During aging in C. elegans, Wnt signaling pathway was identified as another truly aging-dependent pathway with dynamic regulation, which was significantly up-regulated from L4 to D6 and down-regulated from D6 to D15. There were 62 involved genes in Wnt signaling pathway, which were clustered into 4 groups from L4 to D6 (the group from A to D in Figure 4A) and clustered into 3 groups (the group from A to C in Figure 4B). Among these genes appearing in both the group of A in Figure 4A and the group C in Figure 4B, there were 8 dynamically regulated genes containing sma-4, cyd-1, skr-7,9,13,14,F59D6.7 and ZK856.8, which were significantly up-regulated from L4 to D6 and down- regulated from D6 to D15(Table 2).The gene of sma-4 encodes a Smad protein and a homolog of human DPC4, the protein of SMA-4 is similar to members of the vertebrate protein family of Dwarfins. During development, sma-4 functions as part of a DBL-1/SMA-6 TGF-beta-related signaling pathway that controls body size and male tail sensory ray and spicule formation as well as regulates reproductive aging via oocyte and germline quality maintenance [27]. Indeed, TGF-beta signaling pathway was one of our candidate aging-dependent pathways (Figure 5). The Smad-mediated TGF-beta signaling pathway controls numerous cellular responses from cell proliferation, differentiation and extracellular matrix remodelling to embryonic development in species ranging from worms to mammals [20]. It has been shown that TGF-beta superfamily ligands can regulate cellular or physiological processes through non-canonical pathways by the activation of other signaling molecules, e.g. Akt, MAPK, mTOR, and Src independent of Smad proteins [28]. In this study, Wnt signaling pathway would be seen as one of Smad-dependent pathways in TGF-beta family signaling based on the common regulator of SMA-4 (Table 2).

ErbB signaling pathway was the last one of our truly aging-dependent signaling pathways with dynamic expression pattern during aging (Figure 6). The ErbB family of receptor tyrosine kinases (RTKs) couples binding of extracellular growth factor (EGF) ligands to intracellular signaling pathways regulating diverse biologic responses, including proliferation, differentiation, cell motility, and survival [29]. In this pathway, lin-45 was the only one significantly regulated gene, which also appearing in both MAPK signaling pathway and mTOR signaling pathway (Table 2). lin-45 is a component of an EGFR-mediated inductive signaling pathway that causes vulva precursor cells (VPCs) to generate the vulva [30]. Expression of LIN-45(V627E) and LIN-45(ED) was previously identified to occasionally cause the vulva to burst at the L4-to-adult molt, which may reflect defects in either the VPC specification or the resulting vulval cells [31]. Although Hedgehog signaling pathway has been identified as one of the significantly aging-dependent pathways, it was not the true one due to the unsuccessful clustering on samples from D6 to D15 by using involved genes (Additional file 3).

Generally, there were 12 significant genes with dynamic regulation patterns during aging in above five truly significant signaling pathways (Figure 7). They were ced-2, lin-45, ZK856.8, F59D6.7, pdk-1, cyd-1, sma-4, andskr-7, 8, 9, 13, 14, all of which were significantly up-regulated from L4 to D6 and down-regulated from D6 to D15. The information of significances and involved pathways were shown in Table 2. These targets would be regarded as key factors in future studies of genetic mechanisms of aging. Interestingly, 5 members of skr gene family (skr-7, 8, 9, 13, 14) were included. skr-7 encodes a homolog of Skp1 that functions within a particular SCF ubiquitin-ligase (E3) complex by binding both a cullin (a homolog of S. cerevisiae Cdc53) and an F box protein (through direct interaction with an F box motif), which is required for posterior body morphogenesis, embryonic and larval development, and cell proliferation [32, 33]. The most closely related paralogs of skr-7, -8 and -9 were highly similar, but it is not clear whether these genes comprise a functionally redundant set. The function of SKR-13 is not necessarily confined to, or even partially involved with C. elegans ubiquitin-ligase complexes. In two-hybrid assays, SKR-13 does not bind to any known C. elegans cullins (CUL-1 through CUL-6) [33]. The most closely related paralogs of skr-13 in the C. elegans genome are skr-12 and skr-14.

Continued decline in metabolism related pathways

Besides signal transduction, several measurements of metabolic activity have been demonstrated to display age-related changes. For example, about 50% of declines in carbon dioxide generation by populations of live worms grown on plates measured by gas respirometry have been found from day 6 to day 12 [34]. Additionally, about 60% of declines in oxygen consumption of populations of live worms in liquid culture monitored polarographically using Clark type electrodes have also been evaluated from adult day 0 to day 12 [35, 36]. By contrast, it was indicated that a reduced metabolic rate may cause the extended lifespan of worms based on the result of reduced metabolic rate of age-1 and daf-2 mutants compared to wild-type animals [34]. In our study, most of metabolism related pathways were significantly identified with continued decline in expression pattern both from L4 to D6 and from D6 to D15.It has been reported that the couplingof a metabolic shift with endotoxin detoxification results in extreme longevity following mekk-3 knock-down, a mammalian MEKK3-like kinase gene in C. elegans [37]. An age-dependent decreased abundance of proteins involved in mitochondrial function, electron transport chain, citric acid cycle, and fatty acid metabolism as well as increased abundance of proteins involved in glycolysis and oxidative stress response have been revealed by the quantitative proteomics and pathway analysis [38].

Multi-level gene regulatory networks in signaling pathways during aging

There were totally 154 and 158 individual genes in above five truly significant signaling pathways, respectively from L4 to D6 and from D6 to D15. The expression profiles of these genes were shown in Additional file 4. Based on ChIP–seq data analysis to explore the binding site targets of individual factors, 2530, 5974 and 3178 target genes were identified for transcription factors (TFs) of FOS-1, NHR-28 and UNC-62, respectively (Additional file 4). To dissect the association among these regulators or between regulators and target genes in multiple signaling pathways both from L4 to D6 and from D6 to D15, multi-level gene regulatory networks were constructed for these two developmental stages. As a result, NHR-28 was commonly identified at the top of the network, which regulated UNC-62 and FOS-1, the regulator of FOS-1 was regulated by UNC-62 at the same time (Figure 8). HOX co-factor UNC-62 (Homothorax) has been identified as a developmental regulator that binds proximal to age-regulated genes and modulates lifespan by integrating RNAi and genomics approach [39]. C. elegans FOS-1 has been shown to act in uterine and vulval development, which regulates plc-1 expression in the spermatheca to control ovulation [4]. Moreover, 64 interactions between regulatory factors and target genes were identified at the first stage, including 26 genes regulated by NHR-28, 21 genes by UNC-62 and 17 genes by FOS-1. By contrast, 79 interactions between regulatory factors and target genes were identified at the second stage, including 31 genes regulated by NHR-28, 27 genes by UNC-62 and 21 genes by FOS-1 (Additional file 5). 32 common target genes were identified at both stages, but 3 genes only appears at the first stage (ruvb-1, let-60 and siah-1 in Figure 8A) and 10 genes only appears at the second stage (cdc-42, kin-19, mom-5,ife-2, mek-2, C07E3.9, gpa-12, sos-1, let-23and skr-21 in Figure 8B).

In this study, we employed the gene expression datasets of three time points during aging in C. elegans and performed the approach of GSEA on each dataset between two adjacent stages. As a result, multiple genetic pathways and targets were identified as significantly down- or up-regulated. Among them, 5 truly aging-dependent signaling pathways including MAPK signaling pathway, mTOR signaling pathway, Wnt signaling pathway, TGF-beta signaling pathway and ErbB signaling pathway as well as 12 significantly associated genes were identified with dynamic expression pattern during aging. On the other hand, the continued declines in the regulation of several metabolic pathways have been demonstrated to display age-related changes. Furthermore, the reconstructed regulatory networks based on three of aging related ChIP–seq datasets and the expression matrices of 154 involved genes in above signaling pathways provided new insights into aging at the multiple pathways level. In conclusion, the identification of these genetic pathways and target genes by our integrated analysis would support the hyperfunction theory, which is a plausible alternative to the molecular damage theory to explain aging in C. elegans [1]. The combination of multiple genetic pathways and targets needs to taken into consideration in future studies of aging, in which the dynamic regulation would be uncovered.

Gene expression data collection and pre-processing

The gene expression dataset of GSE21784 related to aging was collected from the database of GEO (www.ncbi.nlm.nih.gov/geo/), which was contributed by Youngman MJ et al. [14]. In this study, RNA samples were obtained from wild type C. elegans strain Bristol N2 during aging and then prepared for hybridization to Affymetrix C. elegans Genome Array (GPL200) from synchronized populations of C. elegans at three points: the late larval L4 stage (L4, 1 day pre-adulthood), day 6 of adulthood (D6), and day 15 of adulthood (D15). There are three biological replicates for each time point: 3 for L4 (marked with L4-1, L4-2 and L4-3 in the following analysis), 3 for D6 (D6-1, D6-1 and D6-3), as well as 3 for D15 (D15-1, D15-1 and D15-3). We then divided the datasets into two subsets, including the study of period from L4 to D6 (containing 6 samples: L4-1, L4-2, L4-3, D6-1, D6-1 and D6-3) and the study of period from D6 to D15 (also containing 6 samples: D6-1, D6-1, D6-3, D15-1, D15-1 and D15-3). It aims to uncover the diverse mechanism of aging process.

To assess the influence of preprocessing on the comparison, data preprocessing was performed using software packages developed in version2.6.0 of Bioconductor and R version 2.10.1. Each Affymetrix dataset was background adjusted, normalized and log2 probe-set intensities calculated using the Robust Multichip Averaging (RMA) algorithm in Affy package [41].

The analysis of significantly enriched pathways and genes

The approaches of Gene Set Enrichment Analysis (GSEA) were performed using Category package in version 2.6.0 of Bioconductor to identify the significantly enriched pathways and genes in both above two subsets [42, 43]. In our analysis, the gene sets with less than 10 genes were excluded. The t-statistic mean of the genes was computed in each KEGG (Kyoto Encyclopedia of Genes and Genomes) pathway. Using a permutation test with 1000 times, the cutoff of significance level p value was set as 0.01 for the most significant pathways related to aging. Accordingly, the significant pathways and genes were then identified under the comparison between different time points, including from L4 to D6 and from D6 to D15. The following classification of identified pathways was based on the KEGG pathway maps br08901 of BRITE Functional Hierarchies in KEGG database (http://www.genome.jp/kegg-bin/get_htext?br08901.keg). The annotation of significant C. elegans genes in each pathway was based on the WormBase database, of version WS236 (www.wormbase.org). Next, clustering on groups and genes was performed based on the identified genes' expression in each significant pathway using the method of hierarchical clustering with Euclidean distance. The significance of associated target genes from L4 to D6 or from D6 to D15 was calculated by unpaired t test.

The collection and analysis of aging related ChIP–seq data

In the database of modENCODE, there are almost one hundred available ChIP–seq data associated with individual transcription factors at diverse developmental stages. Recently, Van Nostrand et al. have performed integrative analysis of these data sets to infer gene regulatory interactions [44]. As a result, nine candidate aging regulators (NHR-28, PQM-1, FOS-1, NHR-76, NHR-77, ELT-3, C01B12.2, UNC-62 and SKN-1) were identified, of which three were previously known to have a role in longevity (ELT-3, SKN-1 and UNC-62). However, there are only 3 datasets at the similar stage (at or after L4 stage) of Matthew JY et al.'s gene expression study, including NHR-28 and FOS-1 at L4 stage (L4) as well as UNC-62 at the stage of young adult day 4 (YD). Therefore, these three ChIP–seq data were then downloaded and analyzed to explore the binding site targets of individual factors. Raw sequencing reads were aligned against C. elegans genome (WS232) using bowtie2 [45]. MACS2 was used for peak calling with default parameters [46].

Construction of regulatory networks of signaling pathway genes

The tool of Constructing Multi-level Gene Regulatory Networks (CMGRN) was used to integrate NHR-28 (L4), FOS-1 (L4) and UNC-62 (YD) ChIP-seq information of transcription factor (TF) binding with expression profiling of signaling pathway genes that significantly enriched during aging. CMGRN(http://bioinfo.icts.hkbu.edu.hk/cmgrn/) is based on a methodology combining Bayesian hierarchical model with Gibbs Sampling implementation to identify correct sets of regulator-gene interactions which explore underlying network structures [47]. We aimed to dissect the association between these regulators or between regulators and genes as well as to reconstruct regulatory networks in signaling pathways involved in aging processes from L4 to D6 and from D6 to D15. Firstly, the file of sequence counts of our selected 3 regulators (NHR-28, FOS-1 and UNC-62) on the target genes was used here as the input of ChIP–seq gene counts. Secondly, the target gene symbols containing sequences of each regulator measured by ChIP–seq was used as the input of binding data. The last input file is the gene expression profile. Here, we used the gene expression values of involved genes in the significantly identified pathways by GSEA of each sample, respectively for the stages from L4 to D6 and from D6 to D15. Then the parameters were set as defaults: the minimal value of interaction was selected as 10, the probability cutoff for regulator-gene interaction was selected as 70% and the causal relationship percentage was selected as 20%.