- 1Department of Urology, Zhujiang Hospital of Southern Medical University, Guangzhou, China
- 2Department of Urology, The First Affiliated Hospital of Jinan University, Guangzhou, China
- 3Guangdong Key Laboratory of Urology, Department of Urology, Minimally Invasive Surgery Center, The First Affiliated Hospital of Guangzhou Medical University, Guangzhou Urology Research Institute, Guangzhou, China
Prostate cancer (PCa) is a high morbidity malignancy in males, and biochemical recurrence (BCR) may appear after the surgery. Our study is designed to build up a risk score model using circular RNA sequencing data for PCa. The dataset is from the GEO database, using a cohort of 144 patients in Canada. We removed the low abundance circRNAs (FPKM < 1) and obtained 546 circRNAs for the next step. BCR-related circRNAs were selected by Logistic regression using the “survival” and “survminer” R package. Least absolute shrinkage and selector operation (LASSO) regression with 10-fold cross-validation and penalty was used to construct a risk score model by “glmnet” R software package. In total, eight circRNAs (including circ_30029, circ_117300, circ_176436, circ_112897, circ_112897, circ_178252, circ_115617, circ_14736, and circ_17720) were involved in our risk score model. Further, we employed differentially expressed mRNAs between high and low risk score groups. The following Gene Ontology (GO) analysis were visualized by Omicshare Online tools. As per the GO analysis results, tumor immune microenvironment related pathways are significantly enriched. “CIBERSORT” and “ESTIMATE” R package were used to detect tumor-infiltrating immune cells and compare the level of microenvironment scores between high and low risk score groups. What’s more, we verified two of eight circRNA’s (circ_14736 and circ_17720) circular characteristics and tested their biological function with qPCR and CCK8 in vitro. circ_14736 and circ_17720 were detected in exosomes of PCa patients’ plasma. This is the first bioinformatics study to establish a prognosis model for prostate cancer using circRNA. These circRNAs were associated with CD8+ T cell activities and may serve as a circRNA-based liquid biopsy panel for disease prognosis.
Introduction
Prostate cancer (PCa) is a high morbidity solid tumor most commonly detected in older males in developed countries (Siegel et al., 2020). Radical prostatectomy (RP) offers great therapeutic effect for most cases of localized PCa. Nonetheless, up to 50% patients will experience biochemical recurrence (BCR) after the operation (Freedland et al., 2005; Moris et al., 2020). Research has shown that the prognosis of prostate cancer depends on many factors, such as age, positive margin, pathological stage, Gleason score, PSA, and so on (Kreuz et al., 2020). Due to the heterogeneity of prostate cancer, the existing prognosis judgment system currently fails to meet the needs of every patient. There is an urgent need for a new prostate cancer prognosis prediction model to be constructed.
In the human genome, only 2% of genes are protein coding and up to 98% are non-coding RNAs. The discovery of non-coding RNAs has completely changed our understanding of cancer research (Slack and Chinnaiyan, 2019). The prognosis value of microRNA (miRNA) and long non-coding RNA (lncRNA) are widely reported in cancers (Yang et al., 2019; Mitra et al., 2020; Zhou et al., 2020). As the use of sequencing methods becomes more popular, new kinds of non-coding RNA are emerging. Circular RNA (circRNA) is a new category of non-coding RNA which has a covalent circular structure (Kristensen et al., 2019). By comparing them with other linear RNAs, circRNAs have been shown to exist widely in body fluids and are more difficult to degrade because of their framework (Su et al., 2019). What’s more, circRNAs are ideal biomarkers based on their tissue specific, stable, and easy to detect features (Salzman et al., 2013). The latest reports show that circRNAs play important roles in PCa progression, including proliferation, epithelial-mesenchymal transition (EMT), and apoptosis (Feng et al., 2019; Shen et al., 2020; Shi et al., 2020). Progressions in high-throughput sequencing methods and updated bioinformatics algorithms have facilitated research into circRNAs (Santer et al., 2019). However, there are no reports of prognostic models on circRNAs signature in PCa so far.
In our study, we screened high abundance BCR-related circRNAs from a cohort of 144 PCa patients and developed a new circRNA risk score model to predict the prognosis. The GO analysis of differentially expressed mRNAs between high and low risk score groups shows that the circRNAs may affect BCR using the tumor microenvironment (TME). What’s more, the following CIBERSORT and ESTIMATE analyses proves their effect in TME. In addition, we verified two circRNA’s (circ_14736 and circ_17720) characteristics and tested their biological function with qPCR and CCK8 in vitro. In summary, we are the first in constructing this circRNA risk score model to predict the prognosis of PCa, which may be helpful in guiding therapeutic strategy.
Materials and Methods
Ethics Statement
This study was approved by the Ethics Committee of Zhujiang Hospital, Southern medical University. Informed consent forms were signed by all patients. According to the ethical and legal standards, three blood samples from PCa patients were handled and made anonymous.
Dataset Preparation
circRNA sequencing data and corresponding clinical follow-up information were downloaded from the Gene Expression Omnibus (GEO) database1. We obtained the GSE113124 (Chen et al., 2019) expression profiles, including 144 patients with prostate cancer, for the next analysis by removing the circRNAs of FPKM (Fragments Per Kilobase Million) < 1 from expression dataset. The clinical characteristics of the patients are summarized in Table 1.
Construction of LASSO Model
circRNAs with FPKM > 1 were analyzed the Logistic regression; P < 0.05 was used as the significance cut-off to identify candidate circRNAs associated with BCR and the median of sample expression was used as the cutoff of high and low expression groups. Twenty eight circRNAs were selected, and the least absolute shrinkage and selector operation (LASSO) regression (Ramsay et al., 2018) with 10-fold cross-validation and penalty was used to construct a prognostic prediction model by “glmnet” R software package. Patients of the GSE113124 dataset were randomly assigned to 10 groups; one of these was the test set and the remaining nine groups were the training set. The risk score was calculated by expression profile data and the coefficient of the corresponding circRNAs. The formula is as follows:
where n is the number of circRNAs in the prognostic prediction model, Coef(i) represents the coefficient, and X(i) means the relative circRNAs expression level identified by LASSO regression.
Assessment of the Prognostic Prediction Model
Kaplan-Meier method was employed to draw the survival curves and assess the BCR of high and low risk score groups. Receiver operating characteristic (ROC) curve from “pROC” R package was performed to assess the sensitivity and specificity of the model compared with PSA, Gleason score, and pathologic stage. We identified differentially expressed mRNAs between different risk-score groups with FDR < 0.05 and —log2FoldChange— > 1, and Gene Ontology (GO) pathway analysis of those mRNAs were performed using OmicShare online tool. The background gene list is all Homo sapiens genes with Ensembl release 96 and Ensembl Genomes 43.
Nomogram Construction
The Cox proportional hazard model for 5-year BCR free probability were fit to 144 patients by an “rms” R package, and the variables included PSA, clinical T stage, Gleason score, risk-score, and C-index of the nomogram.
Immune Analysis by CIBERSORT and ESTIMATE
To calculate the tumor-infiltrating immune cells (TIICs) in the tumor samples, we used normalized gene expression profiles and LM22 signature matrix at 1,000 permutations to run the CIBERSORT1, which is a deconvolution algorithm website. Then the result of CIBERSORT was employed to draw a boxplot of immune cells based on risk-score group. For comparing the level of microenvironment scores between high and low risk groups, we obtained all RNA expression profiles by annotating with GRCh37. The result of running the ESTIMATE package (Yoshihara et al., 2013) was conducted to study the difference of microenvironments between the two groups.
Cell Lines and Cell Culture
The prostate cancer cell lines PC3, DU145, C4-2, and LNCaP, were all purchased from Beina Biotechnology Research Institute. PC3 and DU145 were cultured in DMEM medium. C4-2 and LNCaP were cultured in RPMI-1640 medium. Both DMEM and RPMI-1640 medium were supplemented with 10% fetal bovine and 1% Penicillin and Streptomycin. The immortalized prostate epithelial cell line RWPE-1 was recovered from the nitrogen liquid tank in our lab and was maintained in KM medium with 1% recombinant epithelial growth factor. All the cell lines above were cultivated in the incubator with constant 5% CO2 at 37°C.
RNase R Treatment, RT-qPCR, and Agarose Gel Electrophoresis
1 × 106 cells were prepared for total RNA extraction with Trizol reagent. RNA from cytoplasmic and nuclear sections were separated according to the manufacturer of PARISTM Kit (Life technologies). To acquire pure circRNA, total RNA was digested by RNase R (Rev: 20180502; Geneseed) at 37°C for 10 min. Then, RNA was transcribed into cDNA with HiScript II Q RT SuperMix for qPCR (R223-01, Vazyme). SYBR Green real-time PCR Master Mix (QPK-201, TOYOBO) was used to assess the expression level of target circRNAs by Bio-Rad CFX96 PCR machine. β-actin was used as the inner reference. Finally, the PCR product was separated by 1.5% agarose gel electrophoresis to detect the existence and molecule weight of each targeted gene. Convergent and divergent primers of circ_14734, circ_17720, and GAPDH are shown in Supplementary Table 1.
SiRNA Transfection and Cell Proliferation Assay
All small interference RNAs (siRNA) were synthesized by Hippo Biotechnology Company. DU145 and PC3 cells were seeded into the 6-well plates, then the negative control siRNA and siRNAs targeted to objective genes mixed with lipofectamine 3000 (L3000008, ThermoFisher) were added to the prepared cells in the fusion of 60–70%. RT-qPCR was used to detect the knockdown efficiency in DU145 and PC3 cells. Cell proliferation assay was carried out with CCK-8 kit (C0042, Beyotime Biotechnology),and cell viability was evaluated by the absorbance at 450nm using the microplate reader (Bio-Rad iMark). More details can be found in our previous study (Mo et al., 2017). SiRNAs sequence targeted to circ_14734 and circ_17720 are shown in Supplementary Table 1.
Exosome Extraction and Identification
Before exosome isolation, the cultured cells were maintained in exosome free medium for 24 h. In brief, 10 ml culture medium with exosome free FBS was collected and centrifuged at 1,000 rpm for 10 min. Then the supernatant was filtrated with 0.22 μm filter to remove the cell debris. Subsequently, the medium containing exosomes was incubated with one quarter volume of PEG (polyethylene glycol) 8000 buffer overnight at 4°C (Sun et al., 2019). Finally, the mixture was centrifuged at 3,000 rpm for 30 min and re-suspended with 20 μl PBS. The morphology of exosomes was detected by Transmission Electron Microscope (TME), and the particle size distribution was examined via NanoSight. The exosome markers, including HSP70 (Proteintech, No.: 10995-1-AP) and CD63 (Proteintech, No.: 25682-1-AP), were examined by Western blot, and the details of this procedure can be found in our previous study (Mo et al., 2017).
Statistical Analysis
R (version 3.6.1)3 and related packages were applied to all statistical analyses. The Wilcoxon test was used to compare two independent non-parametric samples. The Kaplan-Meier survival curves were verified with the log-rank test. We used mean ± standard deviation to describe the continuous variables in normal distribution while the median (range) was applied to testify the continuous variables in abnormal distribution. P < 0.05 was regarded as statistically significant.
Results
Constructing a CircRNA’s Signature Model to Predict BCR in 144 Patients With Prostate Cancer
Our workflow is showing in Figure 1. By analyzing expression profiles of circRNA in 144 PCa patients, the filter with mean FPKM for circRNAs less than 1 were identified 546 circRNAs from GSE113124. To study the relationship between circRNAs and biochemical recurrence, those circRNAs were subjected to logistic regression with low or high expression levels of each circRNA as the variable. Using a P-value of 0.05 as the cut-off criterion for circRNAs selection, we obtained 28 circRNAs that were significantly linked to biochemical recurrence. The expressions of those 28 circRNAs are shown as a heatmap in Figure 2. Target circRNAs were selected and based on a LASSO regression model. The best model was determined by L1 regularization and 10-fold cross-validation with the AUC (Area Under Curve) as the benchmark (Figures 3A,B). A range of models was constructed for selecting circRNAs. The best model (AUC value = 0.923) incorporating eight circRNAs (circ_30029, circ_117300, circ_176436, circ_112897, circ_178252, circ_115617, circ_14736, and circ_17720) was identified.
Figure 3. Construction and assessment of the eight circRNAs signature associated with biochemical recurrence by LASSO regression. (A) Process of variable selection in LASSO regression with 10-fold cross-validation. (B) Confidence interval in every lambda of LASSO regression. (C–E) circRNA predictor-score analysis of 144 PCa patients. The horizontal axis represents the 144 patients. Heatmap of circRNA expression level. (F) Biochemical recurrence analysis for the classifier. (G) Time-dependent ROC analysis curve for the classifier and clinical factors. (H) Kaplan-Meier estimates of eight circRNAs signature associated with biochemical recurrence.
Non-zero coefficients and the expression of eight circRNAs in the LASSO regression model constructed the following risk-score formula for BCR prediction: Risk score = (-0.0637∗expression level of circ_30029) + (-0.0991 ∗expression level of circ_117300) + (-0.0015∗expression level of circ_176436) + (-0.0185∗expression level of circ_112897) + (-0.0187∗expression level of circ_17720) + (0.0809 ∗expression level of circ_115617) + (0.0027 ∗expression level of circ_14736) + (0.0906 ∗expression level of circ_178252). The essential information and Kaplan-Meier plots of eight circRNAs were shown in Table 2 and Figure 3H. Three circRNAs with positive coefficients (circ_178252, circ_115617, and circ_14736) laid out high expression with an abhorrent outcome. Meanwhile, the negative coefficients in the remnants (circ_30029, circ_117300, circ_176436, circ_112897, and circ_17720) indicated high expression was correlated with less biochemical recurrence.
Table 2. CircRNAs significantly associated with the biochemical recurrence free survival in the test series patients (N = 144).
Figures 3C–E shows the relation among prognostic scores, the biochemical recurrence status, and prostate cancer circRNAs expression in 144 patients ranked by the prognostic score of the eight-circRNAs signature. Patients were divided into high-risk score group and low-risk score group by the median prognostic risk score. The higher risk scores revealed more BCR and less BCR-free time in PCa patients. Kaplan-Meier survival curve shows that patients with high-risk scores have statistically higher BCR than that of the low-risk score group (Figure 3F). As shown in Figure 3G, time-dependent ROC curve analyses were conducted to illustrate the sensitivity and specificity of BCR prediction. The AUC for the eight circRNAs signature prognostic model (AUC = 0.799) was calculated and compared with the AUC of clinical factors (AUC of PSA = 0.557, AUC of Gleason score = 0.626 and AUC of pathological stage = 0.569). Then, multivariate Cox regression analysis was carried out to estimate the prognostic value of the signature on BCR-free survival. With Gleason score, PSA, and pathological stage as covariates, the result shows that the eight-circRNAs signature was an independent prognostic factor for prostate cancer patients (Table 3). What’s more, we constructed a nomogram to predict the 5-year BCR-free survival rate of prostate cancer patients by the risk score and clinical factors, and the C-index of nomogram was 0.796 (Figure 4).
Figure 4. The nomogram to predict likelihood of 5-year BCR-free survival. Nomogram predicting PCa patients’ 5-year BCR-free survival. Each variable is projected to the point scale to obtain a value. The higher the value, the more likely to have BCR.
The Link Between 8 CircRNA’s Signature and Tumor Microenvironment
We employed DESeq2 to work out the differentially expressed mRNAs between the high risk score and low risk score groups. 338 differentially expressed mRNAs were enriched (logFC > 1 or logFC < -1, adjust P < 0.05), including 261 up-regulated mRNAs and 77 down-regulated mRNAs (Figure 5A). We presented the top 20 enriched GO terms (P < 0.05) by GO analyses involving biological process (BP), molecular function (MF), and cellular component (CC) in Figure 5B, which also includes the positive regulation of T cell activation via T cell receptor contact with antigens bound to MHC molecules on antigen presenting cells, MHC class II protein complex assembly, and MHC protein complex assembly.
Figure 5. Tumor immune microenvironment between low-risk and high-risk score group. (A) Volcano Plot shows the differentially expressed mRNAs between high and low risk score groups (The cutoff value is —logFC— > 1 and p < 0.05). (B) GO pathways analysis of differentially expressed mRNAs. (C) The proportion of 22 kinds of immune cells in individual. (D) The proportion of 22 immune cells between high and low risk score groups. The red means high risk score and the green represents low risk score. The parameter p-value was obtained by Wilcoxon test. (E) Kaplan-Meier analysis of BCR and the comparison of different risk-score group based on immune score, stromal score, and estimate score from ESTIMATE package.
Emerging evidence revealed that circRNAs related to the tumor microenvironment, involved in roles such as tumor surveillance, endothelial monolayer permeability, angiogenesis, hypoxia, remodeling of the extracellular matrix (ECM), exosomes, and so on (Ou et al., 2019; Hu et al., 2020; Huang et al., 2020). Based on the results of the GO analysis, we explored our eight circRNAs signature in the tumor immune microenvironment. CIBERSORT was employed to estimate the abundances of 22 immune cell types using mRNA expression data in GSE 113124 (Figure 5C). In Figure 5D, patients in the low-risk group had more CD8+ T cells (P < 0.001), plasma cells (P <0.01), and dendritic cells resting (P <0.01), and patients in the high-risk group had more activated dendritic cells (P <0.05). ESTIMATE R package was used to obtain immune score, stromal score, and ESTIMATE score of patients. According to medians of the three scores, patients were divided into a high score group and low score group. KM curves describe the relationships between biochemical recurrence and immune score, stromal score, and ESTIMATE score (Figure 5E). In ESTIMATE evaluation, it has a consistent trend with CIBERSORT that shows a high-risk score is positively correlated with high immune infiltration, although there is no significant statistical difference.
Experimental Verification of Circ_14734 and Circ_17720
According to a pervious study, circ_14734 and circ_17720 can be detected in PCa patients’ urine (Vo et al., 2019). We wondered if they had potential to be biomarkers. What’s more, we further characterized these two circRNAs by experimental analysis. Firstly, we specially designed the divergent primers targeted to each selected circRNA and identified the circularized site through Sanger sequencing (Figure 6A). Then we examined the basal expression of circ_14734 and circ_17720 through RT-qPCR, and found both to be downregulated in PCa cell lines (Figure 6B). Circ_14734 and circ_17720 were all resistant to RNase R digestion, while the linear mRNA of their parent genes were not (Figure 6C). Additionally, we found that circ_14734 and circ_17720 could be amplified only in cDNA instead of in gDNA (Figure 6D), which indicates that circ_14734 and circ_17720 were not the results of trans-splicing or genomic rearrangements (Jeck and Sharpless, 2014). Moreover, to further obtain the cell distribution of circ_14734 and circ_17720, cytoplasmic and nuclear fraction assay and RT-qPCR were conducted, and the results showed that circ_14734 and circ_17720 were mainly located in the cytoplasm of DU145 and PC3 cells (Figure 6E). Finally, to explore the function of circ_14734 and circ_17720, oligo RNAs containing si-NC and siRNAs specially targeting to the back-splicing junction site of circ_14734 and circ_17720 were transfected into DU145 and PC3 cells. It was confirmed that both first two siRNAs that targeted circ_14734 and circ_17720 efficiently knocked down their expression at around 70% via RT-qPCR. Then, CCK-8 assay was applied to detect the cell viability of NC and knockdown groups. As shown in Figure 6F, knockdown of circ_14734 or circ_17720 both inhibited the proliferation of DU145 and PC3 cells.
Figure 6. Experimental validation of circ_14736 and circ_17720. (A) The back-splicing junction site of circ_14734 and circ_17720 was determined by sanger sequencing. (B) The naive expression of circ_14734 and circ_17720 in prostate cancer cell lines. (C) circ_14734 and circ_17720 could be detected after the treatment of RNase R digestion. (D) the amplification of circ_14734 and circ_17720 by convergent and divergent primers was examined in cDNA and gDNA. (E) The distribution of circ_14734 and circ_17720 by nuclear and cytoplasmic fraction was determined by RT-qPCR. (F) Si-RNAs knock down effect screening. circ_14734 and circ_17720 inhibited DU145 and PC3 cell proliferation via CCK-8 assay. (G,H) The morphology and particle size distribution of Plasma-derived exosomes were verified by TEM and Nano-Sight. (I) Western blot analysis of exosomes isolated from three PCa patients’ blood plasma. The exosomal markers CD63 and HSP70 were detected. J. circ_14734 and circ_17720 were examined in blood plasma exosomes by gel electrophoresis.
Considering circ_14734 and circ_17720 were closely related to immune infiltration, which was regarded as a sort of tumor microenvironment, we began to consider whether there was a potential connection between tumor-cell-derived exosomes and the immune microenvironment. Therefore, we collected exosomes from serum in the manner of precipitation with PEG 8000 and verified the morphology and particle size through TME and NanoSight. The results showed that the exosomes were nearly spherical, and the particle size mainly distributed in around 100 nm (Figures 6G,H). CD63 and HSP70, two exosomal markers, were detected by western blot in three PCa patients’ blood plasma (Figure 6I). Subsequently, RT-qPCR and Agarose gel electrophoresis were performed to testify the existence of circ_14734 and circ_17720, and the results revealed that both circ_14734 and circ_17720 could be detected in exosomes derived from three PCa patients (Figure 6J).
Discussion
PCa remains one of the most common malignant tumors in men (Welch and Albertsen, 2020). Generally, BCR is a useful clinical index to monitor disease progression (Fenton et al., 2018). As reported by the European guidelines, BCR is characterized as two consecutive PSA ≥ 0.2 ng/ml after surgical procedure (Mottet et al., 2020). Moreover, some BCR cases are going to evolve into metastasis, which means a worse prognosis (Dall’Era et al., 2012). Due to the complex molecular mechanism of development and progression, the existing methods of PCa prognosis prediction are still deficient. We need to develop some prediction models for early detection, individualized therapeutic strategy, and distinguishing high-risk patients.
Molecular signatures have been related to predict the prognosis in some kinds of tumors. The risk score model based on non-coding RNAs, such as miRNAs and lncRNAs, has been largely reported in PCa (Pang et al., 2019). Through miRNA analysis in 78 PCa blood samples, Hoey et al. (2019) identified four circulating miRNAs which could accurately stratify patients into high- and low-risk categories after radical prostatectomy. In addition, Shao et al. (2019) built up a prognostic model using seven lncRNA signatures and revealed that the high risk of this model predicted faster BCR of PCa. However, the prognosis value of circRNAs in PCa is still not clear. Compared with other non-coding RNAs, circRNAs have a typical circular structure that resist the degradation of RNase and could be easier detected in body fluid (Su et al., 2019). There is growing research about the dysregulation of circRNAs in cancer. For example, Chen et al. (2020) found that circ-MALAT could directly bind to ribosome and PXA5 mRNA, forming a ternary complex, and functioned as a miRNAs sponge to sustain the self-renewal of hepatocellular cancer stem cells. Additionally, it has been reported that circASPA1 promoted hepatocellular carcinoma (HCC) progression by regulating miR-326/miR-532-5p-MAPK1 signaling, and positively correlated with tumor-associated-macrophage infiltration via the miR-326/miR-532-5p-CSF-1 axis (Hu et al., 2019). Moreover, circRNAs also exert an influence on metabolic reprogramming in malignant tumors. ciRS-122, which was confirmed to be an exosome-transmitted circRNA, could promote glycolysis to strengthen chemoresistance colorectal cancer cells (Wang et al., 2020).
Our study is the first to report a circRNAs signature model in predicting the BCR of PCa. In our work, the mRNAs and circRNAs obtained based on high-throughput sequencing and clinical data of 144 PCa patients are from GSE113124. Firstly, 28 circRNAs whose FPKM was more than 1.0 and were significantly linked to biochemical recurrence in logistic regression were selected. Subsequently, LASSO regression and 10-fold cross-validation targeted eight circRNAs to construct a best risk score model (AUC value = 0.923) for predicting BCR. Our eight circRNAs risk score model (circ_30029, circ_117300, circ_176436, circ_112897, circ_178252, circ_115617, circ_14736, and circ_17720) shows that patients with high-risk scores are statistically significantly higher than the low-risk score group in the rate of BCR by KM survival curve. What’s more, in the ROC and multivariate Cox regression analysis, the prognostic value of our risk score model is better than the existing clinical indexes.
Among the eight circRNAs in our model, only a few have been reported in medical research. Circ_117300 may be specifically found in tuberculosis patients’ plasma (Huang et al., 2018). In addition, circ_176436 inhibited ovarian cancer progression by suppressing miR-518a-5p to induce Fas expression (Zhang et al., 2020). To assess their function in PCa, we gathered the differentially expressed mRNAs between high and low risk score groups. Some T cell and MHC related pathways have been enriched in the GO analysis. As we know, T cells and MHC are two important parts of the immune system, so we hypothesized that these eight circRNAs could affect the tumor microenvironment in prostate cancer. As per the results of our following CIBERSORT and ESTIMATE analyses, we proved our assumption. In the CIBERSORT analysis, CD8+ T cells were markedly decreased in the high-risk score group. CD8+ T cells act as cytotoxic T cells to anti-tumor in most conditions (Qiao et al., 2019; Lynn et al., 2020). For our study, the eight circRNAs may promote BCR by protecting PCa cells as a result of cut down of the CD8+ T cells in the high-risk score group. Coincidentally, we have a consistent trend that the high-risk score is positively correlated with the high immune infiltration in ESTIMATE analysis, although there are no significant statistical differences.
To verify the individual capacity of these circRNAs, we selected circ_14736 and circ_17720 in our experimental validation because they exist in urine. We identified their junction sites by Sanger sequencing, examined the native expression in prostate cell lines, and characterized their circularized structure. What’s more, circ_14736 and circ_17720 showed cell viability inhibition effects in CCK8 assay. Finally, we detected them in the exosomes of PCa patients’ plasma. It gave us some hints that we could develop their prognosis values by blood test. To date, only a minor part of circRNAs’ biological functions have been identified, with most of them playing a miRNA sponge role (Hansen et al., 2013; Memczak et al., 2013). In addition, immunoprecipitation experiments data showed that circRNAs can interact with proteins directly (Ashwal-Fluss et al., 2014; Li Z. et al., 2015; Du et al., 2016; Abdelmohsen et al., 2017; Chen et al., 2018), through roles such as decoy proteins, enhancing protein function, protein scaffolding, and recruiting protein. Interestingly, some circRNAs have a protein or polypeptide translation ability (Legnini et al., 2017; Yang et al., 2018; Zhang et al., 2018), although most circRNAs are regarded as non-coding. In our data, circ_14736 and circ_17720 are mainly located in cytoplasm. What’s more, they both have a binding site with miRNAs and proteins according to the circinteractom database. Furthermore, only circ_17720 revels a protein coding potential (Fickett score: 1.1294, Hexamer score: 0.1549, IRES elements score: 0.90). We will explore the function and mechanism of circ_14736 and circ_17720 in our next study.
Unfortunately, due to the limitation of our clinical samples’ amount and patients’ prognosis information, we cannot verify our model in different datasets. We will monitor any updated circRNAs datasets to improve the accuracy of our model. We are also going to collect more samples to develop their prognosis value in blood.
Conclusion
We constructed an eight circRNAs risk score model to reliably predict the BCR of PCa patients. We found that the BCR predicting effect may be related to the tumor microenvironment. At the same time, we preliminarily verified the function of circ_14736 and circ_17720 in vitro. Further experiments are necessary to clarify their roles in PCa.
Data Availability Statement
The original contributions presented in the study are included in the article/supplementary material, further inquiries can be directed to the corresponding author/s.
Ethics Statement
The studies involving human participants were reviewed and approved by the Ethics Committee of Zhujiang Hospital, Southern medical University. The patients/participants provided their written informed consent to participate in this study.
Author Contributions
JL and XM designed the study. SW and WS analyzed the data. CZ did the experimental validation. TY, WC, GC, ZL, KW, and WZ searched literature. JL and BL wrote the manuscript. All authors read and approved the manuscript.
Funding
This research was supported by grants from the National Natural Science Foundation of China (81773277, 81802567, and 82003271), China postdoctoral science foundation (2019M662979), Science and Technology Program of Guangzhou (201803010014). and Fundamental Research Funds for the Central Universities (21618306).
Conflict of Interest
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
Supplementary Material
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fcell.2020.599494/full#supplementary-material
Abbreviations
AUC, Area Under Curve; BCR, Biochemical recurrence; BP, Biological Process; CC, Cellular Component; cDNA, Complementary DNA; circRNAs, Circular RNAs; DMEM, Dulbecco’s Modified Eagle’s medium; DU145, Duke University 145; ECM, Extracellular matrix; EMT, Epithelial-mesenchymal transition; FBS, Fetal bovine serum; FPKM, Fragments Per Kilobase Million; gDNA, Genomic DNA; GEO, Gene Expression Omnibus; GO, Gene Oncology; HCC, Hepatocellular carcinoma; HGNA, HUGO Gene Nomenclature Committee; KM survival curve, Kaplan-Meier survival curve; LASSO, Least absolute shrinkage and selection operator; LNCaP, Lymph Node Carcinoma of the prostate; lncRNA, Long non-coding RNA; MF, Molecular Function; MHC, Major histocompatibility complex; miRNA, MicroRNA; PCa, Prostate cancer; PSA, Prostate-specific antigen; Rnase R, Ribonuclease R; ROC, Receiver-operating characteristic; RP, Radical prostatectomy; RT-qPCR, Quantitative reverse transcription PCR; siRNA, Small interference RNAs; TIICs, Tumor-infiltrating immune cell; TME, Tumor microenvironment.
Footnotes
References
Abdelmohsen, K., Panda, A. C., Munk, R., Grammatikakis, I., Dudekula, D. B., De, S., et al. (2017). Identification of HuR target circular RNAs uncovers suppression of PABPN1 translation by CircPABPN1. RNA Biol. 14, 361–369. doi: 10.1080/15476286.2017.1279788
Ashwal-Fluss, R., Meyer, M., Pamudurti, N. R., Ivanov, A., Bartok, O., Hanan, M., et al. (2014). circRNA biogenesis competes with pre-mRNA splicing. Mol. Cell 56, 55–66. doi: 10.1016/j.molcel.2014.08.019
Celis, J. E., Gromov, P., Cabezon, T., Moreira, J. M., Friis, E., Jirstrom, K., et al. (2008). 15-prostaglandin dehydrogenase expression alone or in combination with ACSM1 defines a subgroup of the apocrine molecular subtype of breast carcinoma. Mol. Cell. Proteom. 7, 1795–1809. doi: 10.1074/mcp.R800011-MCP200
Chen, L., Kong, R., Wu, C., Wang, S., Liu, Z., Liu, S., et al. (2020). Circ-MALAT1 functions as both an mRNA translation brake and a microRNA sponge to promote self-renewal of hepatocellular cancer stem cells. Adv. Sci. 7:1900949. doi: 10.1002/advs.201900949
Chen, N., Zhao, G., Yan, X., Lv, Z., Yin, H., Zhang, S., et al. (2018). A novel FLI1 exonic circular RNA promotes metastasis in breast cancer by coordinately regulating TET1 and DNMT1. Genome Biol. 19:218. doi: 10.1186/s13059-018-1594-y
Chen, S., Huang, V., Xu, X., Livingstone, J., Soares, F., Jeon, J., et al. (2019). Widespread and functional RNA circularization in localized prostate cancer. Cell 176, 831–843.e22. doi: 10.1016/j.cell.2019.01.025
Dall’Era, M. A., Albertsen, P. C., Bangma, C., Carroll, P. R., Carter, H. B., Cooperberg, M. R., et al. (2012). Active surveillance for prostate cancer: a systematic review of the literature. Eur. Urol. 62, 976–983. doi: 10.1016/j.eururo.2012.05.072
de la Monte, S. M., Tamaki, S., Cantarini, M. C., Ince, N., Wiedmann, M., Carter, J. J., et al. (2006). Aspartyl-(asparaginyl)-beta-hydroxylase regulates hepatocellular carcinoma invasiveness. J. Hepatol. 44, 971–983. doi: 10.1016/j.jhep.2006.01.038
Ding, L., Abebe, T., Beyene, J., Wilke, R. A., Goldberg, A., Woo, J. G., et al. (2013). Rank-based genome-wide analysis reveals the association of ryanodine receptor-2 gene variants with childhood asthma among human populations. Hum. Genom. 7:16. doi: 10.1186/1479-7364-7-16
Du, W. W., Yang, W., Liu, E., Yang, Z., Dhaliwal, P., and Yang, B. B. (2016). Foxo3 circular RNA retards cell cycle progression via forming ternary complexes with p21 and CDK2. Nucleic Acids Res. 44, 2846–2858. doi: 10.1093/nar/gkw027
Feng, Y., Yang, Y., Zhao, X., Fan, Y., Zhou, L., Rong, J., et al. (2019). Circular RNA circ0005276 promotes the proliferation and migration of prostate cancer cells by interacting with FUS to transcriptionally activate XIAP. Cell Death Dis. 10:792. doi: 10.1038/s41419-019-2028-9
Fenton, J. J., Weyrich, M. S., Durbin, S., Liu, Y., Bang, H., and Melnikow, J. (2018). Prostate-specific antigen-based screening for prostate cancer: evidence report and systematic review for the US preventive services task force. JAMA 319, 1914–1931. doi: 10.1001/jama.2018.3712
Freedland, S. J., Humphreys, E. B., Mangold, L. A., Eisenberger, M., Dorey, F. J., Walsh, P. C., et al. (2005). Risk of prostate cancer-specific mortality following biochemical recurrence after radical prostatectomy. JAMA 294, 433–439. doi: 10.1001/jama.294.4.433
Glondu-Lassis, M., Dromard, M., Lacroix-Triki, M., Nirde, P., Puech, C., Knani, D., et al. (2010). PTPL1/PTPN13 regulates breast cancer cell aggressiveness through direct inactivation of Src kinase. Cancer Res. 70, 5116–5126. doi: 10.1158/0008-5472.CAN-09-4368
Hansen, T. B., Jensen, T. I., Clausen, B. H., Bramsen, J. B., Finsen, B., Damgaard, C. K., et al. (2013). Natural RNA circles function as efficient microRNA sponges. Nature 495, 384–388. doi: 10.1038/nature11993
Hoey, C., Ahmed, M., Fotouhi, G. A., Vesprini, D., Huang, X., Commisso, K., et al. (2019). Circulating miRNAs as non-invasive biomarkers to predict aggressive prostate cancer after radical prostatectomy. J. Transl. Med. 17:173. doi: 10.1186/s12967-019-1920-5
Hu, W., Liu, C., Bi, Z. Y., Zhou, Q., Zhang, H., Li, L. L., et al. (2020). Comprehensive landscape of extracellular vesicle-derived RNAs in cancer initiation, progression, metastasis and cancer immunology. Mol. Cancer 19:102. doi: 10.1186/s12943-020-01199-1
Hu, Z. Q., Zhou, S. L., Li, J., Zhou, Z. J., Wang, P. C., Xin, H. Y., et al. (2019). Circular RNA sequencing identifies CircASAP1 as a key regulator in hepatocellular carcinoma metastasis. Hepatology 72, 906–922. doi: 10.1002/hep.31068
Huang, X. Y., Huang, Z. L., Huang, J., Xu, B., Huang, X. Y., Xu, Y. H., et al. (2020). Exosomal circRNA-100338 promotes hepatocellular carcinoma metastasis via enhancing invasiveness and angiogenesis. J. Exp. Clin. Cancer Res. 39:20. doi: 10.1186/s13046-020-1529-9
Huang, Z., Su, R., Yao, F., Peng, Y., Luo, Q., and Li, J. (2018). Circulating circular RNAs hsa_circ_0001204 and hsa_circ_0001747 act as diagnostic biomarkers for active tuberculosis detection. Int. J. Clin. Exp. Pathol. 11, 586–594.
Jamshidi, J., Asnaashari, A., Alipoor, R., Mohammadi, S., Roostaei, S., Samadian, M. M., et al. (2018). ATP2B1 rs2681472 and STK39 rs35929607 polymorphisms and risk of hypertension in Iranian population. Med. J. Islam Repub. Iran. 32:14. doi: 10.14196/mjiri.32.14
Jeck, W. R., and Sharpless, N. E. (2014). Detecting and characterizing circular RNAs. Nat. Biotechnol. 32, 453–461. doi: 10.1038/nbt.2890
Kong, B., Wang, Q., Fung, E., Xue, K., and Tsang, B. K. (2014). p53 is required for cisplatin-induced processing of the mitochondrial fusion protein L-Opa1 that is mediated by the mitochondrial metallopeptidase Oma1 in gynecologic cancers. J. Biol. Chem. 289, 27134–27145. doi: 10.1074/jbc.M114.594812
Kreuz, M., Otto, D. J., Fuessel, S., Blumert, C., Bertram, C., Bartsch, S., et al. (2020). ProstaTrend-A multivariable prognostic RNA expression score for aggressive prostate cancer. Eur. Urol. 78, 452–459. doi: 10.1016/j.eururo.2020.06.001
Kristensen, L. S., Andersen, M. S., Stagsted, L., Ebbesen, K. K., Hansen, T. B., and Kjems, J. (2019). The biogenesis, biology and characterization of circular RNAs. Nat. Rev. Genet. 20, 675–691. doi: 10.1038/s41576-019-0158-7
Laczmanska, I., Karpinski, P., Gil, J., Laczmanski, L., Makowska, I., Bebenek, M., et al. (2017). The PTPN13 Y2081D (T>G) (rs989902) polymorphism is associated with an increased risk of sporadic colorectal cancer. Colorectal. Dis. 19, O272–O278. doi: 10.1111/codi.13727
Legnini, I., Di Timoteo, G., Rossi, F., Morlando, M., Briganti, F., Sthandier, O., et al. (2017). Circ-ZNF609 is a circular rna that can be translated and functions in myogenesis. Mol. Cell 66, 22–37.e9. doi: 10.1016/j.molcel.2017.02.017
Li, W., Ji, W., Li, Z., He, K., Wang, Q., Chen, J., et al. (2015). Genetic association of ACSM1 variation with schizophrenia and major depressive disorder in the Han Chinese population. Am. J. Med. Genet. B Neuropsychiatr. Genet. 168B, 144–149. doi: 10.1002/ajmg.b.32291
Li, Z., Huang, C., Bao, C., Chen, L., Lin, M., Wang, X., et al. (2015). Exon-intron circular RNAs regulate transcription in the nucleus. Nat. Struct. Mol. Biol. 22, 256–264. doi: 10.1038/nsmb.2959
Li, Z., Zhu, W., Xiong, L., Yu, X., Chen, X., and Lin, Q. (2016). Role of high expression levels of STK39 in the growth, migration and invasion of non-small cell type lung cancer cells. Oncotarget 7, 61366–61377. doi: 10.18632/oncotarget.11351
Luu, M., Sabo, E., de la Monte, S. M., Greaves, W., Wang, J., Tavares, R., et al. (2009). Prognostic value of aspartyl (asparaginyl)-beta-hydroxylase/humbug expression in non-small cell lung carcinoma. Hum. Pathol. 40, 639–644. doi: 10.1016/j.humpath.2008.11.001
Lynn, G. M., Sedlik, C., Baharom, F., Zhu, Y., Ramirez-Valdez, R. A., Coble, V. L., et al. (2020). Peptide-TLR-7/8a conjugate vaccines chemically programmed for nanoparticle self-assembly enhance CD8 T-cell immunity to tumor antigens. Nat. Biotechnol. 38, 320–332. doi: 10.1038/s41587-019-0390-x
Memczak, S., Jens, M., Elefsinioti, A., Torti, F., Krueger, J., Rybak, A., et al. (2013). Circular RNAs are a large class of animal RNAs with regulatory potency. Nature 495, 333–338. doi: 10.1038/nature11928
Mitra, R., Adams, C. M., Jiang, W., Greenawalt, E., and Eischen, C. M. (2020). Pan-cancer analysis reveals cooperativity of both strands of microRNA that regulate tumorigenesis and patient survival. Nat. Commun. 11:968. doi: 10.1038/s41467-020-14713-2
Mo, R. J., Lu, J. M., Wan, Y. P., Hua, W., Liang, Y. X., Zhuo, Y. J., et al. (2017). Decreased HoxD10 expression promotes a proliferative and aggressive phenotype in prostate cancer. Curr. Mol. Med. 17, 70–78. doi: 10.2174/1566524017666170220104920
Moris, L., Cumberbatch, M. G., Van den Broeck, T., Gandaglia, G., Fossati, N., Kelly, B., et al. (2020). Benefits and risks of primary treatments for high-risk localized and locally advanced prostate cancer: an international multidisciplinary systematic review. Eur. Urol. 77, 614–627. doi: 10.1016/j.eururo.2020.01.033
Mottet, N., van den Bergh, R. C. N., Briers, E., Cornford, P., De Santis, M., Fanti, S., et al. (2020). EAU - ESTRO - ESUR - SIOG guidelines on prostate cancer - 2020 update. Part 1: screening, diagnosis, and local treatment with curative intent. Eur Urol.
Ou, Z. L., Luo, Z., Wei, W., Liang, S., Gao, T. L., and Lu, Y. B. (2019). Hypoxia-induced shedding of MICA and HIF1A-mediated immune escape of pancreatic cancer cells from NK cells: role of circ_0000977/miR-153 axis. RNA Biol. 16, 1592–1603. doi: 10.1080/15476286.2019.1649585
Pang, X., Zhao, Y., Wang, J., Li, W., Xiang, Q., Zhang, Z., et al. (2019). Competing endogenous RNA and coexpression network analysis for identification of potential biomarkers and therapeutics in association with metastasis risk and progression of prostate cancer. Oxid. Med. Cell. Longev. 2019:8265958. doi: 10.1155/2019/8265958
Qiao, J., Liu, Z., Dong, C., Luan, Y., Zhang, A., Moore, C., et al. (2019). Targeting Tumors with IL-10 prevents dendritic cell-mediated CD8(+) T cell Apoptosis. Cancer Cell 35, 901–915.e4. doi: 10.1016/j.ccell.2019.05.005
Ramsay, I. S., Ma, S., Fisher, M., Loewy, R. L., Ragland, J. D., Niendam, T., et al. (2018). Model selection and prediction of outcomes in recent onset schizophrenia patients who undergo cognitive training. Schizophr. Res. Cogn. 11, 1–5. doi: 10.1016/j.scog.2017.10.001
Ren, A., Sun, S., Li, S., Chen, T., Shu, Y., Du, M., et al. (2019). Genetic variants in SLC22A3 contribute to the susceptibility to colorectal cancer. Int. J. Cancer 145, 154–163. doi: 10.1002/ijc.32079
Salzman, J., Chen, R. E., Olsen, M. N., Wang, P. L., and Brown, P. O. (2013). Cell-type specific features of circular RNA expression. PLoS Genet. 9:e1003777. doi: 10.1371/journal.pgen.1003777
Santer, L., Bar, C., and Thum, T. (2019). Circular RNAs: a novel class of functional RNA molecules with a therapeutic perspective. Mol. Ther. 27, 1350–1363. doi: 10.1016/j.ymthe.2019.07.001
Shao, N., Zhu, Y., Wan, F. N., and Ye, D. W. (2019). Identification of seven long noncoding RNAs signature for prediction of biochemical recurrence in prostate cancer. Asian J. Androl. 21, 618–622. doi: 10.4103/aja.aja_118_18
Shen, Z., Zhou, L., Zhang, C., and Xu, J. (2020). Reduction of circular RNA Foxo3 promotes prostate cancer progression and chemoresistance to docetaxel. Cancer Lett. 468, 88–101. doi: 10.1016/j.canlet.2019.10.006
Shi, J., Liu, C., Chen, C., Guo, K., Tang, Z., Luo, Y., et al. (2020). Circular RNA circMBOAT2 promotes prostate cancer progression via a miR-1271-5p/mTOR axis. Aging 12, 13255–13280. doi: 10.18632/aging.103432
Siegel, R. L., Miller, K. D., and Jemal, A. (2020). Cancer statistics, 2020. CA Cancer J. Clin. 70, 7–30. doi: 10.3322/caac.21590
Slack, F. J., and Chinnaiyan, A. M. (2019). The role of non-coding RNAs in oncology. Cell 179, 1033–1055. doi: 10.1016/j.cell.2019.10.017
Su, M., Xiao, Y., Ma, J., Tang, Y., Tian, B., Zhang, Y., et al. (2019). Circular RNAs in Cancer: emerging functions in hallmarks, stemness, resistance and roles as potential biomarkers. Mol. Cancer 18:90. doi: 10.1186/s12943-019-1002-6
Sun, X., Jung, J. H., Arvola, O., Santoso, M. R., Giffard, R. G., Yang, P. C., et al. (2019). Stem Cell-derived exosomes protect astrocyte cultures from in vitro ischemia and decrease injury as post-stroke intravenous therapy. Front. Cell. Neurosci. 13:394. doi: 10.3389/fncel.2019.00394
Vo, J. N., Cieslik, M., Zhang, Y., Shukla, S., Xiao, L., Zhang, Y., et al. (2019). The landscape of circular RNA in cancer. Cell 176, 869–881.e13. doi: 10.1016/j.cell.2018.12.021
Waerner, T., Gardellin, P., Pfizenmaier, K., Weith, A., and Kraut, N. (2001). Human RERE is localized to nuclear promyelocytic leukemia oncogenic domains and enhances apoptosis. Cell Growth Differ. 12, 201–210.
Wang, X., Zhang, H., Yang, H., Bai, M., Ning, T., Deng, T., et al. (2020). Exosome-delivered circRNA promotes glycolysis to induce chemoresistance through the miR-122-PKM2 axis in colorectal cancer. Mol. Oncol. 14, 539–555. doi: 10.1002/1878-0261.12629
Welch, H. G., and Albertsen, P. C. (2020). Reconsidering prostate cancer mortality - the future of PSA screening. N. Engl. J. Med. 382, 1557–1563. doi: 10.1056/NEJMms1914228
Yang, H., Liu, H., Lin, H. C., Gan, D., Jin, W., Cui, C., et al. (2019). Association of a novel seven-gene expression signature with the disease prognosis in colon cancer patients. Aging 11, 8710–8727. doi: 10.18632/aging.102365
Yang, Y., Gao, X., Zhang, M., Yan, S., Sun, C., Xiao, F., et al. (2018). Novel role of FBXW7 circular RNA in repressing glioma tumorigenesis. J. Natl. Cancer Inst. 110:435. doi: 10.1093/jnci/djx166
Yoshihara, K., Shahmoradgoli, M., Martinez, E., Vegesna, R., Kim, H., Torres-Garcia, W., et al. (2013). Inferring tumour purity and stromal and immune cell admixture from expression data. Nat. Commun. 4:2612. doi: 10.1038/ncomms3612
Zhang, M., Zhao, K., Xu, X., Yang, Y., Yan, S., Wei, P., et al. (2018). A peptide encoded by circular form of LINC-PINT suppresses oncogenic transcriptional elongation in glioblastoma. Nat. Commun. 9:4475. doi: 10.1038/s41467-018-06862-2
Zhang, N., Jin, Y., Hu, Q., Cheng, S., Wang, C., Yang, Z., et al. (2020). Circular RNA hsa_circ_0078607 suppresses ovarian cancer progression by regulating miR-518a-5p/Fas signaling pathway. J. Ovarian Res. 13:64. doi: 10.1186/s13048-020-00664-1
Keywords: prostate cancer, circRNA, biochemical recurrence, immune infiltration, LASSO
Citation: Wang S, Su W, Zhong C, Yang T, Chen W, Chen G, Liu Z, Wu K, Zhong W, Li B, Mao X and Lu J (2020) An Eight-CircRNA Assessment Model for Predicting Biochemical Recurrence in Prostate Cancer. Front. Cell Dev. Biol. 8:599494. doi: 10.3389/fcell.2020.599494
Received: 27 August 2020; Accepted: 16 November 2020;
Published: 10 December 2020.
Edited by:
Kexin Xu, The University of Texas Health Science Center at San Antonio, United StatesReviewed by:
Ping Zhang, Wellcome Trust Centre for Human Genetics (WT), United KingdomAvisek Majumder, University of California, San Francisco, United States
Zhaoyu Liu, Fifth Affiliated Hospital of Guangzhou Medical University, China
Copyright © 2020 Wang, Su, Zhong, Yang, Chen, Chen, Liu, Wu, Zhong, Li, Mao and Lu. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Bingkun Li, cGs5Mzk0MjhAMTYzLmNvbQ==; Xiangming Mao, bXhtNjMxMjIxQDEyNi5jb20=; Jianming Lu, bG91aXNjZmM4QGdtYWlsLmNvbQ==
†These authors have contributed equally to this work