BACKGROUND: Prostate cancer (PCa) is one of most common male neoplasms. TP53 is the tumor suppressor gene with the highest correlation with human tumorigenesis discovered so far. Besides the TP53, immune-related genes attracted much attention since the clinical application of PD-1/PD-L1 (programmed death 1/programmed cell death-ligand 1) related drugs. There is currently a lack of studies that combine TP53 with immune-related genes to analyze the prognosis of prostate cancer patients. METHODS: Differentially expressed genes were filtered out by R package (edgeR) based on the TCGA-PRAD (The Cancer Genome Atlas-Prostate adenocarcinoma) data set. Using the R package (coxph), we distinguished which ones were related to survival prognosis. Constructing high and low risk groups, we used GEO (Gene Expression Omnibus) data set to verify the prediction performance. Subsequently, we explored the functional differences in gene expression between high and low risk groups. RESULTS: A total of six immune-related genes can be seen as prognostic factors in individuals with TP53 mutations. In the high-risk group, genes related to macrophage activation, epithelial cell apoptosis, and inflammation of the skin should be highly expressed. In the low-risk group, highly expressed genes are mainly involved in nucleotide phosphorylation, tRNA metabolism, and mitochondrial metabolism. CONCLUSIONS: Mutations in the TP53 gene can adversely affect the prognosis of prostate cancer and prostate cancer patients with mutations in some immune-related genes together have a worse prognosis. Compared with any other single clinical index, the prognostic score we proposed gave a more accurate forecast. In order to assist clinicians in making predictive assessments, we have also drawn a nomogram of the prognosis of prostate cancer patients. 2021 Translational Andrology and Urology. All rights reserved.
BACKGROUND: Prostate cancer (PCa) is one of most common male neoplasms. TP53 is the tumor suppressor gene with the highest correlation with human tumorigenesis discovered so far. Besides the TP53, immune-related genes attracted much attention since the clinical application of PD-1/PD-L1 (programmed death 1/programmed cell death-ligand 1) related drugs. There is currently a lack of studies that combine TP53 with immune-related genes to analyze the prognosis of prostate cancer patients. METHODS: Differentially expressed genes were filtered out by R package (edgeR) based on the TCGA-PRAD (The Cancer Genome Atlas-Prostate adenocarcinoma) data set. Using the R package (coxph), we distinguished which ones were related to survival prognosis. Constructing high and low risk groups, we used GEO (Gene Expression Omnibus) data set to verify the prediction performance. Subsequently, we explored the functional differences in gene expression between high and low risk groups. RESULTS: A total of six immune-related genes can be seen as prognostic factors in individuals with TP53 mutations. In the high-risk group, genes related to macrophage activation, epithelial cell apoptosis, and inflammation of the skin should be highly expressed. In the low-risk group, highly expressed genes are mainly involved in nucleotide phosphorylation, tRNA metabolism, and mitochondrial metabolism. CONCLUSIONS: Mutations in the TP53 gene can adversely affect the prognosis of prostate cancer and prostate cancer patients with mutations in some immune-related genes together have a worse prognosis. Compared with any other single clinical index, the prognostic score we proposed gave a more accurate forecast. In order to assist clinicians in making predictive assessments, we have also drawn a nomogram of the prognosis of prostate cancer patients. 2021 Translational Andrology and Urology. All rights reserved.
Prostate cancer (PCa) is the most frequently diagnosed male cancer in 105 countries, ranking as the second most frequent cancer and the fifth leading cause of cancer death in men (1-3). Clinically, a great number of patients are diagnosed at an advanced stage, because there is no obvious symptom in the early stage. In recent years, prostate-specific antigen testing has increased the number of men diagnosed with and treated for PCa (2).Genomic analysis is a commonly used method in cancer research worldwide, and PCa is no exception (4-7). A lot of research has been done on over-expressed genes and upregulated signal pathways in PCa, but the exact pathogenesis of the tumor is still unclear (4,5). However, some studies’ results have been translated into clinical treatments that have significantly prolonged the survival of PCa patients (3,6,8,9).Although many factors are considered to influence prognosis, due to the current promotion of individualized treatments, but treatment plans are now more precise, and more and more patients who were considered to have a poor prognosis have benefited significantly (8). Thanks to The Cancer Genome Atlas (TCGA) and GEO databases, we can further understand the genetic differences in PCa patients, such as TP53 (4,7). Compared with the differences in clinical characteristics, genomics has shown a more accurate prediction of survival prognosis in many cancers (4-6), and A great many cancers have significantly improved prognosis since immunotherapy has been applied in clinics, revealing the giant potential of immune factors in the prognosis of cancer patients (3,6,8-10). In this study, we used a database of immune genetics, and Cox proportional hazard regression models to analyze the prognostic risk factors and establish a prognostic model.We present the following article in accordance with the REMARK reporting checklist (available at http://dx.doi.org/10.21037/tau-21-179).
Methods
Data of PCa genes
We searched for PCa gene-expression data in TCGA using the UCSC Xena browser (https:xenabrowser.net) and downloaded both genetic and clinical information. The corresponding somatic mutation data was obtained through analysis of MuTect2. For research we also obtained the human HG38 v22 version of the gene annotation message in the GENCODE database. To test our results, we gathered a treatment-naive cohort of samples from patients with PCa from a GEO set (GSE116918).
Data of immune-related genes (IRGs)
IRGs were searched for in the Immunology Database and Analysis Portal (ImmPort) (https://www.immport.org/home) and a list of them was downloaded. In total, we had 1,811 IRGs.
Statistical analysis
Screening of differentially expressed (DE) IRGs
We first divided all samples into two groups based on the presence or absence of TP53 mutations. For further study, we used the R package edgeR to screen for DE genes between samples with and without TP53 mutations, requiring a fold-change (FC) >1.5 or FC <0.667, and requiring FDR (false discovery rate)-corrected P values <0.05. We excluded non-immune related genes using the list of IRGs.
Survival analysis
For convenience, we performed log2 (FPKM+1) conversion on the expression value of each gene. In order to select the genes tightly associated with survival prognosis, we used the R package coxph to analyze overall survival prognosis for all subjects, requiring P values <0.05.
Immune-related survival risk score
The immune genes related to survival prognosis were used to construct a prognosis score by weighting, and we divided the sample into high- and low-risk groups based on the median risk score, where Exp represents the expression value of the gene, and β represents the weight of the gene expression value:
Immune cell ratio analysis
With the discovery of the tumor microenvironment, the cells surrounding the tumor have become a new research target (10-13). According to research, the infiltration of immune cells into tumors is closely related to the prognosis of patients and the efficacy of immunotherapy (11,13). To confirm the ratio of immune cells in PCa, we used CIBERSORT (https://cibersort.stanford.edu/). In order to make the results more credible, the expressed differences between high- and low-risk groups were determined by t-test, and the P values <0.05 were considered statistically significant.
Tumor mutational burden (TMB)
The TMB, which is regarded as an indicator for immunotherapy, is also used in PCa (14). According to the definition of TMB, we removed synonymous mutations and intron mutations from the PCa mutation data. We calculated the TMB score for both every single sample and the total samples (Table S1):The study was conducted in accordance with the Declaration of Helsinki (as revised in 2013).
Results
Identification of DE immune genes
The mutation data downloaded from TCGA was first used to calculate the overall TMB score, which was 1.16. After analysis, we found a number of mutant genes, among which TP53 had the highest mutation frequency ().
Figure 1
Gene mutation map of prostate cancer.
Gene mutation map of prostate cancer.Using edgeR, a total of 3,101 DE genes (1,089 upregulated; 2,012 downregulated) were screened from samples with and without TP53 mutations ().
Figure 2
Heat map of differentially expressed gene between patients with and without TP53 mutations.
Figure 3
Volcano plot of differentially expressed gene between patients with and without TP53 gene mutations.
Heat map of differentially expressed gene between patients with and without TP53 mutations.Volcano plot of differentially expressed gene between patients with and without TP53 gene mutations.There were 1,811 IRGs downloaded from the ImmPort database, and we found 178 genes (84 upregulated; 94 downregulated) that simultaneously matched the differential expression in the TP53 mutant samples of PCa and immune correlation ().
Figure 4
Venn diagram of differentially expressed (DE) genes and immune-related genes (IRGs).
Figure 5
Heat map of differentially expressed immune-related genes between patients with and without TP53 mutations.
Venn diagram of differentially expressed (DE) genes and immune-related genes (IRGs).Heat map of differentially expressed immune-related genes between patients with and without TP53 mutations.
Enrichment analysis of DE IRGs
The GO(Gene Ontology) term ‘enrichment analysis’ showed that DE IRGs were mainly concentrated in immune response, complement activation, immunoglobulin etc. Moreover, they were also enriched in some pathways, such as neuroactive ligand-receptor interactions, cytokine ligand-receptor interactions, and chemokine signaling pathways (,
Figure S1).
Figure 6
Enrichment analysis of differentially expressed immune-related genes.
Figure 7
Chordal graph of differentially expressed immune-related genes.
Enrichment analysis of differentially expressed immune-related genes.Chordal graph of differentially expressed immune-related genes.Through previous analysis, we already had 178 target genes. Based on the patient’s clinical data, we performed a survival analysis. In the results, there were significant differences in the survival prognosis of six genes (P value <0.05).Three genes (DES, HBEGF and OPRK1) positively correlated with survival prognosis On the contrary, the other three genes (CALCB, OBP2A and UCN3) negatively affected prognosis (,
).
Table 1
Survival analysis of differentially expressed immune-related genes
Gene
Coef
Hazard Ratio
P value
OBP2A
1.44
4.221
0
DES
−0.28
0.756
0.042
CALCB
5.521
249.818
0.001
HBEGF
−0.847
0.429
0.039
UCN3
0.521
1.684
0.017
OPRK1
−1.328
0.265
0.037
Figure 8
Differentially expressed immune-related genes related to survival. (A) The effect of differential expression of CALCB on survival; (B) the effect of differential expression of DES on survival; (C) the effect of differential expression of HBEGF on survival; (D) the effect of differential expression of OBP2A on survival; (E) the effect of differential expression of OPRK1 on survival; (F) the effect of differential expression of UCN3 on survival.
Differentially expressed immune-related genes related to survival. (A) The effect of differential expression of CALCB on survival; (B) the effect of differential expression of DES on survival; (C) the effect of differential expression of HBEGF on survival; (D) the effect of differential expression of OBP2A on survival; (E) the effect of differential expression of OPRK1 on survival; (F) the effect of differential expression of UCN3 on survival.
Risk score for patients with PCa
The six identified genes were used to calculate each patient’s risk score through weighting. Using these scores, we regrouped all patients as high risk or low risk. There was a significant difference in survival prognosis between groups (P value =0.0072). The high-risk group had a worse survival prognosis, and the low-risk group had a better survival prognosis ().
Figure 9
Survival curves of high and low risk score groups.
Figure 10
Risk factor chart for high and low risk score groups. (A,B,C) The definition of the X axis is the serial number of each individual. (A) Scatter plot of individual distribution of high- and low-risk groups; (B) Scatter plot of full sample risk scores from low to high; (C) gene expression heat map of each sample from low risk to high risk.
Survival curves of high and low risk score groups.Risk factor chart for high and low risk score groups. (A,B,C) The definition of the X axis is the serial number of each individual. (A) Scatter plot of individual distribution of high- and low-risk groups; (B) Scatter plot of full sample risk scores from low to high; (C) gene expression heat map of each sample from low risk to high risk.
External dataset validation
The GSE116918 dataset was used to verify the accuracy of our risk score. Prior to validation, we used the optimal threshold (risk score =17.21523) to distinguish between the high and low risk groups. As a result, we found that the high-risk group had a significantly lower survival prognosis than the low-risk group in this dataset (P value =0.037) (,
Figure S2).
Figure 11
Survival curves from GSE116918 dataset.
Survival curves from GSE116918 dataset.
Prognostic performance of PCa risk score
In TCGA’s PCa data, the risk score made good predictions of patient’s 1-, 3-, and 5-year survival. The same conclusion was also obtained with the GSE116818 dataset (,
Figure S3).
Figure 12
ROC curve of 1-, 3-, and 5-year survival rate in TCGA data set.
ROC curve of 1-, 3-, and 5-year survival rate in TCGA data set.
Survival analysis of clinical factors of PCa
Afterwards, we classified all kinds of clinical factors, and found through survival analysis that each clinical factor has no significant difference in survival prognosis, but our risk score can have a good prediction of survival prognosis, which again proves that we have constructed The prognosis score was better than clinical factors. It is of great help to the clinic ().
Figure 13
Survival analysis of clinical factors of prostate cancer data set. (A) Survival difference between TP53 mutation group and control group; (B) survival difference between the 60-year-old group and the 60-year-old group; (C) survival difference between T2 group and T3-4 group; (D) difference in survival between N0 group and N1 group.
Figure 14
Violin diagram of clinical factors of prostate cancer data set. (A) Distribution status and distribution density of TP53 mutation group and control group; (B) the distribution status and density of the group over 60 years old and the group under 60 years old; (C) distribution status and distribution density of T2 group and T3-4 group; (D) distribution status and distribution density of N0 group and N1 group.
Survival analysis of clinical factors of prostate cancer data set. (A) Survival difference between TP53 mutation group and control group; (B) survival difference between the 60-year-old group and the 60-year-old group; (C) survival difference between T2 group and T3-4 group; (D) difference in survival between N0 group and N1 group.Violin diagram of clinical factors of prostate cancer data set. (A) Distribution status and distribution density of TP53 mutation group and control group; (B) the distribution status and density of the group over 60 years old and the group under 60 years old; (C) distribution status and distribution density of T2 group and T3-4 group; (D) distribution status and distribution density of N0 group and N1 group.
GSEA enrichment analysis of high and low risk groups
We performed Gene Set Enrichment Analysis (GSEA) to DE genes in high and low risk groups and found that both groups had some over-expressed genes, which were involved in many biochemical reaction. Genes that are highly expressed in the high-risk group are mainly involved in macrophage activation, epithelial cell apoptosis, and inflammatory response. In the low-risk group, the over-expressed genes mainly participate in nucleotide phosphorylation, tRNA metabolism, and mitochondrial gene expression ().
Figure 15
GSEA enrichment analysis of prostate cancer high and low risk groups. (A) Enrichment plot: macrophage activation; (B) enrichment plot: Positive regulation of tyrosine phosphorylation of STAT protein; (C) enrichment plot: positive regulation of peptidyl tyrosine phosphorylation; (D) enrichment plot: regulation of epithelial cell apoptotic process; (E) enrichment plot: Negative regulation of inflammatory response; (F) enrichment plot: regulation of reactive oxygen species biosynthetic process; (G) enrichment plot: nucleotide phosphorylation; (H) enrichment plot: nucleoside diphosphate kinase activity; (I) enrichment plot: catalytic activity acting on a tRNA; (J) enrichment plot: tRNA metabolic process; (K) enrichment plot: mitochondrial gene expression; (L) enrichment plot: respiratory chain complex IV assembly.
GSEA enrichment analysis of prostate cancer high and low risk groups. (A) Enrichment plot: macrophage activation; (B) enrichment plot: Positive regulation of tyrosine phosphorylation of STAT protein; (C) enrichment plot: positive regulation of peptidyl tyrosine phosphorylation; (D) enrichment plot: regulation of epithelial cell apoptotic process; (E) enrichment plot: Negative regulation of inflammatory response; (F) enrichment plot: regulation of reactive oxygen species biosynthetic process; (G) enrichment plot: nucleotide phosphorylation; (H) enrichment plot: nucleoside diphosphate kinase activity; (I) enrichment plot: catalytic activity acting on a tRNA; (J) enrichment plot: tRNA metabolic process; (K) enrichment plot: mitochondrial gene expression; (L) enrichment plot: respiratory chain complex IV assembly.
Comparison of the proportion of immune infiltrating cells
We analyzed in detail the difference in the proportion of immune cell infiltration between the two groups, but only found that the mast cell active had the significant differences. TIDE, a computational method to model two primary mechanisms of tumor immune evasion, predicted the outcome of melanoma patients treated with first-line anti-PD1 or anti-CTLA4 more accurately than other biomarkers such as PD-L1 level and mutation load (15). So we used the TIDE score (http://tide.dfci.harvard.edu/) to predict the efficacy of immunotherapy in high- and low-risk groups of PCa. However, there was no significant difference between the two groups (P value =0.7) (,
Figure S4).
Figure 16
Difference in the proportion of immune infiltrating cells in high and low risk prostate cancer groups.
Difference in the proportion of immune infiltrating cells in high and low risk prostate cancer groups.
Nomogram of PCa
Through univariate Cox regression analysis for TP53 gene mutation, age, T stage, N stage, and high- and low-risk group in patients with PCa. Excitingly, requiring P value <0.05, we found that the risk group had more significant ability to predict survival prognosis (Table S2).Finally, a nomogram of PCa was constructed, and the survival prognosis of patients was assisted by scoring their various factors ().
Figure 17
nomogram model of prostate cancer patients.
nomogram model of prostate cancer patients.
Discussion
PCa is a high-profile cancer in man with an increasing morbidity and mortality. In this topic, we constructed a prognostic model through the DE immune related genes significantly associated with survival. If we regarded the median score as the demarcation between high and low risk groups, there was no significant difference in two groups from the GSE116918. So we reselected the optimal threshold (risk score =17.21523) to define high and low risk groups, and the result showed that the survival difference between high and low risk groups is significant. Judgment of prognosis is powerlessly only by clinical indicators, nevertheless our prognostic score can offer an excellent assessment for patients.The results of this study have yet to be further confirmed by prospective studies or the inclusion of more experimental subjects. The occurrence and progression of cancer is the combined effect of a variety of gene abnormalities, including the loss or silencing of tumor suppressor genes, and the amplification or activation of proto-oncogenes. The regulation of gene expression in cancer cells is also significantly different from normal cells, such as miRNA and lncRNA. Prostate cancer is heterogeneous, which limits the efficacy of drugs. Therefore, the treatment of prostate should be precise and individualized. Adjuvant treatment plans should be formulated comprehensively based on the results of gene sequencing, combined with molecular subtypes.
Conclusions
The mutation of TP53 and the abnormal expression of some related immune genes have a significant impact on the prognosis of prostate cancer patients. DES, HBEGF, OPRK1, CACLB, OBP2A and UCN3 are prognostic-related genes, and these genes are closely related to the immune response. The comprehensive score of these six genes serves as the criterion for defining the high and low risk groups. The prognosis of the high-risk group is much worse than that of the low-risk group.The article’s supplementary files as
Authors: Anis A Hamid; Kathryn P Gray; Grace Shaw; Laura E MacConaill; Carolyn Evan; Brandon Bernard; Massimo Loda; Niall M Corcoran; Eliezer M Van Allen; Atish D Choudhury; Christopher J Sweeney Journal: Eur Urol Date: 2018-12-12 Impact factor: 20.096
Authors: Freddie Bray; Jacques Ferlay; Isabelle Soerjomataram; Rebecca L Siegel; Lindsey A Torre; Ahmedin Jemal Journal: CA Cancer J Clin Date: 2018-09-12 Impact factor: 508.702
Authors: Lien Spans; Thomas Van den Broeck; Elien Smeets; Stefan Prekovic; Bernard Thienpont; Diether Lambrechts; R Jeffrey Karnes; Nicholas Erho; Mohammed Alshalalfa; Elai Davicioni; Christine Helsen; Thomas Gevaert; Lorenzo Tosco; Karin Haustermans; Evelyne Lerut; Steven Joniau; Frank Claessens Journal: Oncotarget Date: 2016-04-26
Authors: Kasey Jividen; Katarzyna Z Kedzierska; Chun-Song Yang; Karol Szlachta; Aakrosh Ratan; Bryce M Paschal Journal: BMC Cancer Date: 2018-10-10 Impact factor: 4.430
Authors: Shun Wan; Yang He; Bin Zhang; Zhi Yang; Fang-Ming Du; Chun-Peng Zhang; Yu-Qiang Fu; Jun Mi Journal: Front Oncol Date: 2022-04-05 Impact factor: 5.738