Transcriptomics analysis for the identification of potential age-related genes and cells associated with three major urogenital cancers

Age is one of the most important risk factors of the occurrence for tumor patients. The majority of patients with urogenital cancers are the elderly, whose clinical characteristics are greatly affected by age and ageing. Our study aimed to explore age-related genes, cells, and biological changes in three common urogenital cancers via integrative bioinformatics analysis. First, mRNA (count format) and clinical data for bladder cancer, prostate cancer and renal cell carcinoma were downloaded from the Cancer Genome Atlas (TCGA). Through the comparison of clinicopathological characteristics, genes expression and cells infiltration between the old group and the young group, it was found that the clinical characteristics, genes and cells in the tumor microenvironment of different ages were quite different. And 4 key cells, 14 hub genes and some potential pathways were identified and considered as important factors. More importantly, we analyzed the differential landscape of the genes and cells from different perspectives, and confirmed its importance. In conclusion, we identified genes and cell types associated with age-related changes in the tumour microenvironment in urogenital cancer patients. These genes and cell types may play a critical role in the age-associated differences in clinicopathological characteristics among urogenital cancers, thus providing a link between ageing and cancer occurrence. The findings of this study may pave the way for the development of age-tailored approaches to treat cancer and other age-related diseases.


Materials and methods
Data sources and pre-processing. Urogenital cancers are more frequent in elderly patients, and kidney, prostate, bladder tumours are the most common urogenital cancers. The RNA-seq expression profiles (count format) and clinical data of kidney, bladder, and prostate cancers were obtained from the Cancer Genome Atlas (TCGA) using the GDC-client tool. We used data from these three tumours to demonstrate shared characteristics of aging in the TME across urogenital cancers. All statistical analyses were performed using R version 4.0.0 software (R Foundation for Statistical Computing, Vienna, Austria) and GraphPad Prism software (version 8.0; GraphPad Software, La Jolla, CA, USA). After excluding cases with missing age and survival data (4 bladder cancer cases and 3 renal cell carcinoma cases), the patients were divided into young and old groups based on median age. The number of cases included median ages were as follows: bladder cancer, 405 cases (median age, 69 years), prostate cancer, 500 cases (median age, 61 years), and renal cell carcinoma, 884 cases (median age, 60 years).

Influence of age on clinical characteristics of urogenital cancers. Kaplan-Meier analyses were
performed to assess the overall survival of patients in each group and by cancer type. Survival analyses were performed using the log-rank test, performed with GraphPad Prism software. The major clinical characteristics of the patients were also assessed. Identification of age-related genes. Differentially expressed genes (DEGs) between the two groups were identified using edgeR package, and false discovery rate (FDR) were used for adjust to P-value. Genes with |log 2 fold-change (FC)| > 1 and adjust P-value < 0.05 were considered as significant DEGs. Venn diagram analysis (https ://bioin fogp.cnb.csic.es/tools /venny /) was performed to identify the common DEGs among the three urogenital cancer types.

Enrichment analysis of the common DEGs. Gene Ontology (GO) and Kyoto Encyclopaedia of Genes
and Genomes (KEGG) 19,20 enrichment analyses were used to assess the structure, functions, and pathways of the common DEGs. These analyses were conducted using the clusterProfiler R package 21 ; count ≥ 3 and P < 0.05 were used as enrichment cut-offs. The top 10 biological processes, molecular functions, cellular components, and KEGG pathways were identified. The clueGO Cytoscape plugin was used to further analyse biological processes and KEGG pathways.

Protein-protein interaction (PPI) network analysis and hub gene identification. The 188 DEGs
were imported into the STRING database (https ://strin g-db.org/), and PPI analysis was performed using a combined score of ≥ 0.4. There were 139 nodes and 317 edges in the PPI network. The network was reconstructed using Cytoscape, and the cytoHubba plugin was used to identify the top 30 genes in the network using five different algorithms: maximal clique centrality (MCC), density of maximum neighbourhood component (DMNC), maximum neighbourhood component (MNC), edge percolated component (EPC), and the degree method. Genes identified by all of the different algorithms were considered as hub genes.
Relationship between hub genes and clinical characteristics. Univariate Cox regression analysis was used to assess the role of hub genes in patient prognosis, using sex and age as covariates. Receiver operating characteristic (ROC) analysis of hub genes was performed using the pROC R package to assess the differentiation of these genes between the two groups of patients. The relationship between the mRNA levels of the hub genes and patient prognosis in the three urogenital cancers was visualised using the ggplot2 R package. We Co-expression analysis of hub genes. According to the "guilt by association" theory, functionally related genes are often co-expressed 23 . Hence, co-expression analysis can provide insight into the function of poorly annotated genes. The genes co-expressed with the hub genes were identified using PCViz (http://www. pathw aycom mons.org/pcviz ). Then, the chromosomal distribution of hub genes and co-expressed genes was determined using the RIdeogram R package 24 . Metascape (https ://metas cape.org/) was used for further enrichment analysis of hub genes and co-expressed genes 25 .
Changes in cell infiltration patterns. xCell (https ://xcell .ucsf.edu/) predicts the levels of 64 cell types in tissues based on gene expression data 26 . We used the limma R package for differential analysis of cell infiltrates. |log 2 FC| > 1 and P < 0.05 were set as cut-offs, and the results were visualised using a volcano map. Venn diagram analysis was performed to identify cell levels in the three TMEs. Univariate Cox regression analysis was performed to assess the relationship between cell levels and patient prognosis, and the results were visualised using ggplot2. The association between cell infiltration levels and clinical characteristics in bladder cancer patients was analysed using WGCNA R package.
The relationship between key cells and hub genes. The corrplot R package was used to perform a Pearson correlation analysis of hub genes and key cells, and thereby acquire insight into their possible biological links.

Results
Influence of age on the clinical characteristics of urogenital cancers. Kaplan-Meier analysis revealed that in bladder and kidney cancer, elderly patients had a significantly worse prognosis than younger patients; however, no significant difference in prognosis was observed between the two age groups among prostate cancer patients ( Fig. 1A-C). Notably, only 10 of the 500 prostate cancer patients died, likely contributing to the lack of a significant difference in prognosis between the two groups. The proportion of T and M pathological stage patients was higher in the elderly group for all three cancer types, whereas no significant difference was observed in the proportions of young and elderly patients with N stage disease ( Fig. 1D-F). High-grade bladder tumours and high Gleason stage prostate tumours were more frequent in elderly patients than in young patients.  . GO enrichment analysis revealed that the most common biological processes among the DEGs were spinal cord development, striated muscle contraction, antimicrobial humoral response, immune response, and steroid metabolic process (Fig. 3A). Contractile fibre, myofibril, myosin filament, sarcomere, Z disc, I band, myosin complex, blood microparticle, and muscle myosin complex were the most significantly enriched cellular components (Fig. 3B). The most common molecular functions among the DEGs were monooxygenase activity, hormone activity, oxidoreductase activity, substrate-specific channel activity, steroid hydroxylase activity, structural constituent of muscle, oxygen binding, FMN binding, intracellular calcium channel activity, and intracellular chloride channel activity (Fig. 3C). KEGG pathway analysis revealed significant enrichment of the retinol and cytochrome P450-mediated metabolism of xenobiotic compounds (Fig. 3D).
Since there were only a few enriched pathways, we used clueGO to further analyse enriched pathways and biological processes. We identified musculoskeletal movement, sensory perception of mechanical stimuli, oxidoreductase activity, antimicrobial humoral immune response, ventral spinal cord interneuron fate commitment, acute phase response, dorsal spinal cord development, and intracellular cilium activated chloride channel activity as significantly enriched processes (Fig. 3E).
PPI network and age-related hub genes. To investigate the interrelations of the 188 age-related DEGs and identify hub genes, we conducted a PPI analysis using different algorithms (MCC, DMNC, MNC, EPC, and    Elderly patients had a higher mutation frequency compared with young cancer patients. Notably, the mutation frequency of MYH7, MHY6, and MYH7B was considerably higher than that of the other hub genes, exceeding 2.5% (Fig. 4C). Meanwhile, ROC analysis results of these hub genes, namely Area Under Curve (AUC), are shown in Table 2, which showed significant differentiation in gene expression between the two groups. We also assessed the relationship between the 14 hub genes and patient prognosis, and found that the prognostic value of the different hub genes varied according to cancer type (Fig. 4B). Importantly, MYH6 expression levels were significantly associated with age in bladder cancer, while MYBPH exhibited a negative correlation with age (Fig. 4D).
Co-expression analysis of hub genes. Co-expression analysis is often used to infer the function of poorly annotated genes. We identified 79 genes that were co-expressed with the 14 hub genes and constructed a PPI network (Fig. 5A). Here, genes with similar expressions and functions are closer to each other, thus forming a classic bicentric form of 'hairball' . The genes in the upper small hairball structure are mainly involved in the inflammation and immune response of the tissue, while the genes in the lower hairball are mainly involved in the composition and movement of muscle tissue. Chromosome distribution analysis of the co-expressed genes revealed that 21 of the 79 (22.5%) co-expressed genes were located on chromosome 17 (Fig. 5B). Hence, Chromosome 17 may have many age-related genes, which are involved in the regulation of age-related microenvironment changes and cancer initiation. Intriguingly, gene mutations on chromosome 17 have been associated with senile frontotemporal dementia and senile diseases, such as Parkinson's disease 27 . Metascape analysis of the coexpressed genes revealed enrichment in muscle/ motor-related pathways, such as actin, myofibril, and skeletal muscle pathways, as well as in various cancer-related pathways, including interstitial migration and epithelial cell development (Fig. 5C,D). These pathways may be involved in the age-mediated epigenetic regulation of urogenital cancers.
Age-related changes in cell infiltrates. Next, we analysed the levels of 48 cell types in the TME in old and young urogenital cancer patients (Fig. 6A-C). The levels of preadipocytes, CD8 + effector memory T (Tem), CD4 + T, and CD4 + Tem cells in the TME differed significantly between the two age groups in urogenital cancers (Fig. 6D). Preadipocytes levels were higher in the TME of elderly patients than young patients, whereas the levels of the other three immune cell types were lower. Interestingly, preadipocyte levels in the TME were associated with poor patient survival, whereas the levels of CD8 + Tem, CD4 + T, and CD4 + Tem cells were associated with improved survival (Fig. 6E). Consistently, preadipocyte levels in the TME were strongly correlated with the clinical TNM stage in bladder cancer patients (Fig. 6F).

The relationship between age-related cells and hub genes. Correlation analysis of age-related cells
and hub genes revealed a strong association between preadipocyte levels in the TME and the expression levels of ORM1, ORM2 and HP (Fig. 7). Additionally, there was a strong positive correlation between the levels of CD4 + T and CD4 + Tem cells in the TME.

Discussion
The incidence of kidney, prostate, and bladder cancer increases with age, and 82% of newly diagnosed urogenital cancer patients are aged 60 years and older 28 . Although the strong links between ageing and tumorigenesis have become evident, the age-related changes in the TME in urogenital cancer patients remain elusive. Advances in sequencing technologies and bioinformatics tools have enabled unbiased characterisation of age-related changes in the microenvironment in urogenital cancer. In this study, we identified 14 hub genes and 4 cell types enriched in the ageing TME. The genes and cell types identified in this study may mediate age-induced genomic and clinicopathological alterations in urogenital cancers. Age is one of the major factors affecting the occurrence of patients with urogenital cancers. In this study, we found that elderly patients with bladder or renal cell cancer had poorer overall survival and a higher clinical stage than younger patients, consistent with previous findings 29,30 . However, we found no significant difference in prognosis between the two age groups in prostate cancer. It is worth noting that only 10 of 500 prostate cancer patients had died, making it statistically challenging to assess differences in patient prognosis. A study of 259 prostate cancer patients showed that older patients had poorer cancer-specific survival than younger patients 31 .
In this study, we identified 14 hub genes associated with age-related changes in the TME (Table 1), which may contribute to differences in the overall survival and clinicopathological characteristics of patients in different age groups. ORM1 and ORM2 encode endoplasmic reticulum membrane proteins regulating lipid homeostasis and protein quality control 32 . HP expression levels were strongly associated with those of ORM1 and ORM2. Interestingly, all three genes have been implicated in various inflammatory diseases, such as sarcoidosis and chronic obstructive pulmonary disease 33,34 . MYH6, MYH7, and MYH7B are myosin heavy chain proteins, while MYLPF is a myosin light chain protein and MYBPH a myosin-binding protein. Mutations in these proteins have been linked to various muscle-related disorders and heart diseases 35 . CSRP3 is an autophagy-regulating protein essential for degrading muscle-related components and maintaining muscle structure and function 36 . TCAP is a giant elastic protein with kinase activity that extends half the length of a sarcomere and maintains cellular structure 37 . TNNT3 is a fast skeletal muscle troponin T (TnT) initiating muscle contractions 38 . TRIM63 regulates proteasome degradation of cardiac troponin I/TNNI3 and other sarcomeric-associated proteins 39 . INSM1 and NEUROG3 are essential markers for pancreatic neuroendocrine tumours and have been implicated in gastrointestinal tumours, diarrhoea, and malabsorption 40,41 . Hence, these hub genes are mainly involved in muscle structure, metabolism, and inflammation/immune. Interestingly, premature aging and muscle tissue expression disorders were often detected in immunodeficiency mice, and MYH7 and TRIM63 gene expression disorders were also detected 42,43 . In addition, HP and ORM1 were related to immune response, and the expression was decreased in immunodeficiency virus (HIV) infected patients 44,45 . These evidences seem to suggest a link between muscle tissue disorders/immunodeficiency and ageing/cancer. Their relationship with ageing and carcinogenesis required further investigation.
Enrichment analysis of genes co-expressed with the 14 hub genes revealed their role in muscle structure and contraction regulation, ion channel (Cl − , Ca 2+ ) regulation, inflammation, antibacterial humoral immune responses, and metabolism, in addition to various cancer-related pathways. Cancer is widely considered as an age-related disease. Cancer and ageing share various biological processes, including changes in intracellular communication, protein stability, and metabolism, as well as mitochondrial dysfunction 46 . In this study, we found that changes in muscle activity and the immune system were common to both ageing and cancer. Myocardial and skeletal muscle motor capabilities decrease with ageing, and the accumulation of mutations can cause organ dysfunction in older individuals 47 . Furthermore, both cancer and ageing are characterised by disruption of metabolic homeostasis, including enhanced glycolysis and impaired oxidative phosphorylation 48,49 . Changes in the levels of steroids, including adrenocorticosteroid, androgen, and oestrogen, can promote cancer development and progression 50 . Disruption of immune homeostasis is also shared by ageing and cancer, and numerous studies    13,14 . However, more studies are required to elucidate the mechanistic links between cancer and ageing.
In this study, we identified preadipocytes, CD4 + Tem cells, CD4 + T cells, and CD8 + Tem cells as cell markers of age-related changes in the TME. Ageing significantly affects the TME, promoting tumour progression and metastasis 13 . Interestingly, the levels of tumour-infiltrating fibroblasts in elderly prostate patients were found to be particularly low. Preadipocytes levels were strikingly high in the TME of elderly patients. The increase in adipose tissues is one of the hallmarks of ageing, likely reflecting a link between ageing and cancer 46,51 . Preadipocytes are closely related to macrophages and are often dedifferentiated during ageing, switching to a pro-inflammatory, tissue-remodelling, senescent-like state 51 . A recent study has shown that preadipocyte dedifferentiation is regulated by the JAK pathway 52 . The levels of tumour-infiltrating CD4 + Tem cells, CD4 + T cells, and CD8 + Tem cells  www.nature.com/scientificreports/ were lower in our elderly patients than in the young patients. Immune dysfunction is a common feature of ageing; this phenomenon is known as immunosenescence and has been linked to impaired immune surveillance and increased risk of cancer development 14 . Immunosenescence is associated with extensive alterations in the immune system, and CD4 + and CD8 + T cells are particularly susceptible to immunosenescence 53 . Meanwhile, immunodeficiency and immunosenescence are both manifested as decreased or absent immune cell content. CD4+ T cells and CD8+ T cells were significantly reduced in common immunodeficiency diseases 54 . The incidence of multiple cancers is increased in patients with acquired immune deficiency, which also reflects the influence of immune deficiency on the promotion of tumour incidence to some extent 55 . In addition, patients with immunodeficiency are more likely to develop flora disorder and urinary tract infections, which may also increase the incidence of corresponding tumours 56,57 . Therefore, the infiltration level of immune cells is related to the tumorigenesis, and modulating the infiltration of these four cells in TME may provide a novel platform to treat cancer. This study is not without limitations. First of all, this study was only analysed in urogenital cancers and may not be applicable to all tumours. Secondly, matched adjacent-normal expression data is very rarely available in tumour datasets so this is probably just a fact of life for this kind of study but it does mean that it is not possible to correct for general effects of aging unrelated to cancer. It is notable for example that most of the DEGs identified were not obviously related to a possible age-related tumour response e.g. immunosenesence, but rather to pathways one might expect to be affected by normal aging processes e.g. muscle. Additionally, muscle tissue changes and immunosenesence, as common features of elderly patients, may be associated with tumorigenesis. However, whether muscle tissue changes and immunosenesence are a cause or an effect of tumorigenesis need a further study.
In conclusion, we identified 14 genes and 4 cell types associated with age-related changes in the TME in urogenital cancer patients. These genes and cell types may play an important role in the age-related differences in the clinicopathological characteristics of urogenital cancer patients, thus linking ageing with cancer development