Predicting Prostate Adenocarcinoma Patients’ Survival and Immune Signature: A Novel Risk Model Based on Telomere-Related Genes

doi:10.21203/rs.3.rs-3428019/v1

Download PDF

Research Article

Predicting Prostate Adenocarcinoma Patients’ Survival and Immune Signature: A Novel Risk Model Based on Telomere-Related Genes

https://doi.org/10.21203/rs.3.rs-3428019/v1

This work is licensed under a CC BY 4.0 License

You are reading this latest preprint version

Background

Telomere-related genes (TRGs) play an essential role in the carcinogenesis and progression of prostate adenocarcinoma(PRAD). However, the prognostic value of TRGs remains unclear in PRAD.

Methods

We conducted a study using The Cancer Genome Atlas-Prostate Adenocarcinoma (TCGA-PRAD) dataset as the training group and the Memorial Sloan-Kettering Cancer Center (MSKCC) and Gene Expression Omnibus (GEO) datasets as the validation group. We developed a risk model and a nomogram to predict survival rates in patients with PRAD. The expression of model genes and their possible regulatory mechanisms were then analyzed. Furthermore, we explored the relationship between the risk model and immune cell infiltration, chemotherapy drug sensitivity, and specific signaling pathways using the CIBERSORT algorithm, the Genomics of Drug Sensitivity in Cancer (GDSC) database, and motif enrichment.

Results

The developed risk model was based on seven key TRGs (HELLS, TOP3A, SRC, LARP7, BUB3, THRSP, and GTF2H4). Moreover, this model was an independent prognostic factor for patients with PRAD and was significantly associated with T and N stages. Among seven TRGs, TOP3A and BUB3 were not only significantly positively correlated with the clinical T and N phases of PRAD, but also expression was significantly increased in PRAD tissues compared to adjacent normal tissues. The model was a good predictor of 1-, 3-, and 5-year survival, and patients in the high-risk group had significantly shorter overall survival than those in the low-risk group. The integrated nomogram can be a good predictor of 3- and 5-year survival in patients with PRAD. Finally, compared to the low-risk group, the high-risk group had a higher response to chemotherapy and immunosuppression, which provided potential guidance to treatment options for patients in the high-risk group.

Conclusion

In summary, a new risk model based on TRGs was successfully developed in PRAD. This risk model is valuable for guiding the selection of immunotherapy and chemotherapy in the clinical treatment of patients with PRAD.

Telomere

risk model

immune microenvironment

prostate adenocarcinoma

immune checkpoint inhibitors

Prostate adenocarcinoma (PRAD) is a prevalent male malignancy, displaying diverse clinical presentations, treatment responses, and recurrence patterns due to its molecular heterogeneity [1, 2]. Precise treatment plans are challenging to formulate for individual patients. Therefore, there is a crucial need to establish reliable prognostic models [3]. Existing models, such as the Gleason scoring system, ISUP prostate cancer grading system, and AJCC TNM staging system, integrate patient characteristics but have limitations [4–6]. Advancements in sequencing technology offer an opportunity to identify critical genes associated with PRAD prognosis, enabling more accurate risk assessment and novel molecular subtyping methods surpassing current staging approaches [7].

Telomeres, which are unique structures located at the ends of chromosomes, consist of highly repetitive DNA sequences and play a crucial role in maintaining the integrity of the genome. In normal cells, telomeres gradually shorten with each cell division, leading to cell cycle arrest or apoptosis [8]. However, tumor cells employ unique mechanisms like telomerase activation and telomerase-independent alternative lengthening telomeres for maintenance [9, 10]. Studies have demonstrated that telomere length changes and maintenance have been linked to tumor survival, proliferation, metastasis, and poor clinical outcomes in various cancers, including prostate cancer [11–13]. Prostate cancer is characterized by significant telomere shortening associated with genomic instability, chromosomal abnormalities, and oncogene activation. Heaphy et al. found that men with more pronounced intercellular telomere length changes in PRAD tissue or shorter telomeres in cancer-associated stromal cells had a higher risk of disease progression or death [14, 15]. Recent studies also have examined the effects of core components of telomerase, including the telomerase RNA gene TERC and the telomerase reverse transcriptase, on PRAD. Baena et al. observed that increased TERC expression was associated with a poorer prognosis in patients with PRAD [16]. Furthermore, Poos et al. found that the expression levels of PITX1, a transcription factor that regulates TERC, increased with increasing levels of Ki67 in PRAD, indicating a poorer clinical outcome [17, 18].

However, alterations in telomeres, present in both precancerous lesions and normal cells, impact cancer-specific telomere length. Reduced expression of telomerase-associated proteins hinders their suitability as ideal prognostic markers. Studies have demonstrated a significant association between telomere-related genes (TRGs) and PRAD prognosis [19, 20]. We investigated this correlation by analyzing extensive gene expression profiles in patients with PRAD and developed a predictive risk model based on TRGs. This model aims to establish a precise and personalized strategy for managing prostate cancer. Additionally, we assessed the clinical significance and immunotherapeutic capacity of the TRG-based risk model and essential genes to facilitate the development of more accurate and individualized treatment approaches for patients with PRAD.

2.1 Data acquisition and preprocessing

Clinical information and genomic data from 502 PRAD tumor samples and 52 standard prostate samples were downloaded from the Cancer Genome Atlas (TCGA) database (https://tcga-data.nci.nih.gov/tcga/) [18]. Transcriptomic data and survival information for 140 patients with PRAD were downloaded from the Memorial Sloan Kettering Cancer Centre project (MSKCC) using cBioPortal (http://www.cbioportal.org/). We obtained the GSE70768 dataset from the NCBI GEO (http://www.ncbi.nlm.nih.gov/geo) public database [21]. This dataset contained data from 111 patients with PRAD, including their complete expression profiles and survival information. The platform number for the annotations in the matrix file was GPL10558. Traditional screening methods were used for genetic screening [22, 23]. We obtained 1,323 TRGs (relevance score > 1) from the GeneCards database (https://www.genecards.org). All of the data were available for free online. A schematic overview of the study design with respect to the above data is shown in Fig. 1.

2.2 Genes function enrichment analysis

Gene Set Enrichment Analysis (GSEA) was performed using the “clusterProfiler” package in R [24]. Gene ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analyses were used for the evaluation of relevant category functions. GO terms and KEGG pathways with p < 0.05 were considered significantly enriched. In addition, we performed functional annotation of key genomes through the Metascape database. Specific genes then underwent GO and KEGG pathway analyses. Min overlap ≥ 3 and p ≤ 0.01 were considered statistically significant.

2.3 Model development and prognosis

We developed a prognostic model comprising seven genes by employing LASSO regression and multivariate Cox regression of prognosis-associated TRGs using the ‘glmnet’ package in R. After including the expression values of each gene, a risk score formula was constructed for each patient. The formula was weighted with the estimated regression coefficients of each gene using LASSO regression analysis. Based on the median risk score value, we categorized the patients into low- and high-risk groups [25]. Kaplan-Meier and log-rank statistical methods were used to compare the difference in survival between the two groups. Additionally, we determined the prognostic accuracy of the risk score model using LASSO regression analysis, stratified analysis, and receiver operating characteristic (ROC) analysis.

2.4 Nomogram construction

We established a prognostic nomogram using regression analysis. To express the interrelationship between the variables in the risk model, we constructed line segments with scales on the same plane based on the level of gene expression and clinical symptoms at a specific ratio. By constructing a regression model, we assigned a score to each value level of each influencing factor based on their contribution to the outcome variable (as indicated by the magnitude of the regression coefficient). Subsequently, we added up the individual scores to obtain a total score, which allowed us to determine the predictive value.

2.5 Drug sensitivity evaluation

Based on the Genomics Database of Cancer Drug Sensitivity In Cancer (GDSC, https://www.cancerrxgene.org/), the largest pharmacogenomic database available, we predicted the chemotherapeutic sensitivity of each tumor sample using the “pRRophetic” R package. We used regression to derive IC₅₀ values for each unique chemotherapeutic drug. Additionally, we conducted 10 cross-validations on the GDSC training set to assess the accuracy of regression and prediction. All parameters, including “combats,” were set to their default levels to eliminate batch effects and calculate the mean of duplicate gene expression.

2.6 Immune cell infiltration analysis

We analyzed RNA-seq data from various subgroups of patients with PRAD using the CIBERSORT algorithm and determined the proportions of 22 immune infiltrating cells. A p-value < 0.05 was considered significant for gene expression and immune cell content in the Pearson correlation analyses [26].

2.7 Immunohistochemical (IHC) staining

Clinical tissues from patients with PRAD were fixed and sectioned. The sections were heated, dewaxed, immersed in 100% xylene, and dehydrated using a graded alcohol series. Hydrogen peroxide (0.3%) was added. And the sections were then boiled in citrate buffer (pH 6.0) for 30 minutes using an autoclave. Sealing was achieved by adding 10% bovine serum. Subsequently, 100 µL of primary antibody, diluted in 1% bovine serum albumin, was added to each section. The sections were then incubated overnight in a humidified chamber at 4°C with the following primary antibodies: HELLS antibodies (1:300, PA5-64099, Thermo Fisher), TOP3A antibodies (1:200, PA5-121538, Thermo Fisher), SRC antibodies (1:200, ab109381, abcam), LARP7 antibodies (1:200, 17067-1-AP, proteintach), THRSP antibodies (1:400, 13054-1-AP, proteintach), GTF2H4 antibodies (1:200, PA5-103346, Thermo Fisher), and BUB3 (1:200, ab133699, abcam ). The corresponding secondary antibodies were added to tissue sections. The sections were observed under a microscope after 3,3′-diaminobenzidine development and hematoxylin staining.

2.8 Gene set variance analysis(GSVA)

To identify differentially expressed genes between the high- and low-risk groups, we performed an enrichment analysis of the gene set using the PRAD expression profile (GSEA; http://www.broadinstitute.org/gsea) [27]. The “c2.cp. Kegg.v7.4.entrez.gmt” gene set was downloaded from the MSigDB database and applied filters to include gene sets with sizes ranging from 15 to 500 genes. GSEA was then performed using the “clusterProfiler” package in R. GSVA was performed using the GSVA package. The “limma” algorithm was used to determine pathway differences between the high- and low-risk groups. By intersecting the results obtained from the transcriptome-based GSEA and proteome-based GSVA, we identified differential pathways between the high- and low-risk groups.

2.9 Regulatory network analysis

Transcription factors are proteins that bind to specific DNA sequences. The R package “RcisTarget” was used to predict the transcription factors associated with the model genes. The area under the curve (AUC) was calculated for each motif-motif set pair to evaluate the overexpression of each motif within a gene set. Subsequently, we calculated the normalized enrichment score (NES) for each motif based on the distribution of AUCs for all motifs in the genome.

2.10 Statistical analysis

All statistical analyses were performed using the R software (version 4.0.2). Survival curves were generated using the Kaplan-Meier method, and the differences between the curves were assessed using the log-rank test. Furthermore, a multivariate analysis was conducted using the Cox proportional hazards model. Both single- and multi-factor Cox analyses were performed by using the R survival package, and forest plots were generated to visualize the results. Two-group comparisons were assessed using the Wilcoxon test function, whereas multi-group comparisons were assessed using the Kruskal-Wallis test function. Correlation analysis was performed using the Pearson method. The ssGSEA scores for the samples were calculated using the ssGSEA function in the “GSVA” package. P < 0.05 was considered significant.

3.1 Screening and functional enrichment analysis of TRGs associated with PRAD prognosis

The univariate Cox regression analysis revealed 136 TRGs associated with PRAD prognosis (p < 0.05) (Supplementary file 1). We also performed a functional enrichment analysis of these TRGs. The GO enrichment analysis revealed that the majority of these genes were related to telomere regulation, chromosomal area, DNA-dependent ATPase activity, and other related pathways (Fig. 2A). In addition, the KEGG enrichment analysis revealed that these genes were considerably enriched in pathways such as homologous recombination and cell cycle (Fig. 2B). Furthermore, the Metascape enrichment analysis revealed that these genes were predominantly involved in biological activities such as DNA metabolism, catalytic activity, DNA activity, and DNA metabolism (Fig. 2C). To investigate the potential relationships between the TRGs, we performed a Protein-protein interaction (PPI) network analysis using the Cytoscape software (Supplementary Figure S1). The results imply that the interactions between the TRGs are complex, closely related, and closely associated with DNA-related biological processes, suggesting that 136 TRGs may play an essential role in the progression and prognosis of PRAD.

3.2 Construction and validation of the TRG-based risk model

We used 136 TRGs as candidate genes for modeling and selected TCGA-PRAD cohorts with complete expression profiles and survival data as the training set. Detailed clinical characteristics of the 500 patients with PRAD in the TCGA-PRAD cohort are shown in supplementary table 1. The MSKCC and GSE70768 patient cohorts were used as external validation sets. Additionally, clinical data were collected from patients with PRAD, and the final candidate genes were screened using the LASSO regression algorithm (Figs. 3A–3C). The formula for the risk model was: Risk Score = (HELLS × 0.00716375741617321) + (TOP3A × 0.0301970796315794) + (SRC × 0.0327209383443109) + (LARP7 × 0.0408035537429572) + (BUB3 × 0.11506165146503) + (THRSP × 0.158045396171819) + (GTF2H4 × 0.426874220611804). The patients with PRAD were divided into high- and low-risk groups based on the median risk score value. In the training cohort, the high-risk group had significantly lower OS than the low-risk group (Fig. 3D). Similarly, we validated the stability of the risk model by calculating survival differences between patients with PRAD in the GEO and MSKCC cohorts (Figs. 3E–3F). To further assess the performance of the risk model, we generated time-dependent ROC curves. The risk model’s 1-, 3-, and 5-year survival AUC values were 0.99, 0.88, and 0.81, respectively (Fig. 3G). We also performed ROC curve analysis on external validation sets to validate the accuracy of the model. The results indicated that the risk model can be useful in predicting prognosis (Figs. 3H–3I).

3.3 Differential expression and regulatory networks of the model genes

We performed a correlation analysis between clinical data and model genes and found that HELLS, TOP3A, BUB3, and GTF2H4 were significantly positively correlated with clinical T and N stages of PRAD (Fig. 4A–4B). To evaluate the expression levels of the seven critical genes used to construct the risk model in PRAD and normal tissues, we downloaded RNA sequencing data from TCGA and GEO databases for differential analysis. The results of TCGA datasets revealed that the mRNA expression of HELLS, TOP3A, LARP7, BUB3 and THRSP was higher in PRAD compared to normal. However, only the expression of BUB3 was found to be higher in PRAD (Fig. 4C–4D). To further investigate the expression of the model genes in clinical samples, we collected tissues from patients with PRAD and performed IHC staining on cancerous and adjacent normal tissues. The IHC representative results revealed that TOP3A, SRC and BUB3 exhibited stronger staining in PRAD tissues compared to adjacent normal tissues (Fig. 4E-4F). The staining intensity of other genes did not show significant differences between cancerous and adjacent normal tissues.

Next, we investigated the specific regulatory mechanisms of the model genes in the risk model. We identified several transcription factors and common regulatory mechanisms associated with these genes. Using cumulative recovery curves, we conducted an enrichment analysis of these transcription factors (Fig. 4G). The results show that the motif annotation with the highest NES was cisbp_M0497. This motif was overrepresented in the BUB3, GTF2H4, LARP7, and TOP3A genes (Supplementary Figure S2). Further details can be found in the specification document (Supplementary file 2).

3.4 Exploration of the prognostic value and clinical relevance of the risk model

We validated the prognostic value of the risk model in each patient with PRAD. The univariate and multivariate analyses revealed that the risk score was also an independent prognostic factor for PRAD (Fig. 5A). Additionally, the AUC value of the risk score was greater than that of the conventional prognostic scoring system (AUC: risk score = 0.78), further indicating that the risk model has superior prognostic value over other clinical characteristics (Fig. 5B). To investigate further the correlation between the risk model and clinical indicators, we presented the risk score values of all samples as box plots (Fig. 5C). We discovered that the risk score values were significantly correlated between T and N clinical groups but not between age and survival status groups.

3.5 Construction and validation of the risk score-related nomogram model

To further improve the predictive power and applicability of the risk model, a predictive nomogram (Fig. 6A) was developed by integrating the age, T category, N category, and risk score of patients with PRAD. The nomogram can systematically predict the 3- and 5-year overall survival (OS) of patients with PRAD. The calibration plots demonstrated that the nomogram model had high accuracy and congruence between actual and predicted results (Fig. 6B). The 3- and 5-year survival AUC values of the predictive nomogram were 0.78 and 0.70, respectively (Fig. 6C). Comparing the 3-and 5-year survival AUC values of the nomogram with other clinical factors, we found that the nomogram had a substantially higher AUC value, indicating that our nomogram has good discrimination (Figs. 6D–6E). The decision curve analysis indicated that the nomogram has greater clinical applicability (Supplementary Figure S3). These findings suggest that the nomogram model is suitable for predicting the 3-year and 5-year survival of patients with PRAD.

3.6 Exploration of signaling mechanisms related to the risk model

To elucidate the potential impact of the risk score on the molecular progression of PRAD, we analyzed the distinct signaling pathways implicated in the high- and low-risk groups. The GSEA analysis revealed significant enrichment in several related pathways, such as the GO enrichment pathways of DNA recombination and endoplasmic reticulum mannose trimming (Fig. 7A), as well as the KEGG enrichment pathways of fundamental basal transcription factors and base excision repair (Fig. 7B). Additionally, the GSVA analysis indicated that differential pathways in both groups were predominantly enriched in E2F targets, DNA repair, G2M checkpoint, and other signaling pathways, suggesting that alterations in these signaling pathways in the high- and low-risk groups affected prognosis (Fig. 7C). Our findings suggest that these signaling pathways are potential therapeutic targets for developing more effective treatments for high-risk patients with PRAD.

3.7 Analysis of the tumor immune microenvironment in the high- and low-risk groups

The tumor immune microenvironment considerably impacts tumor diagnosis, survival outcomes, and clinical treatment sensitivity[28]. Consequently, we analyzed the relationship between risk score and tumor immune infiltration (Supplementary Figure S4). The immune cell content of each patient with PRAD is shown in Fig. 8A. The results revealed that the low-risk group samples had significantly higher levels of T-cell CD4 memory quiescent, dendritic cell (DC) quiescent, and mast cell quiescent and significantly lower levels of regulatory T cells (Tregs) compared to the high-risk group. Additionally, we performed a Pearson correlation analysis on 22 immune cells, revealing multiple pairs of significant correlations between several immune cells and risk scores (Fig. 8B). Furthermore, risk scores were positively correlated with Tregs and macrophage M2 but negatively correlated with mast cell resting (Fig. 8C). We also investigated the correlation between the seven model genes and the immune cells (Fig. 8D). The expression levels of the model genes showed significant correlations with various immune cells. For example, BUB3 was positively correlated with T cell CD4 memory quiescence, macrophage M1, and DC quiescence but negatively correlated with natural killer cell quiescence, DC activated, mast cell activated, and neutrophils. Additionally, we performed a correlation analysis between these model genes and immune ssGSEA scores (the set of immune-related genes with correlation scores > 7 in the GeneCards database) (Fig. 8E) and found a positive correlation between HELLS expression and immune scores.

3.8 Correlation between risk score and immune-related genes

Immune checkpoint inhibitors (ICIs) have been used in recent years for the treatment of several cancers. We explored the differences in the effectiveness of treatment with ICIs between the high- and low-risk groups to inform the clinical management of PRAD. We extracted several sets of immune-related genes, including immunomodulators, from the TISIDB database to further investigate the association between risk scores and immune checkpoints, immunomodulators, chemokines, and cell receptors (Fig. 9A–9D)(Supplementary Figure S5). We found that nearly all immune checkpoint genes (CTLA-4, TGFB1, CD274, LAG-3, TIGIT, and BTLA) were considerably expressed in the high-risk group relative to the low-risk group, which may partly explain the favorable prognosis associated with the low-risk group and suggest that patients with high risks may benefit from treatment with ICIs. In addition, we analyzed the effect of immunotherapy in the high- and low-risk groups using the tumor immune dysfunction and exclusion (TIDE) algorithm. Significant differences in CD274, microsatellite instability (MSI), and TIDE scores were observed between the high- and low-risk groups (Fig. 9E–9G). The higher expression of CD274 and TIDE scores and lower MSI scores in the high-risk group further indicated a poor prognosis. However, although the effectiveness of ICIs is unknown, high-risk patients may benefit more than low-risk patients.

3.9 Sensitivity of chemotherapeutic agents in PRAD

Several studies have demonstrated that conventional endocrine therapy does not improve OS in patients with PRAD. However, early administration of neoadjuvant chemotherapy combined with endocrine therapy can significantly improve survival [29]. For example, treatment regimens with docetaxel combined with endocrine therapy have yielded definitive clinical results and are recommended as part of the standard of care for men who are physically strong and beginning long-term hormone therapy. Using the “pRRophetic” R package, we predicted the chemotherapy sensitivity of each tumor sample, offering additional drug options and potential regimens for neoadjuvant chemotherapy. The results revealed a significant correlation between risk scores and drug sensitivity to gemcitabine, doxorubicin, axitinib, bleomycin, cytarabine, doxorubicin, and vincristine but not to cisplatin (Fig. 10A).

Telomeres play a crucial role in the development of PRAD. In the early stages of PRAD, there is a noticeable shortening in the length of telomeres. This shortening of telomeres leads to genomic instability, resulting in the inactivation of tumor-suppressor genes and the production of oncogenes involved in the initiation and progression of PRAD [30]. The length of telomeres in PRAD tissues can be a potential prognostic marker [31]. PRAD cells activate telomerase to maintain the shortened telomere length at a level that supports unlimited replication. In contrast, normal prostate cells have undetectable levels of telomerase activity. Telomere stability is closely related to the prognosis of PRAD. In castration-resistant prostate cancer, telomere dysfunction leads to more invasive cancer cells [32]. Studies have demonstrated a significant association between certain TRGs affecting telomere alterations and PRAD prognosis. To investigate the predictive value of TRGs, we constructed a risk model based on TRGs using a public database. Subsequently, we validated the predictive ability and clinical relevance of this risk model. Through our study, we aim to offer novel insights that may contribute to identifying potential treatment options for PRAD.

We used LASSO regression and Cox survival analyses to construct a risk model based on seven TRGs: HELLS, SRC, LARP7, BUB3, THRSP, and GTF2H4. We found that HELLS, TOP3A, BUB3, and GTF2H4 showed significant correlations with PRAD’s clinical T and N stages. Furthermore, immunohistochemical analysis of clinical PRAD tissues demonstrated elevated staining levels of TOP3A, SRC, and BUB3 compared to adjacent normal tissues. These results suggest that these three significantly different genes may be potential therapeutic targets. Further investigations are required to delve. Notably, the risk model exhibited robust predictive capabilities, as evidenced by OS and ROC curve analyses. The risk model exhibited powerful predictive capabilities, as evidenced by overall survival OS and ROC curve analyses. Additionally, it proves to be an independent prognostic factor for patients with PRAD and has shown superiority compared to traditional predictive scoring systems. We found a significant correlation between the risk model and the T/N tumor staging. Through GSEA and GSVA, we identified pathways enriched in both the high- and low-risk groups, including DNA recombination, endoplasmic reticulum mannose modification, fundamental transcription factors, and base excision repair pathways. However, there were differences in specific ways, such as DNA repair and the G2M checkpoint. These differential signaling pathways may be potential therapeutic targets that could lead to the development of more effective treatments for high-risk patients with PRAD. Moreover, our integrated nomogram demonstrated high accuracy in predicting the 3- and 5-year survival rates of patients with PRAD. Although risk models have shown excellent predictive power in patients with PRAD, comparing our risk model with other published prognostic models is essential to determine its superiority. In clinical practice, surgeons typically assess prognosis and guide treatment based on the TNM stage, biopsy Gleason score, and pretreatment prostate-specific antigen levels. Therefore, to develop a more accurate nomogram, we need to collect additional samples and incorporate Gleason score staging and prostate-specific antigen levels into our model [33].

Typically, a local immunosuppressive zone develops within the PRAD tumor, severely suppressing immune cell function [34]. However, DCs, the most potent antigen-presenting cells, can overcome the immunosuppressive effect and directly present the ingested antigen to T lymphocytes. DCs induce cell-mediated immune responses and exhibit antitumor effects on cytotoxic T cells [35, 36].Increased CD4 + or CD8 + cells correlate with poorer biochemical recurrence and PRAD-specific survival [37]. However, memory CD4 + T cells can exert antitumor immune effects via multiple mechanisms, including direct regulation of granzyme and perforin expression and indirect regulation of CD 8 + memory T cell production [38]. As intrinsic immune cells, most cells play a regulatory role in PRAD. Increased numbers of tumor-infiltrating mast cells have been associated with the progression and prognosis of PRAD [2]. However, the specific role of mast cells remains controversial. Mast cells can activate metalloproteinase 9, which is crucial in early PRAD development and can promote tumor vascular production [39]. The samples in the low-risk group had significantly higher levels of T cell CD4 memory resting, DCs resting, and mast cells resting compared to the high-risk group. Additionally, the low-risk group had significantly lower levels of Tregs compared to the high-risk group. These results suggest that the high-risk group had a poorer prognosis than the low-risk group. Tregs play an essential role in maintaining immune tolerance and homeostasis. Tregs can prevent the activation of CD8 + T cells by inhibiting DC expansion and immunogenicity, leading to poor immunotherapeutic outcomes in PRAD [40]. Tregs, in addition to their exert antitumor immune effects, also induce the development of macrophages toward an M2-like tumor-promoting phenotype, thereby enhancing the proliferation of PRAD cells. The differences in immune cell infiltration between the low-risk and high-risk groups further indicate a more favorable prognosis for the low-risk group. Tumor-associated macrophages are a dynamic multifunctional cell population with an M2-like phenotype that mediates immunosuppression through the production of interleukin 10, transforming growth factor-β, prostaglandin E2, and matrix metalloproteinase-7. The M2-like phenotype is associated with aggressive tumors characterized by higher Gleason scores, increased metastasis rates, and lower cancer-specific survival rates [41]. In the present study, risk scores showed a significant positive correlation with Tregs and M2 macrophages while displaying a significant negative correlation with mast cell quiescence. These findings suggest that higher risk scores are associated with unfavorable prognosis in patients with PRAD, and this is consistent with the results of our training and validation sets.

Peng et al. demonstrated that lower TIDE scores were associated with a more favorable response to immunotherapy [42]. Tumors with elevated MSI scores tend to have a better prognosis, which may be related to the immunological properties of the tumor [43–45]. Elevated CD274 expression may indicate that tumor cells play a role in evading immune surveillance and attack [46]. In our study, the high-risk group exhibited higher CD274 expression levels and higher TIDE scores while demonstrating lower expression levels of MSI. These findings further support the notion of a poorer prognosis for the high-risk group. In recent years, ICIs have shown remarkable effectiveness in treating many cancers. ICIs block the loci of PD-1, PD-L1, and CTLA-4, which are involved in inhibitory signaling pathways in T cells. The use of anti-PD-L1 monoclonal antibodies has shown excellent effectiveness in renal and uroepithelial cell carcinoma[47, 48]. Upon analying immunomodulators and immune checkpoints in both the high-risk and low-risk cohorts, noteworthy disparities were identified in the expression levels of numerous immune-related genes. Notably, most immune genes showed higher levels in the high-risk group compared to the low-risk group. Thus, high-risk patients may benefit more from ICIs than low-risk patients, but the effectiveness is unclear. Studies have shown that conventional antitumor chemotherapeutic agents, in addition to directly inhibiting tumor growth, can induce immunogenic cell death (ICD) [49, 50]. In the present study, a significant correlation was found between risk scores and drug sensitivity to gemcitabine, doxorubicin, axitinib, bleomycin, cytarabine, doxorubicin, and vincristine, with generally lower IC50 values in the high-risk group than in the low-risk group. The developed risk model shows promise as a predictor of chemotherapeutic response, potentially aiding in the identification of the most appropriate chemotherapeutic regimen for each patient with PRAD. And we believe that ICD-inducing chemotherapeutic agents combined with ICIs will be a promising treatment option for patients with advanced PRAD who have high-risk scores.

Although this is the first risk model based on TRGs in PRAD, it still has some notable limitations. First, the accuracy and applicability of the risk model need to be validated using additional clinical data. We discussed the differences in model gene expression levels between primary prostate cancer tissues and normal paracancerous tissues and compared them with the mRNA expression of model genes in the TCGA database. However, we needed more clinical data and samples to validate the accuracy of the results. The infiltration levels of T cell CD4 memory resting, DCs resting, and mast cells resting were significantly higher in the low-risk group compared to the high-risk group. Additionally, the low-risk group exhibited negative correlations with Tregs and macrophage M2 expression levels. In addition to providing a clinical method for determining the prognosis of patients with PRAD, these findings further validate the exceptional prognosis of the low-risk group and the impact of the immune environment. Although the high-risk group had a poor prognosis, there was high sensitivity to the chemotherapeutic agents and ICIs, suggesting that chemotherapy combined with ICSs may be an effective treatment option for patients in the high-risk group.

In summary, we constructed a risk model based on seven TRGs to predict the prognosis of patients with PRAD. This model showed commendable performance in independently evaluating the prognosis of OS at 1, 3, and 5 years. In addition, we further determined the relationship between the immune environment, ICIs, chemotherapy sensitivity and risk models. Compared to the low-risk group, the high-risk group had a higher response to chemotherapy and immunosuppression, which provided potential guidance to treatment options for patients in the high-risk group.

AUC：Area under the curve

DC：Dendritic cell

GDSC：Genomics of Drug Sensitivity in Cancer

GEO：Gene Expression Omnibus

GO：Gene ontology

GSEA：Gene Set Enrichment Analysis

GSVA：Gene set variance analysis

ICD：Immunogenic cell death

ICIs：Immune checkpoint inhibitors

IHC：Immunohistochemical

KEGG：Kyoto Encyclopedia of Genes and Genome

MSI：Microsatellite instability

MSKCC：Memorial Sloan-Kettering Cancer Center

NES：Normalized enrichment score

OS：Overall survival

PPI：Protein-protein interaction

PRAD：Prostate adenocarcinoma

ROC：Receiver operating characteristic

TCGA-PRAD：Cancer Genome Atlas-Prostate Adenocarcinoma

TIDE：Tumor immune dysfunction and exclusion

TRGs：Telomere-related genes

Acknowledgments

We would like to acknowledge the Guangzhou Laboratory of Regenerative Medicine and Health, Guangdong Province，for granting us access to their state-of-the-art laboratory facilities. Without their generous provision of equipment and resources, this study would not have been possible.

Founding

This research was supported by the National Natural Science Foundation of China (Grant No. [82072808]). The authors gratefully acknowledge the financial support provided by the National Natural Science Foundation of China. The funding from the National Natural Science Foundation of China played a crucial role in conducting this study.

Competing Interests

The authors have no relevant financial or non-financial interests to disclose.

Author Contribution

All authors contributed to the study conception and design. Conceptualization, methodology, and writing, JZ and JC; data curation and formal analysis, JZ; data processing, JZ ,and JC; investigation, JC and YL; visualization, JZ and HL; validation, HL; conceptualization, reviewing, editing, and supervision, XJ. and WD. Specifically, JZ and JC contributed equally to this work. All authors have read and approved the final manuscript.

Data availability statement

The raw data used were obtained from publicly available databases in this study. The specific databases and their accession numbers or URLs are provided in this paper's 'Methods' section. The data were accessed on 2023-3 and were available until then. Researchers can access the same data by visiting the respective databases and following the instructions provided by the database administrators. Patients with primary prostate cancer (n = 3) and paracancerous tissue samples from the same patients (n = 3) were obtained from the Department of Urology, the Fifth Affiliated Hospital of Guangzhou Medical University. Further inquiries can be directed to the corresponding authors.

Ethics statement

The Ethics Committee of the Fifth Affiliated Hospital of Guangzhou Medical University has approved the analysis of diagnostic residual tissues for research purposes (Approval No. KY01-2020-03-01). Patients provided written informed consent before agreeing to provide biospecimens for research use at this institute.

Consent to participate

Informed consent was obtained from all individual participants included in the study.

Consent to publish

The authors affirm that human research participants provided informed consent for publication of the images in Figure(3) 3e and 3f.

Siegel RL, Miller KD, Fuchs HE, Jemal A. Cancer statistics, 2022. CA Cancer J Clin. 2022;72(1):7-33. https://doi.org/10.3322/caac.21708
Sridaran D, Chouhan S, Mahajan K, et al. Inhibiting ACK1-mediated phosphorylation of C-terminal Src kinase counteracts prostate cancer immune checkpoint blockade resistance. Nat Commun. 2022;13(1):6929. https://doi.org/10.1038/s41467-022-34724-5
Tzelepi V, Grypari IM, Logotheti S, et al. Contemporary Grading of Prostate Cancer: The Impact of Grading Criteria and the Significance of the Amount of Intraductal Carcinoma. Cancers (Basel). 2021;13(21). https://doi.org/10.3390/cancers13215454
Buyyounouski MK, Choyke PL, Mckenney JK, et al. Prostate cancer - major changes in the American Joint Committee on Cancer eighth edition cancer staging manual. CA Cancer J Clin. 2017;67(3):245-253. https://doi.org/10.3322/caac.21391
Paner GP, Stadler WM, Hansel DE, Montironi R, Lin DW, Amin MB. Updates in the Eighth Edition of the Tumor-Node-Metastasis Staging Classification for Urologic Cancers. Eur Urol. 2018;73(4):560-569. https://doi.org/10.1016/j.eururo.2017.12.018
Wenzel M, Wurnschimmel C, Chierigo F, et al. Pattern of Biopsy Gleason Grade Group 5 (4 + 5 vs 5 + 4 vs 5 + 5) Predicts Survival After Radical Prostatectomy or External Beam Radiation Therapy. Eur Urol Focus. 2022;8(3):710-717. https://doi.org/10.1016/j.euf.2021.04.011
Basourakos SP, Tzeng M, Lewicki PJ, et al. Tissue-Based Biomarkers for the Risk Stratification of Men With Clinically Localized Prostate Cancer. Front Oncol. 2021;11:676716. https://doi.org/10.3389/fonc.2021.676716
Shay JW, Wright WE. Telomeres and telomerase: three decades of progress. Nat Rev Genet. 2019;20(5):299-309. https://doi.org/10.1038/s41576-019-0099-1
Meeser A, Bartenhagen C, Werr L, et al. Reliable assessment of telomere maintenance mechanisms in neuroblastoma. Cell Biosci. 2022;12(1):160. https://doi.org/10.1186/s13578-022-00896-2
Subasri M, Shooshtari P, Watson AJ, Betts DH. Analysis of TERT Isoforms across TCGA, GTEx and CCLE Datasets. Cancers (Basel). 2021;13(8). https://doi.org/10.3390/cancers13081853
Livingstone J, Shiah YJ, Yamaguchi TN, et al. The telomere length landscape of prostate cancer. Nat Commun. 2021;12(1):6893. https://doi.org/10.1038/s41467-021-27223-6
Sharma S, Chowdhury S. Emerging mechanisms of telomerase reactivation in cancer. Trends Cancer. 2022;8(8):632-641. https://doi.org/10.1016/j.trecan.2022.03.005
Shay JW. Role of Telomeres and Telomerase in Aging and Cancer. Cancer Discov. 2016;6(6):584-93. https://doi.org/10.1158/2159-8290.CD-16-0062
Heaphy CM, Joshu CE, Barber JR, et al. The prostate tissue-based telomere biomarker as a prognostic tool for metastasis and death from prostate cancer after prostatectomy. J Pathol Clin Res. 2022;8(5):481-491. https://doi.org/10.1002/cjp2.288
Heaphy CM, Yoon GS, Peskoe SB, et al. Prostate cancer cell telomere length variability and stromal cell telomere length as prognostic markers for metastasis and death. Cancer Discov. 2013;3(10):1130-41. https://doi.org/10.1158/2159-8290.CD-13-0135
Baena-Del VJ, Zheng Q, Esopi DM, et al. MYC drives overexpression of telomerase RNA (hTR/TERC) in prostate cancer. J Pathol. 2018;244(1):11-24. https://doi.org/10.1002/path.4980
Poos AM, Schroeder C, Jaishankar N, et al. PITX1 Is a Regulator of TERT Expression in Prostate Cancer with Prognostic Power. Cancers (Basel). 2022;14(5). https://doi.org/10.3390/cancers14051267
Simmonds P, Loomis E, Curry E. DNA methylation-based chromatin compartments and ChIP-seq profiles reveal transcriptional drivers of prostate carcinogenesis. Genome Med. 2017;9(1):54. https://doi.org/10.1186/s13073-017-0443-z
Briggs FB, Goldstein BA, Mccauley JL, et al. Variation within DNA repair pathway genes and risk of multiple sclerosis. Am J Epidemiol. 2010;172(2):217-24. https://doi.org/10.1093/aje/kwq086
Joyce BT, Zheng Y, Nannini D, et al. DNA Methylation of Telomere-Related Genes and Cancer Risk. Cancer Prev Res (Phila). 2018;11(8):511-522. https://doi.org/10.1158/1940-6207.CAPR-17-0413
Chen L, Zheng Y, Jiang C, Yang C, Zhang L, Liang C. The established chemokine-related prognostic gene signature in prostate cancer: Implications for anti-androgen and immunotherapies. Front Immunol. 2022;13:1009634. https://doi.org/10.3389/fimmu.2022.1009634
Dai W, Zheng P, Luo D, et al. LPIN1 Is a Regulatory Factor Associated With Immune Response and Inflammation in Sepsis. Front Immunol. 2022;13:820164. https://doi.org/10.3389/fimmu.2022.820164
Feng D, Li D, Shi X, et al. A gene prognostic index from cellular senescence predicting metastasis and radioresistance for prostate cancer. J Transl Med. 2022;20(1):252. https://doi.org/10.1186/s12967-022-03459-8
Yu G, Wang L, Han Y, He Q. clusterProfiler: an R package for comparing biological themes among gene clusters. Omics : a journal of integrative biology. 2012;16(5)
Du J, Hou K, Mi S, et al. Malignant Evaluation and Clinical Prognostic Values of m6A RNA Methylation Regulators in Glioblastoma. Front Oncol. 2020;10:208. https://doi.org/10.3389/fonc.2020.00208
Zeng D, Ye Z, Shen R, et al. IOBR: Multi-Omics Immuno-Oncology Biological Research to Decode Tumor Microenvironment and Signatures. Front Immunol. 2021;12:687975. https://doi.org/10.3389/fimmu.2021.687975
Subramanian A, Tamayo P, Mootha VK, et al. Gene set enrichment analysis: a knowledge-based approach for interpreting genome-wide expression profiles. Proc Natl Acad Sci U S A. 2005;102(43):15545-50. https://doi.org/10.1073/pnas.0506580102
Douglas H, Lisa MC. Accessories to the Crime: Functions of Cells Recruited to the Tumor Microenvironment. Cancer Cell. 2012;21(3)
Quinn DI, Sandler HM, Horvath LG, Goldkorn A, Eastham JA. The evolution of chemotherapy for the treatment of prostate cancer. Ann Oncol. 2017;28(11):2658-2669. https://doi.org/10.1093/annonc/mdx348
Storchova Z, Kloosterman WP. The genomic characteristics and cellular origin of chromothripsis. Curr Opin Cell Biol. 2016;40:106-113. https://doi.org/10.1016/j.ceb.2016.03.003
Fordyce CA, Heaphy CM, Joste NE, Smith AY, Hunt WC, Griffith JK. Association between cancer-free survival and telomere DNA content in prostate tumors. J Urol. 2005;173(2):610-4. https://doi.org/10.1097/01.ju.0000143195.49685.ce
Zhou J, Richardson M, Reddy V, et al. Structural and functional association of androgen receptor with telomeres in prostate cancer cells. Aging (Albany NY). 2013;5(1):3-17. https://doi.org/10.18632/aging.100524
Van Poppel H, Albreht T, Basu P, Hogenhout R, Collen S, Roobol M. Serum PSA-based early detection of prostate cancer in Europe and globally: past, present and future. Nat Rev Urol. 2022;19(9):562-572. https://doi.org/10.1038/s41585-022-00638-6
Hegde PS, Karanikas V, Evers S. The Where, the When, and the How of Immune Monitoring for Cancer Immunotherapies in the Era of Checkpoint Inhibition. Clinical cancer research : an official journal of the American Association for Cancer Research. 2016;22(8)
Obermajer N, Urban J, Wieckowski E, et al. Promoting the accumulation of tumor-specific T cells in tumor tissues by dendritic cell vaccines and chemokine-modulating agents. Nat Protoc. 2018;13(2):335-357. https://doi.org/10.1038/nprot.2017.130
Saxena M, Bhardwaj N. Re-Emergence of Dendritic Cell Vaccines for Cancer Treatment. Trends Cancer. 2018;4(2):119-137. https://doi.org/10.1016/j.trecan.2017.12.007
Ness N, Andersen S, Valkov A, et al. Infiltration of CD8+ lymphocytes is an independent prognostic factor of biochemical failure-free survival in prostate cancer. Prostate. 2014;74(14):1452-61. https://doi.org/10.1002/pros.22862
Kwek SS, Lewis J, Zhang L, et al. Preexisting Levels of CD4 T Cells Expressing PD-1 Are Related to Overall Survival in Prostate Cancer Patients Treated with Ipilimumab. Cancer Immunol Res. 2015;3(9):1008-16. https://doi.org/10.1158/2326-6066.CIR-14-0227
Pittoni P, Tripodo C, Piconese S, et al. Mast cell targeting hampers prostate adenocarcinoma development but promotes the occurrence of highly malignant neuroendocrine cancers. Cancer Res. 2011;71(18):5987-97. https://doi.org/10.1158/0008-5472.CAN-11-1637
Wing K, Sakaguchi S. Regulatory T cells exert checks and balances on self tolerance and autoimmunity. Nat Immunol. 2010;11(1):7-13. https://doi.org/10.1038/ni.1818
Chanmee T, Ontong P, Konno K, Itano N. Tumor-associated macrophages as major players in the tumor microenvironment. Cancers (Basel). 2014;6(3):1670-90. https://doi.org/10.3390/cancers6031670
Jiang P, Gu S, Pan D, et al. Signatures of T cell dysfunction and exclusion predict cancer immunotherapy response. Nat Med. 2018;24(10):1550-1558. https://doi.org/10.1038/s41591-018-0136-1
Petrelli F, Ghidini M, Ghidini A, Tomasello G. Outcomes Following Immune Checkpoint Inhibitor Treatment of Patients With Microsatellite Instability-High Cancers: A Systematic Review and Meta-analysis. JAMA Oncol. 2020;6(7):1068-1071. https://doi.org/10.1001/jamaoncol.2020.1046
Roudko V, Bozkus CC, Orfanelli T, et al. Shared Immunogenic Poly-Epitope Frameshift Mutations in Microsatellite Unstable Tumors. Cell. 2020;183(6):1634-1649.e17. https://doi.org/10.1016/j.cell.2020.11.004
Schrock AB, Ouyang C, Sandhu J, et al. Tumor mutational burden is predictive of response to immune checkpoint inhibitors in MSI-high metastatic colorectal cancer. Ann Oncol. 2019;30(7):1096-1103. https://doi.org/10.1093/annonc/mdz134
Clark CA, Gupta HB, Curiel TJ. Tumor cell-intrinsic CD274/PD-L1: A novel metabolic balancing act with clinical potential. Autophagy. 2017;13(5):987-988. https://doi.org/10.1080/15548627.2017.1280223
Isaacsson VP, Antonarakis ES. PD-1/PD-L1 pathway inhibitors in advanced prostate cancer. Expert Rev Clin Pharmacol. 2018;11(5)
Santoni M, Massari F, Cheng L, et al. PD-L1 inhibitors for the treatment of prostate cancer. Curr Drug Targets. 2020
Galluzzi L, Humeau J, Buqué A, Zitvogel L, Kroemer G. Immunostimulation with chemotherapy in the era of immune checkpoint inhibitors. Nature reviews. Clinical oncology. 2020;17(12)
Yajie S, Xiangru F, Chao W, Jonathan FL, Honglin J, Jianxun D. Role of nanoparticle-mediated immunogenic cell death in cancer immunotherapy. Asian J Pharm Sci. 2020;16(prepublish)

Supplementary Figure S6 is not available with this version.

Supplementary figure S6: CeRNA network

No competing interests reported.

SupplementaryfigureS1.pdf
Supplementary figure S1. PPI graph of TRPGs, PPI: Protein-protein interaction
SupplementaryfigureS2.pdf
Supplementary figure S2. DCA analysis curves for nomogram model
SupplementaryfigureS3.pdf
Supplementary figure S3: Estimated scores of 22 immune cell subtypes
SupplementaryfigureS4.pdf
Supplementary figure S4: Differential expression of MHC-related genes
SupplementaryfigureS5.pdf
Supplementary figure S5: Enrichment of motifs and transcription factors
Supplementaryfile1.pdf
Supplementaryfile2GWASresult.pdf
Supplementarytable1.pdf

Download PDF

Editorial decision: Revision requested
11 Dec, 2023
Reviews received at journal
08 Dec, 2023
Reviewers agreed at journal
29 Nov, 2023
Reviewers agreed at journal
28 Nov, 2023
Reviewers invited by journal
22 Nov, 2023
Editor assigned by journal
21 Nov, 2023
Submission checks completed at journal
20 Nov, 2023
First submitted to journal
10 Oct, 2023

You are reading this latest preprint version

Predicting Prostate Adenocarcinoma Patients’ Survival and Immune Signature: A Novel Risk Model Based on Telomere-Related Genes

Status:

Version 1

Abstract

Background

Methods

Results

Conclusion

Figures

1 Introduction

2 Materials and Methods

2.1 Data acquisition and preprocessing

2.2 Genes function enrichment analysis

2.3 Model development and prognosis

2.4 Nomogram construction

2.5 Drug sensitivity evaluation

2.6 Immune cell infiltration analysis

2.7 Immunohistochemical (IHC) staining

2.8 Gene set variance analysis(GSVA)

2.9 Regulatory network analysis

2.10 Statistical analysis

3 Results

3.1 Screening and functional enrichment analysis of TRGs associated with PRAD prognosis

3.2 Construction and validation of the TRG-based risk model

3.3 Differential expression and regulatory networks of the model genes

3.4 Exploration of the prognostic value and clinical relevance of the risk model

3.5 Construction and validation of the risk score-related nomogram model

3.6 Exploration of signaling mechanisms related to the risk model

3.7 Analysis of the tumor immune microenvironment in the high- and low-risk groups

3.8 Correlation between risk score and immune-related genes

3.9 Sensitivity of chemotherapeutic agents in PRAD

4 Discussion

5 Conclusion

Abbreviations

Declarations

References

Supplementary Figure

Additional Declarations

Supplementary Files

Status:

Version 1