Capture of gene–gene relationships from transcriptome dataTo investigate age-related changes in gene co-expression, we made use of RNA-seq data collected postmortem by the Genotype-Tissue Expression (GTEx) consortium14 from 30 different human tissues, spanning 948 donors aged between 20 years and 79 years. We first aimed to illustrate age-related changes in gene expression coordination within and between modules. For this, we hand-picked five sets of genes, constructed based on Gene Ontology (GO) term membership. This selection includes both cell-type-specific (for example, antigen binding) and ubiquitous (for example, mitochondrial respiratory chain) gene sets as well as both broad regulatory groups (for example, extracellular matrix (ECM) components) and protein complexes (for example, polymerase II (Pol-II) core complex). We computed all pairwise Pearson correlations in samples collected from either young (20–29 years) or old (60–69 years) donors. We focused our analyses on two tissues with differences in cellular composition and function—Brain and Blood—and on data from all tissues pooled together (Methods). We observed strong within-module correlations across ages and tissues (Fig. 1a, original data, diagonal blocks, and Fig. 1b, colored), in line with a modular organization of gene expression. Additionally, we observed between-module correlations (Fig. 1a, original data, off-diagonal blocks, and Fig. 1b, beige), particularly between genes encoding for members of the mitochondrial respiratory chain and the RNA Pol-II core complex (Fig. 1a,b), representative of regulatory crosstalk between functional modules. Both within-module and between-module correlations differed between tissues and age groups. Notably, we observed age-related differences in gene co-expression both specific to individual tissues (Fig. 1a,b, antigen binding – ECM in blood) and shared across tissues (Fig. 1a,b, mitochondrial respiratory chain in the cross-tissue analysis). These examples confirm the notion that gene co-expression relationships change with age. However, reduced co-expression between individual gene pairs does not necessarily result from reduced expression coordination. Instead, it may also result from adaptive changes requiring a different regulation of the respective genes.Fig. 1: Representative gene–gene relationship changes with age, captured by pairwise correlation and a GRN-based approach.a, Heatmap of pairwise Pearson correlation coefficients for selected cellular functions (modules) and tissues. Correlations were computed based on the original expression data (Original) or expression data reconstructed based on model predictions from c (Reconstructed), in young (20–29 years) or old (60–69 years) samples from Brain, Blood or pooled GTEx tissues (Cross-tissue). b,d, Quantification of within-module (colored) and between-module (beige) correlations for selected modules, in young (opaque) and old (transparent) samples. Correlations were computed based on the original expression data (b) or expression data reconstructed based on model predictions from c. d, Between-module correlations are presented for selected module pairs. A two-sided Mann–Whitney test was used to identify modules with age-related correlation changes: MRC (P = 0.0034 in the Original Brain data, P = 0.000094 in the Reconstructed Brain data, P = 0.018 in the Reconstructed Blood data, P = 0.047 in the Original Cross-tissue data, P = 0.00033 in the Reconstructed Cross-tissue data); MRC – Pol-II (P = 0.055 in the Original Brain data, P = 0.030 in the Reconstructed Brain data, P = 0.00083 in the Original Cross-tissue data, P = 0.0068 in the Reconstructed Cross-tissue data); Antigen binding – ECM (P = 0.023 in the Original Brain data, P = 4.7 × 10−10 in the Original Blood data, P = 0.00028 in the Reconstructed Blood data, P = 0.072 in the Original Cross-tissue data). c, Methodological approach used to capture gene–gene relationships in our regulatory model. Through a combination of regularized linear regression (LASSO) and stability selection, we identified stable predictors for each gene in the transcriptome—that is, genes whose expression pattern across the training data (cancer cell line transcriptomic data) is informative of the expression pattern of the target gene. A linear model was then fit to explain the expression pattern of the target gene (a) based on the pattern of the stable predictors (b and c). The weights of this linear model can then be used in other datasets to reconstruct the expression pattern of the target gene a based on the expression pattern of the stable predictors b and c observed in those datasets. Although the scheme shows two predictors, b and c, for target gene a, the number of predictor genes is not limited; rather, the optimal number of predictors is determined individually for each target gene (Methods). For illustrative purposes, gene sets were restricted to the following functions: respiratory chain complex members I–IV (MRC, GO:0045271, GO:0005749, GO:0005750 and GO:0005751, orange), components of collagen-containing ECM (GO:0062023, purple), Pol-II core complex members (GO:0005665, pink), nucleosome members (GO:0000786, dark green) and peptide antigen binding partners (GO:0042605, brown). ***P < 0.001; **P < 0.01; *P < 0.05; .P < 0.1. cor., correlation; Mito., mitochondrial.Therefore, we next aimed to systematically map out robust regulatory expression relationships between genes that better reflect functional requirements for coordinated expression than simple pairwise co-expression. To capture robust regulatory relationships that explain expression coordination across different tissues, we inferred gene–gene relationships using transcriptomic data collected from 1,443 cancer cell lines15,16. The benefit of using in vitro gene expression data for model training, as opposed to using tissue data, is the fact that such data do not result from a mix of different cell types and are not affected by changes in cell type composition. As opposed to single-cell data, bulk in vitro data provide greater coverage and sensitivity. To derive gene expression relationships from large transcriptomic datasets, we devised a procedure that selects a small set of genes whose expression patterns are, as a linear combination, most predictive for the expression of a given target gene (as in refs. 11,17,18). The goal of this model is not to predict the global average expression level of a gene but the deviation from this global average in a given individual—that is, predicting inter-individual variation in gene expression. The network resulting from this procedure is much sparser than computing all pairwise correlations between all genes: only a very small fraction (0.04%) of gene pairs has non-zero relationships (ignoring directionality). These gene–gene relationships do not exclusively represent causal relationships between regulators and their targets. Instead, they represent gene pairs that are robustly co-expressed across a wide range of human tissues and cell types16.Next, we assessed the capacity of our model to capture known direct regulatory relationships between regulatory TFs and their target genes. We collected gene sets corresponding to TF targets from the Dorothea database19 and tissue-specific TF–target relationships inferred from the GTEx data18. For each of those sets, we computed edge connectivity in our network (Methods) and observed that the connectivity of genes regulated by the same TF was consistently and significantly higher than the background network connectivity (Extended Data Fig. 1a). This was true for both the cross-tissue (Dorothea sets) and tissue-specific (GTEx tissue sets) TF target gene sets. This result supports the notion that the connectivity patterns in our network capture generic and tissue-specific regulatory programs.Using the gene–gene relationships captured by our regulatory model, we reconstructed the expression pattern of each gene across samples of the same tissue (Fig. 1c and Methods). Thus, the reconstructed data consist of the expected (predicted) pattern for each gene given the expression pattern observed for its regulatory neighbors (that is, directly connected genes in the network, corresponding to the most strongly co-regulated gene pairs). First, we observed that tissue-specific differences in the gene–gene correlation patterns observed in the original data (Fig. 1a) were largely preserved in the reconstructed data (Extended Data Fig. 3)— for instance, the weaker correlation between a subset of the members of the mitochondrial respiratory chain in Blood but not Brain. This was true not only for within-module correlations but also for correlations between different modules, as is the case between genes encoding for peptide antigen binding partners and ECM components. This observation confirms that our network captures regulatory neighborhoods that are conserved across cell types and tissues and, at the same time, enables prediction of tissue-specific expression levels. Second, we observed that age-related co-expression changes observed in the original data (Fig. 1a,b) were also captured in the reconstructed data (Fig. 1d and Extended Data Fig. 3). This observation is important because it supports the notion that co-expression changes (especially co-expression reduction) do not necessarily result from ‘sloppy regulation’ but may reflect regulated, physiologically plausible cellular responses.The results shown in Fig. 1 are limited to a small subset of hand-picked genes and tissues to illustrate aspects of age-related changes in gene expression coordination, but they were not intended to replace a systematic, transcriptome-wide analysis. Therefore, we next tested the predictive power of our network model across all genes and tissues present in the GTEx dataset. Inter-individual differences in gene expression were predicted better than random guessing for the vast majority of genes and tissues (Extended Data Fig. 1b, blue distributions, and ref. 16). Thus, the network model captures gene–gene relationships that are conserved across a large diversity of cellular states. However, we note that the capacity of our model to correctly capture relevant regulatory neighborhoods is dependent on the gene and dataset at hand16. To gain a better understanding of the limitations of our model, we identified genes whose observed expression profiles differ strongly from the predicted based on their regulatory neighborhood (Methods). Our analyses revealed that most of these poorly predicted genes fell into one of two groups. The first group consisted of lowly expressed genes with few regulatory neighbors, which themselves were also lowly expressed (Extended Data Fig. 4). Thus, these genes correspond to network regions that were essentially inactive in the respective tissue or cell type, and, therefore, the observed expression variation reflected mostly background noise. The second group corresponded to highly expressed genes with extremely low variance and ubiquitous expression across tissues (Extended Data Fig. 4), suggesting that these are housekeeping genes with essentially constant expression. Because our model predicts expression variation across samples (that is, relative differences in expression between different samples), a constant expression results in virtually no relative differences between samples and is, thus, hard to model. Additionally, we functionally characterized genes based on the agreement/disagreement between observed and predicted expression levels across tissues (Methods). We consistently found genes involved in general cellular functions, such as RNA binding, translation and oxidative phosphorylation, to be well predicted across tissues, whereas genes encoding for signaling receptors were poorly predicted (Extended Data Fig. 4f). In subsequent analyses, we excluded genes that were not predictable by our model in a given tissue—that is, we used our network model only in a gene–tissue context where it is actually applicable.Next, we asked how much the model predictions would improve if tissue-specific networks were available. We, therefore, retrieved publicly available scRNA-seq data collected from peripheral blood mononuclear cells (PBMCs) from healthy donors20 and trained a blood-specific model using the same strategy as before. Compared to the cross-tissue (cancer cell line–based) network, the blood network performed worse on all GTEx tissues, except for GTEx blood, where both networks performed equally well. (Methods and Extended Data Fig. 1b). Further comparison of these two networks revealed a much higher coverage of the cross-tissue network (16,624 predictable genes) compared to the blood-specific network (4,141 predictable genes). The expression pattern of genes involved in blood-specific regulatory functions was better predicted by the blood-specific network, whereas the cross-tissue network performed better for genes involved in basic cellular functions (Extended Data Fig. 2). We concluded from this analysis that most gene–gene relationships captured by the cross-tissue network are conserved across cell types, and the potential benefit of having tissue-specific or cell-type-specific networks is partially offset by the larger training dataset available for the cross-tissue network.Age-related changes in gene expression coordinationNext, we systematically investigated changes in gene regulatory programs across the entire transcriptome and across human tissues. Here, we aimed to separate those changes in regulatory programs from changes in gene expression levels. As described above, our network model consists of gene–gene relationships that are invariant across a large diversity of cellular states. Thus, our notion was that changes in expression coordination would modify the ability of our network model to correctly predict the expression of a gene. Whenever the regulation of a gene ‘aligns’ with the model structure, it will be predictable, whereas regulatory inputs deviating from the model structure would reduce its predictability.To quantify such changes in regulatory programs, we established a metric of (de)regulation, based on the extent of agreement between the expression profile of a given gene in the data and the reconstructed profile according to its regulatory neighborhood (that is, comparing observed versus predicted expression levels). To quantify this agreement, we used Spearman’s rho, computed between the expected and observed expression of a gene across samples. We refer to this metric as the ‘predictability’ of a gene in each group of samples (Fig. 2a). Our choice of Spearman over Pearson correlation was based on its lower sensitivity to outliers, which avoids that a few individuals with outlying expression of the target gene dictate the overall predictability. To capture age-related changes in predictability, we made use of the wide age range of the donors in the GTEx dataset. We restricted our analysis to the 20 tissues with highest sample number and split the data into six similarly sized age groups, defined as age decades (20–29 up to 70–79; Methods). Splitting the samples by decade resulted in 30–62 samples per tissue and age split (Supplementary Table 1). We then computed the predictability of each gene in every tissue–age split and restricted our analysis to genes showing sufficiently high average predictability across age groups (Methods). This resulted in the selection of 3,291–5,830 genes per tissue (Supplementary Table 3). We then regressed predictability of these genes as a function of age using a linear model (Predictability ~ Age; Fig. 2a). This procedure resulted in one slope for each gene in each tissue, quantifying the change of its predictability with age. To determine the statistical significance of the resulting predictability slopes, we compared the P value distribution for the obtained regression slope with the null P value distribution, obtained after repeatedly shuffling the age groups (Methods). We observed large differences in the predictability signatures of different tissues, with some tissues showing an enrichment in small P values but not others (Extended Data Fig. 5). Notably, we did not observe a clear global trend toward decreasing (or increasing) predictability with age. Although some tissues (for example, Blood and Thyroid) had more genes with decreasing predictability, other tissues had more genes with increasing predictability (for example, Breast – mammary tissue, Nerve and Esophagus – mucosa). We focused our subsequent analysis on eight tissues with substantially more genes showing an age-dependent predictability change than expected by chance—that is, those tissues with an inflation of small P values (Adipose – visceral, Artery – tibial, Blood, Brain, Breast – mammary tissue, Esophagus – mucosa, Testis and Thyroid; Fig. 2b and Supplementary Table 1). First, we analyzed predictability changes of 370 genes that could be analyzed across all of those eight tissues, because their average predictability was sufficiently high in all of them (Fig. 2c). When comparing the predictability slopes of those 370 genes, we observed that many of them showed similar predictability changes across multiple tissues. This was especially true for genes with age-related decrease in predictability (Fig. 2c, left). When grouping tissues according to their age-related predictability change profiles, we observed a separation of Testis, Breast and Brain from the remaining tissues, accompanied by a more predominant age-related increase in predictability in these tissues. The distinct predictability changes in these tissues are consistent with earlier analysis of the GTEx cohort revealing regulatory programs distinguishing Brain and Testis from the other tissues21.Fig. 2: Transcriptome-wide predictability changes with age across tissues.a, Computational approach used to identify predictability changes with age. Predictability is quantified as the Spearman correlation between the observed expression patterns (in the original data) and the predicted expression patterns (in the reconstructed data), corresponding to the expected pattern given the expression of regulatory neighbors. A high correlation indicates that the expression pattern of a gene fits the regulatory relationships captured by the model (top left), whereas a low correlation indicates the opposite. Predictability is quantified for groups of samples at six different age groups: the decades spanning 20–29 to 70–79. For each gene, predictability is modeled as a linear function of age. b, Distribution of P values for the regression of predictability values within each age group against the mean age of the age group. P values correspond to the two-sided t-test on regression coefficients, without multiple testing correction. Red line: P value distribution obtained from the real data. Gray background: average P value distribution across 100 permutations of the age groups. Black lines: five individual permutations randomly picked from the background. The dashed vertical line indicates the highest P value among the genes considered statistically significant in each tissue (orange and blue bars in e). The number of genes included in each tissue-specific analysis can be found in Supplementary Table 3. c, Heatmap of the predictability slopes across all 370 genes, independently of significance level, ordered by increasing average predictability across tissues. Only the 370 genes for which the regression analysis was performed in all eight tissues were included—that is, the genes with a high average predictability in all eight tissues (Supplementary Table 3). d, Hallmark gene sets enriched in age-related gene–gene relationship changes, captured by GSEA. The heatmap shows all hallmarks with statistically significant (FDR < 0.05) enrichment in at least one tissue. e, Number of genes with predictability increase (blue) and decrease (orange) among the top 100 most significant genes per tissue. ***FDR < 0.001; **FDR < 0.01; *FDR < 0.05; .FDR < 0.1. FDR, false discovery rate.We next sought to understand which cellular functions were most affected by the observed age-related predictability changes, either commonly across tissues or specifically for each tissue. To this end, we performed gene set enrichment analysis22 (GSEA) on the age-related predictability slopes for all genes included in the analysis in each tissue, to identify gene sets with significant enrichment in age-related predictability changes in each tissue (Methods). We focused our analysis on gene sets (Molecular Signatures Database (MSigDB) hallmarks, Fig. 2d, or GO terms, Extended Data Fig. 7) showing a significant enrichment in at least one tissue. This analysis revealed that genes involved in oxidative phosphorylation were among the most affected by predictability changes across several tissues, showing a predictability decrease in all tissues apart from Brain and Testis (Fig. 2d and Extended Data Fig. 7). Targets of the E2F TF family, known to regulate the expression of genes involved in the transition between the G1 and S phases of the cell cycle, and targets of the oncogene MYC, also involved in cell cycle regulation, apoptosis and differentiation, showed the same trend toward decreased predictability in all tissues apart from Brain and Testis. These results suggest age-related changes in gene–gene relationships of cell proliferation and mitochondrial genes.Other functions differently affected across tissues included metabolism (mTORC1 signaling and glycolysis; Fig. 2d), development (adipogenesis, myogenesis and epithelial–mesenchymal transition, Fig. 2d; axonogenesis and neuron development, Extended Data Fig. 7), cell growth (mTORC1 signaling and TGFβ signaling; Fig. 2d) and immune response (TNF/NF-κB signaling and IL-6/JAK/STAT3 signaling; Fig. 2d). Of note, gene sets with tissue-specific functions also showed age-related predictability changes in the respective tissue. Among these, we highlight the increased predictability of hypoxia genes in Artery and of coagulation genes in Blood (Fig. 2d) and the decreased predictability of synaptic genes in Brain (Extended Data Fig. 7).To further corroborate those results, we selected in each tissue the 100 genes with the most significant predictability changes with age (lowest regression P values), independent of the direction of that change. We termed these genes ‘high-confidence predictability hits’. For those genes, we also observed considerable differences between tissues, with Blood, Thyroid, Artery and Adipose Tissue showing mostly decreases in predictability with age, and Testis, Breast, Brain and Esophagus showing mostly increases in predictability (Fig. 2e). To exclude that these differences between tissues are driven by differences in sample numbers, we repeated our analysis with similar sample numbers (n = 30) across tissues and age groups (Extended Data Fig. 5 and Supplementary Table 2), which removed differences between tissues due to different sample sizes. This resulted in statistical significance (Extended Data Fig. 5b) and slope directions (Extended Data Fig. 5d) consistent with those obtained with larger sample sizes. Finally, we also repeated our analysis excluding the age group 70–79, as the smaller number of samples in this group might have influenced the predictability values and skewed the regression slope (Supplementary Table 1). However, we once again observed consistent statistical significance (Extended Data Fig. 6b) and preference toward positive or negative slopes (Extended Data Fig. 6d) upon exclusion of this age group. Taken together, these results suggest that predictability is affected by age in some tissues more than others and that the direction of this effect varies between tissues.To confirm these findings in an external dataset unaffected by tissue composition changes, we computed age-related predictability changes in scRNA-seq data from PBMCs of 982 human donors20. We computed age-related slopes across all cell types as well as per cell type (Methods) and compared to the age-related slopes previously obtained in GTEx blood. This comparison revealed an agreement between high-confidence predictability hits found in GTEx blood and those found in T cells (Extended Data Fig. 8), suggesting that the tissue-level age trends may be driven by T cell age trends instead of age-related changes in cell type composition.Regulatory relationship changes with ageNext, we sought to analyze the causes underlying age-associated predictability changes. By definition, low predictability can result from changes in the structure of the regulatory model itself—that is, if the regulatory relationships captured by our model are not met (Fig. 3a, ‘Correlation loss’). However, other factors can modulate predictability. One such factor is the average expression level of the gene at hand: if average gene expression levels are low, the quantification becomes noisier (fewer reads per gene), resulting in lower predictability (Fig. 3a, ‘Low expression’). Another factor is gene expression variance, as a small range of expression values impairs correlation quantification and, thus, predictability (Fig. 3a, ‘Low variance’).Fig. 3: Factors underlying age-related predictability changes.a, Possible factors explaining observations of low predictability and resulting comparisons of observed values (in the original data, x axis) versus predicted values (data reconstructed with our regulatory model, y axis). Top left, low average expression of the target gene results in noisy quantifications. Top right, low variance of the target gene results in almost constant expression values. Correlation approaches thus capture noisy fluctuations around the mean. Bottom left, loss of correlation with the regulatory neighborhood results in a failure of the model to correctly capture the expression pattern of the target gene. b, Expression slope across age (x axis) against variance slope across age (y axis) for background genes (small light gray points). Top 100 significant genes with predictability increase (large blue points) or decrease (large orange points) with age. Age-related variance changes are cropped to remove outliers (see Extended Data Fig. 9a for full data). c, Number of genes with predictability decrease (orange) or increase (blue) in each of the three scenarios tested: ‘Expression change’, ‘Variance change’ or ‘Total’ (total number of predictability hits, for reference).First, we quantified average expression level and variance changes with age. We observed a clear trend toward variance increase with age across all tissues except Artery and Thyroid (Fig. 3b, y axis), in line with earlier reports of increased inter-individual transcriptomic variability with age23,24. The direction of average expression changes with age varied across tissues, with Artery, Brain and Breast showing mostly expression decrease with age, Thyroid showing mostly expression increase with age and the remaining tissues showing similar expression changes in both directions (Fig. 3b, x axis). We found the quantification of average expression and variance changes with age to be influenced by the 70–79 age group in some tissues (Fig. 3b versus Extended Data Fig. 9b). For this reason, the process described below for the identification of genes with predictability, average expression or variance changes is limited to genes showing the same trend when including and excluding the 70–79 age group (Methods).To quantify the contribution of average expression and variance changes to age-associated predictability changes, we compared age-related changes in each of these factors with age-related predictability changes. We applied a cutoff on the absolute values of the expression fold changes or variance slopes with age (0.001) and then counted the number of high-confidence predictability hits (either increase or decrease in predictability) that also showed age-related changes in expression levels (Fig. 3c, ‘Expression changes’) and in expression variance (Fig. 3c, ‘Variance changes’).We observed relatively few genes whose age-associated expression changes were linked to predictability changes (Fig. 3c, ‘Expression change’), suggesting that age-related differences in average expression levels do not explain the predictability changes captured by our approach. Genes with increased variance with age almost always increased their predictability, and, conversely, genes with decreased variance mostly decreased their predictability with age (Fig. 3c, ‘Variance change’), suggesting that age-related changes in gene expression variance influence changes in predictability.Taken together, our analyses suggest that predictability changes captured by our regulatory model are rarely explained by age-related changes in average expression levels. Instead, predictability changes are, in part, explained by changes in inter-individual variability.We next assessed whether correlations to neighbors involved in the same cellular process (MSigDB hallmark) and to neighbors in different cellular processes had a similar impact on the observed changes in predictability. To this end, we grouped the neighbors of each predictability hit according to their functional annotation and summarized age-related correlation changes (Methods) across all neighbors in the same process. This analysis was performed separately for targets with increased and decreased predictability. We then compared the total contribution of within-set relationships (that is, age-related correlation changes between a hit and a neighbor in the same process) and between-set relationships (that is, age-related correlation changes between a hit and a neighbor in a different process). For decreasing predictability hits, we observed that, across all tissues except Adipose, between-set relationships had more impact on predictability changes than within-set relationships (Fig. 4a, bottom). For all tissues except Breast, Esophagus and Testis, between-set relationships had more impact on age-related predictability decreases compared to increases (Fig. 4a, top versus bottom). To take into account differences in the number of neighbors between gene sets, we repeated our analysis averaging the contributions across all regulatory neighbors of the same gene set (weighted average according to the weights in the regulatory model) and found a similar trend toward within-set relationships contributing more strongly toward increasing predictability (Extended Data Fig. 10a). Our results are consistent with a weakening coordination of different functional modules with age, coupled with the strengthening of specific functional responses.Fig. 4: Within-set and between-set regulatory relationships altered with age.a, Contribution of within-set and between-set correlation changes to age-related predictability changes, represented as the weighted sum of age-related correlation changes between genes with increasing (top) or decreasing (bottom) predictability with age and their regulatory neighbors within (full color) or between (transparent) the same gene set. Correlation changes with age are quantified as the slope of the correlation ~ Age regression. These values are weighted with the coefficients of the regulatory model, so that correlation changes in strongly connected regulatory neighbors (higher coefficient) are prioritized. b,e, Subnetwork of the neighborhood of genes with age-related predictability changes. Nodes represent genes, colored by predictability slope in the respective tissue. Connections between nodes represent gene–gene relationships captured by the regulatory model. b, Neighborhood of LAMTOR5, colored by predictability slope with age in Artery – tibial. e, Neighborhood of ITIH2, colored by predicatability slope with age in Blood. c,f, Correlation between expression of genes with predictability changes (LAMTOR5 in c and ITIH2 in f) and regulatory neighbors across age groups. d,g, Comparison of expression of genes LAMTOR5 (d) and ITIH2 (g) in the original data (observed expression, x axis) and reconstructed expression based on the regulatory neighborhood (predicted expression, y axis). The linear regression fit to the trend is shown along with the 95% confidence interval bands. pred., predicted.To further explore the link between predictability and within-set/between-set relationship changes, we focused on the subnetworks created by the 20 strongest (highest absolute value of age slopes) predictability hits in each tissue and their immediate regulatory neighbors. We observed both disconnected subnetworks, where one or very few genes showed strong predictability changes with age (exemplified in Extended Data Fig. 10b), and larger subnetworks connecting several genes with strong predictability changes (exemplified in Fig. 4b,e). We then selected for further analysis two of the obtained subnetworks, showing distinct age-related behaviors and each captured in a different tissue.Between-module gene–gene relationship changes with ageOne of the subnetworks captured by our analysis (Fig. 4b) was centred on mTOR signaling, which is a well-known contributor to age-related phenotypes and is involved in lifespan extension25,26. This subnetwork included the predictability hit LAMTOR5, a gene encoding for a member of the Ragulator complex, located in the lysosomal membrane and involved in mTORC1 activation upon nutrient sensing. In nutrient-rich conditions, mTORC1 is recruited to the lysosomal membrane27, resulting in the promotion of cellular growth and proliferation, coupled to translation and biosynthesis28. In our analysis, LAMTOR5 showed a loss of correlation with most of its regulatory neighbors in Artery – tibial: cytochrome C oxidases COX14 and COX6B1, NADH dehydrogenase complex member NDUFB3, vacuolar ATPase cytosolic (V1) domain subunit ATP6V1F and oligosaccharyl transferase complex subunit OST4 (Fig. 4c). This loss of correlation took place in parallel with a decrease in predictability with age (Fig. 4d).Interestingly, also in Artery – tibial, we observed age-related changes in the predictability of another vacuolar ATPase (v-ATPase) subunit (in the transmembrane domain; Extended Data Fig. 10b), also linked to another Ragulator complex member, LAMTOR1. v-ATPases are located on the membrane of multiple organelles, including the lysosome, where they pump protons into the lumen and acidify it29. The observed relationships between Ragulator subunits (LAMTOR1 and LAMTOR5) and v-ATPase subunits in both the transmembrane and cytosolic domains (ATP6V1F, ATP6V0D1, ATP6V0C, ATP6V1H, ATP6V1D and ATP6V1A) are consistent with the fact that v-ATPases are required for mTORC1 activation by amino acids30. Furthermore, v-ATPases have been described to physically interact with the Ragulator complex through LAMTOR1 (ref. 30), in line with the relationship observed between this gene and multiple of the v-ATPase subunits.In addition, we observed tight gene–gene relationships between Ragulator complex members and electron transport chain members (NDUFB3, COX14 and COX6B1 as direct neighbors of LAMTOR5; UQCRFS1, COX6A1, COX8A, COX17 and COX7A2 as indirect neighbors). These observations are in line with previous work showing that mTORC1 promotes the expression of mitochondrial regulators, such as TFAM and electron transport chain components31. These gene–gene relationships were also altered with age, as exemplified by the age-related loss of correlation between LAMTOR5 and its direct electron transport chain neighbors NDUFB3, COX14 and COX6B1 (Fig. 4c).Within-module gene–gene relationship changes with ageAnother subnetwork captured by our analysis was composed of serum-specific proteins, most of which showed increased predictability with age in Blood (Fig. 4e). This subnetwork included serum albumin ALB; apolipoproteins APOA1, APOA2, APOB and APOH, involved in lipid transport between tissues; members of the coagulation cascade F2 (thrombin), FGA (a subunit of fibrinogen) and F11 (coagulation factor XI); inter-alpha-trypsin inhibitor chains ITIH2 and AMBP; and serine protease HABP2, involved in hyaluronic acid binding. Serum albumin has been reported to have anticoagulant action32, and high-density lipoproteins (HDLs)— composed of APOA1 and, to a lesser extent, APOA2 and APOE—have been shown to modulate platelet reactivity and protect against thrombosis33. Additionally, inter-alpha inhibitor activity is required for the formation of the protective hyaluronan coat surrounding some cell types34. Hyaluronan has been shown to activate thrombin by blocking its inhibitor antithrombin35, supporting the relationships observed between these genes and suggesting a central role of the coagulation cascade in this subnetwork.We observed an age-related increase of ITIH2 predictability (Fig. 4g), in parallel with an increased correlation to its regulatory neighbors (Fig. 4f), and increased variance across individuals. Inter-alpha inhibitors suppress pro-inflammatory responses and have been shown to improve the outcome of ischemic stroke36, which may point toward an adaptive response to an age-related increase in inflammation. This is in line with previous reports of increased blood coagulation potential during aging37, including in centenarians38.Taken together, our analyses uncovered gene–gene relationships that seem to be affected by the aging process. The weakening of gene–gene relationships among mTOR regulators, v-ATPase subunits and electron transport chain members may suggest an age-related de-coupling of these processes, which should be coordinated to ensure appropriate response to cellular growth cues. On the other hand, the strengthening of gene–gene relationships between diverse serum proteins may reflect an age-related increase in the coordination of the coagulation cascade and an adaptive response to increased inflammation.
