A two‑circRNA signature predicts tumour recurrence in clinical non‑functioning pituitary adenoma.

Clinical non‑functioning pituitary adenoma (NFPA) accounts for >30% of all pituitary adenomas, and the recurrence rate is notably high. The ability to predict tumour recurrence during initial surgery will aid in determining if adjunctive therapy is required to reduce recurrence. With the aim of developing a circular RNA (circRNA) signature to improve prognosis prediction in NFPA, the present study examined the circRNA expression profiles in 73 patients with NFPA from Beijing Tiantan Hospital using high‑throughput RNA chip technology. The dataset was randomly separated into a training group and a test group using an R program. In the training group (n=37), a Cox proportional hazards regression model was used to analyse the genes associated with the recurrence and progression‑free survival (PFS) of patients with NFPA. Meanwhile, a random survival forest algorithm, Kaplan‑Meier and receiver operating characteristic curve (ROC) analyses were used to determine the multi‑circRNA signature with the largest area under the ROC curve (AUROC) and verify its efficacy in the test group (n=36). In the training and test groups, the signatures of two circRNAs (hsa_circ_0000066 and hsa_circ_0069707) were specifically associated with the PFS of patients with NFPA (log‑rank P<0.05). Furthermore, the two‑circRNA signature had a high prediction accuracy for tumour recurrence, with an AUROC of 0.87 and 0.67 in the training and test groups, respectively; and the discriminative power of the signature was greater compared with that of age. The present study is the first to suggest a circRNA signature with a clinical application value for predicting recurrence/progression in patients with NFPA.

Introduction

Pituitary adenomas (PAs) are clinically common intracranial neoplasms (1). Non-functioning pituitary adenomas (NFPAs) account for nearly half of all PAs (2). The majority of NFPAs are detected due to neurological manifestations, including headaches, vision and visual field disorders or incidental imaging examinations (3). The first-line treatment is surgical resection, but the clinical efficacy depends on the tumour size in addition to the surgeon's ability (4,5) as numerous macroadenomas invade important surrounding structures, including the cavernous sinus, internal carotid artery and optic chiasm (6). This invasion makes it difficult to achieve total resection. Following partial tumour resection, any residual tumour are a potential source of regeneration (7). Even certain tumour types that are completely resected during the initial surgery still have a risk of regrowth (8). For non-functioning pituitary macroadenoma tumour types, the recurrence rate may range from 12 to 58% (9). A previous study demonstrated that NFPA recurrence and regrowth was likely associated with the biomarkers high-mobility group protein HMG-1/HMG-Y and mouse double minute 2 homolog (10). In a separate study, fascin actin-bundling protein 1 was suggested to increase the risk of PA recurrence (11). However, the molecular mechanisms of tumour regrowth and recurrence are unclear and require further research.

Postoperative radiation is considered an effective adjunctive therapy that is able to reduce recurrence. However, it is often difficult to determine whether radiation therapy is a reasonable therapeutic method for a patient, as it may result in progressive hypopituitarism and other long-term complications (12). Therefore, radiotherapy exclusively benefits patients with a high risk of tumour recurrence following surgery, and its use is controversial in patients without definite evidence of tumour regrowth (13). A method for predicting recurrence during the initial surgery will help to determine if adjunctive radiation therapy is required.

Circular RNAs (circRNAs) are covalently closed RNAs that have been identified to regulate mammalian gene expression and may be novel markers for cancer diagnosis (14–18). The levels of hsa_circ_100269 are decreased in gastric cancer (GC) tissues, and a signalling pathway involving hsa_circ_100269 and hsa_miR_630 serves a critical role in GC cell growth (19). hsa_circ_001988 is downregulated in colorectal cancer, and its expression is associated with differentiation and perineural invasion (20). The expression levels of hsa_circ_0001649 are downregulated in hepatocellular carcinoma and colorectal cancer, indicating that hsa_circ_0001649 may be a novel biomarker (21,22). Additionally, circRNAs serve a strong regulatory function in carcinoma; for example, they are able to function as efficient microRNA (miRNA/miR) sponges (15,23,24). circRNAs sponge miRNAs through complementary sequences to reduce the function of the microRNAs (25). For instance, the circRNA cerebellar degeneration-related protein 1, which contains miR-7 binding sites, is able to reduce the growth of human cells by inhibiting miR-7 binding (26,27); and circRNA homeodomain interacting protein kinase 3 binds directly to miR-124 (28,29), thus regulating cancer development. In hepatocellular carcinoma, hsa_circ_0067934 accelerates tumour growth and metastasis by modulating the miR-1324/Frizzled class receptor 5/Wnt/β-catenin axis (30). In addition, hsa_circ_100290 sponges the miR-29 family in oral cancer (31). Accordingly, circRNAs may be novel potential cancer biomarkers and therapeutic targets. Thus, the identification of pituitary-associated circRNAs, in addition to their clinical roles and molecular mechanisms, is crucial for understanding the development and progression of NFPA.

In the present study, tumour regrowth from residual cells and tumour reappearance following total resection were collectively referred to as tumour recurrence. The aim was to develop a circRNA signature that was associated with recurrence and could be used for treatment guidance in patients with NFPA.

Materials and methods

Patients and samples

A total of 73 specimens were obtained from patients who underwent surgery at Beijing Tiantan Hospital (Beijing, China) between October 2007 and July 2014. The patients included 34 males and 39 females, with a median age of 55 years (range 25–73 years) and a median follow-up of 78 months (range 36–121 months). Tumors with Hardy-Wilson classification grade IV and/or Knosp classification grades III and IV (32,33) were defined as invasive. Table I includes the clinical and pathological features of the adenomas. All tissue samples were immediately snap-frozen and stored in liquid nitrogen or at −80°C. The tumour type was determined based on pathological diagnosis and the preoperative clinical and biochemical examination results. This study was ethically approved by the Tiantan Hospital Ethics Committee of Beijing Tiantan Hospital, and written informed consent was obtained from each patient prior to the study.

Table I.

Summary of patient demographics and clinical characteristics.

Characteristic	Training set (n)	Testing set (n)	Total (n)
Sex
Female	19	20	39
Male	18	16	34
Age, years
≤55	19	20	39
>55	18	16	34
Vision and visual field disorders			0
Yes	18	28	46
No	19	8	27
Headache
Yes	18	18	36
No	19	18	37
Hypogonadism
Yes	2	8	10
No	35	28	63
Knosp classification
I	6	5	11
II	9	9	18
III	8	8	16
IV	14	14	28
Hardy classification
I	0	2	2
II	15	2	17
III	7	10	17
IV	15	0	15
Pathological type
Gonadotrophin	12	9	21
Mult	4	1	5
Null cell	21	25	46
Silent GH	0	1	1
Invasion
Yes	25	23	48
No	12	13	25

Mult, clinical non-functioning pituitary adenomas with multiple hormones positive immunostaining; GH, growth hormones.

Total RNA isolation and circRNA microarray

Total RNA from a tumour was isolated and purified using a phenol-free mirVana™ miRNA Isolation kit (cat. no. AM1561; Ambion; Thermo Fisher Scientific, Inc., Waltham, MA, USA), according to the manufacturer's protocol, and the purity and quantity were measured using a Nanodrop spectrophotometer. Subsequently, these 73 total RNA samples (1.8 < Optical Density 260 < 2.1) were used to produce fluorescence-labelled cRNA targets for the SBC human ceRNA array V1.0 (4×180 K). An Agilent Microarray scanner (Agilent Technologies, Inc., Santa Clara, CA, US) was used to scan the slides and hybridize the labelled cRNA targets. Feature Extraction software 10.7 (Agilent Technologies, Inc.) was used to extract the data. The R 3.2.3 software (algorithm ‘mean’ and ‘var’; www.r-project.org/) was used to normalize the raw data and for data processing.

Selection of circRNA subsets associated with progression-free survival (PFS)

The present study used the R 3.2.3 (www.r-project.org/) program algorithm ‘sample’ to randomly divide the dataset into training (n=37) and test groups (n=36). The association between the continuous expression level of each circRNA and patient PFS in the training set was determined using univariable Cox regression analysis: A gene with a univariate Cox coefficient >0 was considered to be a risk factor, <0 was considered to be a protective factor and -log10 P>1.30 (-log100.05) was regarded as significant.

Selection of circRNA prognostic signature based on a risk score model

Next, a model for selecting circRNA signatures with prognostic indicators was developed as follows:

Risk score = ∑Ni = 1 (Explg*Coef) (30,34).

Where N is the number of predictive circRNAs, Explg indicates the expression of circRNAs and Coef is the evaluated regression coefficient of the circRNAs determined from the univariable Cox regression analysis.

circRNAs were screened using a random survival forest-variable hunting (RSFVH) algorithm, and then a total of 9 circRNAs were selected (30) as fewer circRNAs provide a greater value. All 9 of the circRNAs were analysed. The total of all 9 node permutations and combinations was 29−1=511 and the sensitivity and specificity of the survival prediction of the risk score of the 29−1=511 signatures in the training dataset were compared with time-dependent receiver operating characteristic (ROC) curves, and the area under the curve (AUROC) values were obtained. In the training test, the marginal value was obtained from the median risk score. According to the median risk score, the patients with NFPA were separated into a high-risk group and a low-risk group in each dataset.

Statistical analysis

The equality of survival distributions in different groups for the training and test NFPA cohorts was assessed using Kaplan-Meier survival analyses, and two-sided log-rank tests were used to evaluate the statistical significance. Furthermore, the association with clinical features was analysed via χ2 tests. The risk grade and other clinical features in the available data when independent were tested using multivariate Cox regression analysis and data stratification analysis. P<0.05 was considered to indicate a statistically significant difference. The R 3.2.3 program was used for the analysis of all data and to generate plots; and the pROC, survival, ggplot2 and random forest SRC packages were downloaded from Bioconductor (www.r-project.org) (35,36). The analysis sequence used for development of a risk scoring model and to verify the validity of signatures is presented in Fig. 1.

Figure 1.

Schedule of the present study. The order of analyses conducted to develop the risk score model and validate the efficiency of the signature to predict outcomes. CircRNA, circular RNA; ROC, receiver operating characteristic; AUC, area under curve.

Functional analysis of circRNAs with prognostic value

Pearson's correlation coefficients were used to determine the co-expression association between prognostic circRNAs and protein-coding genes. Gene ontology (GO) (37,38) and Kyoto Encyclopaedia Gene and Genome (KEGG) (39–41) enrichment analyses using the Cytoscape plugin ClueGo (version 3.2.3) were performed to determine the co-expressed protein-coding genes and prognostic circRNAs for predicting the biological function of prognostic circRNAs (42). Cytoscape is a functional annotation tool used for assessing the over-representation of a category of interest gene concentration. Enrichment analyses were performed using functional annotation maps and functional annotation clustering options and were restricted to the KEGG pathway and GO terminology in the ‘biological process’ category. Functional annotations with P<0.05 were considered to indicate a statistically significant difference.

Results

Selection of a prognostic circRNA signature with a high predictive power for recurrence progression

The present study examined the pituitary tumour tissues of 73 patients using an Agilent microarray. From this analysis, a total of 88,750 circRNAs were identified. Next, the patients were randomly divided into two groups: A training set (n=37) and a test set (n=36). The training set was used to select the circRNA prognostic signature. The test set was used to evaluate the discriminative power of the circRNA signature, which was obtained independently. In this model selection section, all analyses were based on the training set.

First, univariate Cox proportional hazards regression analysis of circRNA expression using recurrence as the dependent variable was performed. A total of 7,481 circRNAs (data not shown) were determined to be correlated (P<0.05) with patient recurrence/progression and a volcano plot was used to present the univariate Cox results (Fig. 2A). Next, a RSFVH algorithm was used to select 9 circRNAs that were most closely associated with progression/relapse in the set of 7,481 circRNAs, according to the replacement vital score of the RSFVH algorithm (Fig. 2B).

Figure 2.

Identifying the circRNA signature in the training dataset. (A) Univariate Cox proportional hazards regression analysis of circRNA expression profiling data in the training dataset. (B) Procedure of the random forest supervised classification algorithm. (C) Procedure for identifying the final signature. The accuracies of all 511 signatures were calculated, and the 9 highest accuracies for k=1–9 are presented in the plot. (D) Receiver operating characteristic curve for the circRNA signature prediction model in the training dataset. CircRNA, circular RNA; AUC, area under the curve; PFS, progression-free survival.

The application of time-dependent ROC curves and the comparison of the sensitivity and specificity of the risk score survival prediction may result in better predictive characteristics. The risk score is based on the risk score model of 29−1=511 combinations of the circRNA training set (data not shown) for each circRNA of each patient in the training dataset. All risk scores for the circRNA signatures were calculated using this method. Then, the circRNA combination of hsa_circ_0000066 and hsa_circ_0069707, which had the largest AUROC, was selected (Fig. 2C and D; Table II). A risk score for the combination of hsa_circ_0000066 and hsa_circ_0069707 was obtained as follows: Risk score=(2.38 × expression value of RP11-702B10.1) + (8.27 × expression value of hsa_circ_0069707). The AUROC for the circRNA signature in the model was 0.87 in the training set. Thus, the optimal survival prediction value was obtained based on the expression-based risk score.

Table II.

Identities of circRNAs in the prognostic expression signature and their univariable cox association with prognosis.

							Chromosome location
Database ID	Best transcript	Gene symbol	Genomic length	Coefficient[a]	P-value[a]	Gene expression level association with poor prognosis	(GRCh38/hg38)
hsa_circ_0000066[b]	NM_018150	RNF220	132	2.38	<0.001	High	chr1:44890310-44890442:-1
hsa_circ_0069707[b]	NM_015030	FRYL	402	8.27	<0.001	High	chr4:48712535-48782316:-1

Derived from the univariable Cox regression analysis in the training set.

CirBase database. CircRNA, circular RNA; RNF220, ring finger protein 200; FRYL, FRY like transcription coactivator.

Evaluation of the circRNA signature for predicting the recurrence and PFS of patients with NFPA

Subsequently, the predictive power of the circRNA signature was assessed. The test set (n=36) was utilized for an independent evaluation, and the training set (n=37) was evaluated in parallel for comparison.

In the training set, based on the median value of the risk score of the circRNA signature being used as the cutoff value, the patients were separated into a high-risk group (n=18) and a low-risk group (n=19). The PFS time was significantly shorter in the high-risk patients compared with the low-risk patients (median survival: 56 months vs. >80.2 months; P<0.001; Fig. 3A). The rate of no recurrence or progression in high-risk patients was <1% and that in low-risk patients was >90%.

Figure 3.

CircRNA signature to predict recurrence in and progression-free survival of patients with non-functioning pituitary adenomas. Kaplan-Meier survival curves of patients sorted into high- and low-risk groups using the circRNA signature in the (A) training and (B) test datasets. P-values were calculated using a log-rank test. The predictive power of the circRNA signature compared with that of age based on ROC curves in (C) training and (D) test datasets. CircRNA, circular RNA; ROC, receiver operating curves; AUROC, area under the ROC; CI, confidence interval.

In the independent test set, the patients were separated into two groups (n=18). Verification of the recurrence prediction ability of the circRNA signature was performed using Kaplan-Meier curves for the high-risk and low-risk groups in the test dataset; the data were plotted and are presented in Fig. 3B, revealing a significant difference in PFS between the two groups (median PFS time: 70.2 months vs. >60 months; log-rank sum test, P=0.047). The rate of no recurrence/progression was ~16.6% in the high-risk group and 50% in the low-risk group.

CircRNA signature is independent of clinicopathological features

The clinical significance of the circRNA signature in NFPAs was obtained by associating the signature with a range of clinical pathological parameters from the dataset (n=73). As presented in Table III, except for age, vision and visual field disorders, and headache, there was no significant association between the circRNA markers and clinical pathological variables, including sex and histological grade (Table III).

Table III.

Association of the circRNA signature with clinicopathological characteristics in patients with pituitary adenoma.

	Training group			Test group			Entire group
Variables	Low risk	High risk	P-value	Low risk	High risk	P-value	Low risk	High risk	P-value
Sex			0.62			>0.99			0.73
Female	11	8		10	10		21	18
Male	8	10		8	8		16	18
Age, years			0.41			0.09			0.05
≤55	8	11		7	13		15	24
>55	11	7		11	5		22	12
Vision and visual field disorders			>0.99			0.05			0.13
Yes	13	13		11	17		24	30
No	6	5		7	1		13	6
Headache			0.03			0.74			0.29
Yes	13	5		8	10		21	15
No	6	13		10	8		16	21
Hypogonadism			>0.99			0.30			0.34
Yes	1	1		2	5		3	6
No	18	17		15	13		33	30
Knosp classification			0.49			0.10			0.93
I	2	4		4	1		6	5
II	6	3		2	7		8	10
III	5	3		3	5		8	8
IV	6	8		9	5		15	13
Hardy classification			0.52			0.43			0.27
II	6	9		10	13		16	22
III	4	3		3	3		7	6
IV	9	6		5	2		14	8
Pathological type			0.48			0.47			0.49
Gonadotroph	4	8		4	5		8	13
Mult	2	2		1	0		3	2
Null cell	13	8		13	12		26	20
Silent GH				0	1		0	1
Invasion			>0.99			>0.99			0.93
Yes	13	12		12	11		25	23
No	6	6		6	7		12	13

circRNA, circular RNA; Mult, Clinical non-functioning pituitary adenomas with multiple hormones positive immunostaining; GH, growth hormones.

To assess whether the ability of the circRNA markers in predicting the progression of recurrence was independent of other clinical features, multivariable Cox regression analysis was performed using the NFPA dataset and a signature-based risk score based on two circRNAs and other clinical features. The analysis indicated that the circRNA marker risk score was independent of the clinical characteristics in predicting survival, but that the signature was significantly associated with survival (high-risk group vs. low-risk group; hazard ratio=10.07, 95% confidence interval, 2.92–34.77; P<0.001, n=73; Table IV).

Table IV.

Univariable and multivariable Cox regression analysis of the circRNA signature and survival of patients with non-functioning pituitary adenoma in the training, test group and entire group.

		Training set (n=37)				Test set (n=36)				Entire dataset (n=73)
			95% CI of HR				95% CI of HR				95% CI of HR
Variables	Comparison	HR	Lower	Upper	P-value	HR	Lower	Upper	P-value	HR	Lower	Upper	P-value
Univariable analysis
Signature	High-risk vs. low risk	2.36	1.59	3.51	<0.001	1.04	0.88	1.22	0.67	9.89	2.92	33.54	<0.001
Sex	Male vs. female	1.45	0.42	4.99	0.55	0.58	0.17	1.93	0.38	0.96	0.42	2.18	0.93
Age, years	>55 vs. ≤55	0.42	0.11	1.59	0.20	0.78	0.23	2.60	0.68	0.58	0.24	1.40	0.23
Knosp classification	I–II vs. III–IV	0.97	0.56	1.67	0.91	1.20	0.69	2.09	0.52	1.08	0.74	1.58	0.68
Multivariable analysis
Signature	High-risk vs. low risk	3.42	1.82	6.41	<0.001	1.08	0.88	1.33	0.45	10.07	2.92	34.77	<0.001
Sex	Male vs. female	0.40	0.03	4.69	0.47	0.50	0.13	1.84	0.29	0.80	0.34	1.87	0.61
Age, years	>55 vs. ≤55	3.79	0.52	27.75	0.19	0.94	0.27	3.30	0.92	0.80	0.32	2.01	0.64
Knosp classification	I–II vs. III–IV	1.77	0.66	4.79	0.26	1.29	0.71	2.37	0.40	1.14	0.75	1.73	0.54

CircRNA, circular RNA; HR, hazard ratio; CI, confidence interval.

Comparing the predictive power of the circRNA signature with age

In order to compare the sensitivity and specificity of age and the circRNA-specific survival predictions, ROC analysis was performed as a larger AUROC usually indicates a better predictive model. In the dataset (n=37), the prognostic power of the circRNA marker was significantly stronger compared with age for the training and test groups (AUROC=0.87 and 0.69 for signature, respectively, vs. AUROC=0.65 and 0.60 for age, respectively; Fig. 3C and D). These data further confirm that the circRNA marker is a novel predictive marker with a high accuracy and thus has major clinical value.

Functional characterization of the circRNAs in the prognostic signature

Pearson's correlation coefficient was used for 73 patients to calculate the co-expression association between the two circRNAs (hsa_circ_0000066 and hsa_circ_0069707) and protein-coding genes to determine the potential biological effects of these molecular markers. Of the two circRNAs, at least one may be associated with the expression levels of 623 protein-coding genes, and the interaction network was determined with Cytoscape (Pearson's correlation coefficient >0.50; P<0.05; Figs. 4 and 5). Co-expressed protein-coding genes with the whole human genome as the background were obtained via GO and KEGG pathway function enrichment analysis. A total of 623 protein-coding genes were enriched in 57 GO/KEGG terms according to GO functional annotation (data not shown). These important GO/KEGG terms were used in the ClueGo plugin in Cytoscape to organize an interaction network with similar functionality. Several groups of functionally associated GO terms were analysed in the present study and should be further investigated in the future. Based on the present analysis, the two circRNAs likely to influence protein-coding gene interactions and major biological processes, including vesicle transport and the cellular response to unfolded proteins, may be involved in tumorigenesis (Fig. 6).

Figure 4.

Co-expression of the selected circular RNA hsa_circ_0069707 with other functional genes visualized with Cytoscape (Pearson's correlation coefficient >0.50, P<0.05).

Figure 5.

Co-expression of the selected circular RNA hsa_circ_0000066 with other function genes visualized with Cytoscape (Pearson's correlation coefficient >0.50, P<0.05).

Figure 6.

Functional enrichment of the protein-coding genes co-expressed with the prediction circular RNAs using ClueGo. GO, Gene Ontology; KEGG, Kyoto Encyclopedia of Genes and Genomes; IRE1, inositol-requiring enzyme 1; ER, endoplasmic reticulum.

Discussion

NFPA does not induce endocrine symptoms via hormone hypersecretion, and thus differs from functioning PA. Surgical resection is the preferred treatment (3–5). Postoperative tumour recurrence is difficult to detect; regrowth is discovered only when the tumour is large enough to be observed in regular follow-up imaging examination or to cause notable symptoms resulting from a tumour mass effect (43). This situation may exacerbate the patient prognosis. Considering the severity of NFPA recurrence, numerous studies are in progress to identify novel predictors that may identify early recurrence (7). The purpose of the present study was mainly to separate patients into high-risk groups or low-risk group in order to identify the most effective and least toxic treatment method.

CircRNAs are closed-loop RNAs that were originally thought to be by-products of splicing or other mRNA processing errors. However, current evidence demonstrates that circRNA may be involved in disease development and even in the development of a tumour (44). The level of hsa_circRNA_100855 was revealed to be significantly (P<0.001) higher in laryngeal squamous cell cancer tissues compared with in the corresponding adjacent non-neoplastic tissues, and the levels of hsa_circRNA_104912 were significantly (P<0.001) lower (25). Filamin binding LIM protein 1 (FBLIM1) expression was regulated by circFBLIM1, which may serve as a competing endogenous RNA through the sponging of miR-346 in hepatocellular cancer (45). Zhang et al (46) established a four-circRNA-based classifier (hsa_circRNA_101308, hsacircRNA_104423, hsa_circRNA_104916 and hsa_circRNA_100269) to predict early recurrence for patients with stage III gastric cancer following radical surgery. Therefore, these results suggest that circRNAs may serve as novel diagnostic markers and treatment targets.

To the best of our knowledge, the present study is the first to identify a group of circRNAs expressed in NFPA. In the present study, Cox's regression and RSFVH algorithm were used to select 9 circRNAs that were most closely associated with progression or relapse in the set of 7,481 circRNAs. Then, a risk score survival prediction method that exhibited good predictive characteristics was used to combine the circRNAs. Finally, a risk score for the combination of hsa_circ_0000066 and hsa_circ_0069707 was obtained, which had the largest AUROC curve with the highest predictive power. The two-circRNA-based classifier was able to separate patients with NFPA into low-risk or high-risk early recurrence groups. Disease recurrence may be predicted more accurately by directly constructing an early recurrence model rather than using traditional categorical indicators. According to the present analysis, hsa_circ_0000066b is able to bind to hsa_circ_0069707, and they work together to affect the PFS time of patients by modulating the response of transport vesicles and cells to unfolded proteins. Thus, hsa_circ_0000066b and hsa_circ_0069707 serve an important role in NFPA recurrence. Bioinformatics may analyse and infer only the functions of these circRNAs; thus, it remains necessary to confirm the biological effects of these two circRNAs in tumourigenesis in experimental studies.

In addition to the limited availability of NFPA sequencing data, the present study has a number of limitations that need to be considered. First, the prognostic circRNAs confirmed here are likely not the only circRNA candidates associated with NFPA PFS as only a fraction of human circRNAs (88,750 out of 140,000+) were included in the present analysis. Therefore, the present results should be further validated through prospective and multi-centre studies. Secondly, the present study lacks information on the mechanisms by which the two circRNAs affect the prognosis of patients with NFPA. Further functional experimental studies of primary cells or the 293 cell line should be performed to determine whether these circRNA directly affect NFPA progression. Lastly, although the present results were outlined in the test dataset set as much as possible based on data availability, this marker has not yet been prospectively tested in clinical trials. However, despite these deficiencies, the association between this circRNA signature and PFS in the present dataset suggests that it is a potent prognostic marker for NFPA.

In conclusion, to the best of our knowledge, the present study is the first circRNA signature identified that predicts tumour recurrence in patients with NFPA with a high prediction accuracy and thus may be used for treatment guidance and prognosis estimation.