A circadian rhythm-related gene signature for predicting relapse risk and immunotherapeutic effect in prostate adenocarcinoma

Abstract

Prostate adenocarcinoma (PRAD) represents the most common male carcinoma in developed countries, its high relapse risk contributes to the second-leading cause of cancer-related deaths. Therefore, it is required to develop an effective signature for predicting the relapse risk of PRAD. To identify a circadian rhythm- (CR-) related predictive signature, we analyzed RNA-seq data of patients with prostate adenocarcinoma (PRAD) from the TCGA and GEO cohort. Seven circadian rhythm- (CR-) related genes (FBXL22, MTA1, TP53, RORC, DRD4, PPARGC1A, ZFHX3) were eventually identified to develop a CR-related signature. AUCs for 3-year overall survival were 0.852, 0.856 and 0.944 in the training set, validation set and an external independent test set (GSE70768), respectively. Kaplan-Meier curve analysis showed that the high-risk group has a reduced relapse-free survival (RFS) than the low-risk group in the training set, validation set, and test set, respectively (P < 0.05). We constructed a prognostic nomogram combining the CR-related signature with T staging to precisely estimate relapse risk of PRAD patients. Finally, we observed that the CR-related gene signature was associated with tumor mutation burden, multiple immune checkpoint molecules and microsatellite instability, and thus could predict response to immune checkpoint inhibitors in PRAD. Conclusively, we developed a circadian rhythm-related gene signature for predicting RFS and immunotherapy efficacy in prostate adenocarcinoma.

Keywords: FBXL22; circadian rhythm; immune checkpoint inhibitor; prognosis; prostate adenocarcinoma.

Figure 1

Establishment of a circadian rhythm- (CR-) related gene signature in prostate adenocarcinoma (PRAD). (A) Enrichment score (ES) of CR-related gene set was significantly enriched in normal tissue than tumor tissue. (B) High-ES patients had an improved relapse-free survival (RFS) than their counterparts. (C) CR-related gene set was positively enriched in normal tissue compared with tumor tissue. (D) 13 CR-related genes were qualified in univariate Cox regression analysis (FBXL22, FBXL6, MTA1, NR2F6, TP53, BTRC, GHRL, RORC, CIPC, DRD4, MTOR, PPARGC1A, ZFHX3; P<0.05). (E, F) 13 qualified genes were further filtered using LASSO regression analysis to eliminate multicollinearity and seven eligible genes (FBXL22, MTA1, TP53, RORC, DRD4, PPARGC1A, ZFHX3) were eventually acquired for establishment of a CR-related gene signature for predicting RFS in PRAD patients.

Figure 2

Evaluation of the performance of the gene signature. (A) AUC for predicting three-year RFS is 0.852 in the training set (n=236). (B) AUC for predicting three-year RFS is 0.856 in the validation set (n=102). (C) Low-risk patients had an improved RFS than high-risk patients in the training set. (D) Low-risk patients had an improved RFS than high-risk patients in the validation set.

Figure 3

The discriminative power of the CR-related signature. (A) Principle component analysis (PCA) demonstrated that low-risk group was apparently distinct from high-risk group in Dim 1 in the training set. (B) As risk score increased, chance of tumor recurrence increased in the training set. (C) Principle component analysis (PCA) showed that the low-risk group was apparently distinct from the high-risk group in Dim 1 in the validation set. (D) As risk score increased, chance of tumor recurrence increased in the validation set.

Figure 4

Comparison of the CR-related gene signature with other indicators for RFS in PRAD. (A) CR-related gene signature showed an improved predictive performance than other clinical indicators, with an AUC of 0.831. (B) CR-related gene signature (Signature 1) showed an improved predictive performance than other four reported gene signatures. (C) CR-related gene signature also showed an ideal predictive performance in an independent cohort of PRAD (GSE70768), with an AUC of 0.944. (D) Low-risk group had an improved RFS than high-risk group in an independent cohort of PRAD (GSE70768) (with a cutoff point of median value; log-rank test, P = 0.003).

Figure 5

Functional enrichment analysis for CR-related gene signature. (A) Enriched biological processes included cytoplasmic translation, oxidative phosphorylation, mitochondrial respiratory chain complex assembly. (B) Enriched cell components included mitochondrial protein-containing complex, ribosome, inner mitochondrial membrane protein complex. (C) Enriched molecular functions included structural constituent of ribosome, oxidoreduction-driven active transmembrane transporter activity and electron transfer activity. (D) Enriched KEGG pathways included ribosome, diabetic cardiomyopathy, chemical carcinogenesis (E, F) Gene set enrichment analysis (GSEA) showed the top 10 KEGG signaling pathway including AMPK signaling pathway, dopaminergic synapse, central carbon metabolism in cancer, NF-kappa B signaling pathway and others.

Figure 6

Clinical application of CR-related gene signature for PRAD patients. (A) Multivariate Cox regression was used to investigate the ability of age, clinical T stage, clinical N stage, Gleason score, prostate specific antigen (PSA) and the CR-related gene signature (risk score) to predict RFS. Risk score remained to be a valid predictor (P<0.05). (B) A nomogram combing risk score and clinical T stage was constructed to predict 1-, 3- and 5-year RFS for individual PRAD patient. (C) ROC analysis for nomogram showed an impressive predictive performance, with AUC of 0.944. (D) Calibration curve showed agreement between actual and predicted RFS, indicating an ideal predictive capability.

Figure 7

Association of the CR-related gene signature with tumor infiltrating immune cells in PRAD. (A) Comparison of tumor infiltrating immune cells between low- and high-risk groups demonstrated that there existed a significant difference in the abundance of B cell naive, B cell plasm, T cell CD4 memory, T cell follicular helper, and monocyte (P < 0.05). (B–D) Risk score was inversely correlated with B cell naive, B cell plasm, T cell CD4 memory, T cell follicular helper, and monocyte (R < 0, P < 0.05). (E) Risk score was positively correlated with T cell CD4 memory (R < 0, P < 0.05). (F) Risk score was inversely correlated with monocyte (R < 0, P < 0.05).

Figure 8

The profiling of somatic nucleotide variation for PRAD patients of different risk groups. (A, B) The waterfall plot of nucleotide variation rate in high-risk group and low-risk group. (C, D) The bar plot and box plot analyses of risk scores in the high-mutation group and low-mutation group. (E) The correlation of tumor mutation burden and risk score. *P <0.05, ** P <0.01, ***P<0.001.

Figure 9

The effects of CR-related risk score on response to ICIs. (A) MSI levels was elevated in the high-risk patients than in low-risk patients (P <0.001). (B) MSI levels were positively correlated with risk score (P < 0.05, R = 0.37). (C) Risk score was significantly increased in MSI-H group than MSI-L/MSS group (data of MSI generated from PreMSIm R package; P < 0.001). (D) Risk score was significantly higher in MSI-H group than MSI-L/MSS group (data of MSI generated from UCSCXenaShiny R package; P < 0.05). (E) Comparison of multiple immune checkpoint mRNA (PDCD1, LAG3, CD40, CTLA4, PDCD1LG2) in the high-risk patients and the low-risk patients. *P <0.05, ** P <0.01, ***P<0.001.