A Novel Angiogenesis-Related Gene Signature to Predict Biochemical Recurrence of Patients with Prostate Cancer following Radical Therapy

Abstract

Background: Prostate cancer (PCa) is one of the most common malignancies in males; we aim to establish a novel angiogenesis-related gene signature for biochemical recurrence (BCR) prediction in PCa patients following radical therapy.

Methods: Gene expression profiles and corresponding clinicopathological data were acquired from The Cancer Genome Atlas (TCGA) and Gene Expression Omnibus (GEO) database. We quantified the levels of various cancer hallmarks and identified angiogenesis as the primary risk factor for BCR. Then machine learning methods combined with Cox regression analysis were used to screen prognostic genes and construct an angiogenesis-related gene signature, which was further verified in external cohorts. Furthermore, estimation of immune cell abundance and prediction of drug responses were also conducted to detect potential mechanisms.

Results: Angiogenesis was regarded as the leading risk factor for BCR in PCa patients (HR = 1.58, 95% CI: 1.38-1.81), and a novel prognostic signature based on three genes (NRP1, JAG2, and VCAN) was developed in the training cohort and successfully validated in another three independent cohorts. In all datasets, this signature was identified as a prognostic factor with promising and robust predictive values. Moreover, it also predicted higher infiltration of regulatory T cells and M2-polarized macrophages and increased drug sensitivity of angiogenesis inhibitors in high-risk patients.

Conclusions: The angiogenesis-related three-gene signature may serve as an independent prognostic factor for BCR, which would contribute to risk stratification and personalized management of PCa patients after radical therapy in clinical practice.

Figure 1

Flowchart of this study. Four stages are included in our study. Firstly, angiogenesis was identified as the primary risk factor for BCR by single sample gene set enrichment analysis and survival analysis. Secondly, random survival forest and cox regression analyses screened prognostic genes and constructed an angiogenesis-related gene signature. Thirdly, the predictive ability of the novel signature was further verified in another three validation cohorts. Moreover, functional enrichment and immune cell infiltration analyses and drug sensitivity prediction were also conducted to perform further investigation.

Figure 2

Identification of angiogenesis as the leading risk factor in BCR of PCa patients. (a) Heatmap shows the ssGSEA score of prognostic hallmark activities in four datasets. Angiogenesis ranks first among prognostic factors leading to BCR (b) and its influence on BCR was further validated through meta-analysis (c). (d) Kaplan–Meier (K-M) survival curves of BCR survival between high- and low-score groups stratified by the angiogenesis score. (e) Percentile chart of BCR distribution between high- and low-score groups. (f) Violin plot illustrated different stromal, immune, and estimate scores between groups. (g) Box plot of normalized angiogenesis-related score in BCR or non-BCR patients. BCR: biochemical recurrence. ^∗∗P < 0.01; ^∗∗∗P < 0.001.

Figure 3

Construction of the angiogenesis-related gene signature. (a) Box plot shows the aberrant expression levels of angiogenesis-related genes in BCR in comparison with non-BCR patients. (b) The interaction of 12 angiogenesis-related prognostic genes in PCa. Expression differences between BCR and non-BCR tissues were depicted in different colors. Genes upregulated in BCR: red; genes with no significant alterations: green. Risk factors are colored in blue while protective factors are depicted in yellow. The lines connecting genes represented their interaction with each other. The size of each circle represents the prognosis effect of each gene scaled by P value. (c) Variable importance plot of the random survival forest analysis comparing rankings between minimal depth and variable of importance (VIMP). The VIMP rank is reported on the x-axis. The minimal depth (rank order) is on the y-axis. The vertical line divides variables with positive VIMP (left) from those with negative VIMP (right; unimportant). The horizontal line indicates the minimal depth threshold: important variables are below the line. The variables on the diagonal red line are those ranked equally by the two methods. (d) Venn diagram shows commonly selected six genes by univariate cox regression and random survival forest analyses. (e) Bar graph displays coefficients of each gene included in the signature.

Figure 4

Angiogenesis-related gene signature serves as an independent prognostic factor with promising and robust predictive values. Forest plot illustrates the risk score calculated by angiogenesis-related gene signature is an independent risk factor for BCR in the training (a), validation I (d), validation II (g), and validation III (j) cohorts, respectively. K-M survival curves show poor prognosis in high-risk patients divided by the gene signature in four cohorts (b), (e), (h), and (k). ROC curves illustrate the promising and stable predictive ability of the gene signature in four cohorts (c), (f), (i), and (l). This novel signature can also be a prognostic factor for metastasis (m), and it is suitable for risk stratification and metastasis prediction in the validation I cohort (n), (o).

Figure 5

Functional enrichment analysis, immune cell estimation, and drug sensitivity prediction in high- and low-risk groups. Lollipop plot shows pathways significantly enriched in high-risk groups in both training and whole validation cohorts (a). GSEA plots of enriched hallmark activities in the training (b) and whole validation cohorts (c). Box plots show the abundance of immune cells between high- and low-risk groups in the training (d) and whole validation cohorts (e). Bee graph of the predictive IC₅₀ values for three anti-angiogenesis agents in the training (f) and whole validation cohorts (g). IC₅₀: half-maximal inhibitory concentration. ^∗P < 0.05; ^∗∗P < 0.01; ^∗∗∗P < 0.001.