The Prognostic Value of the Circulating Tumor Cell-Based Four mRNA Scoring System: A New Non-Invasive Setting for the Management of Bladder Cancer

Simple Summary Bladder cancer with similar diagnosis based on traditional classification exhibits different behaviors and therapeutic outcomes. Thus, circulating tumor cells (CTCs) represent a more accurate approach to investigate bladder cancer features. Our results demonstrate that risk score based on EGFR, TRPM4, TWIST1, and ZEB1 four-gene signature in CTCs is markedly and undoubtedly associated with recurrence, suggesting an innovative and non-invasive strategy to manage both non muscle invasive and muscle invasive bladder cancer progression without the necessity of repetitive and onerous cystoscopies. Abstract Bladder cancer (BC) is one of the most expensive lifetime cancers to treat because of the high recurrence rate, repeated surgeries, and long-term cystoscopy monitoring and treatment. The lack of an accurate classification system predicting the risk of recurrence or progression leads to the search for new biomarkers and strategies. Our pilot study aimed to identify a prognostic gene signature in circulating tumor cells (CTCs) isolated by ScreenCell devices from muscle invasive and non-muscle invasive BC patients. Through the PubMed database and Cancer Genome Atlas dataset, a panel of 15 genes modulated in BC with respect to normal tissues was selected. Their expression was evaluated in CTCs and thanks to the univariate and multivariate Cox regression analysis, EGFR, TRPM4, TWIST1, and ZEB1 were recognized as prognostic biomarkers. Thereafter, by using the risk score model, we demonstrated that this 4-gene signature significantly grouped patients into high- and low-risk in terms of recurrence free survival (HR = 2.704, 95% CI = 1.010–7.313, Log-rank p < 0.050). Overall, we identified a new prognostic signature that directly impacted the prediction of recurrence, improving the choice of the best treatment for BC patients.


Introduction
Bladder cancer (BC) is the tenth most common cancer worldwide and its incidence is gradually increasing, particularly in developed countries [1]. The American Cancer Society estimates that for 2022, there will be about 81,180 new cases and 17,100 deaths for BC in biomarkers' and 'overall survival' resulted in 141 hits and in 86 hits for 'genitourinary cancer' with 'transient receptor potential channel'. All research articles were written in English, dating back to the maximum year of publication of 2010, and the evaluated biomarkers involved in bladder cancer growth and progression were included. In vitro and in vivo studies concerning human samples were also included. Studies with limited information or with an inaccurate number of samples were excluded. The validity of the articles, after applying the inclusion and exclusion criteria, were assessed by three researchers (C.A., M.B.M., and F.M.).
The Cancer Genome Atlas (TCGA) and Genotype-Tissue Expression (GTEx) datasets were screened and the accessible transcriptomic data of patients diagnosed with BC were obtained (http://www.gepia.2cancer-pku.cn, accessed on 22 October 2018) with OS available for 404 patients.

Study Patients
This translational study included 60 patients with NMIBC or MIBC, hospitalized before transurethral resection of bladder (TURB), at the Urologic and Andrologic Clinics, University of Perugia, Perugia, Italy. Moreover, to improve the clarity/soundness of the results, the exclusion criteria were: another histological type of cancer, had received neoadjuvant or adjuvant chemotherapy, had a medical history of the same or another type of cancer, and had received a T0 stage classification. Enrolled patients had a clinical and radiological evidence of BC and did not undergo systemic/locally therapy before and after the TURB and the blood sample collection. The diagnosis was confirmed by histological analysis. The patient data regarding the age, tumor stage, histological type, regional lymph node involvement, and number of metastatic sites were also tabulated and statistically analyzed (Table 1). Local was used to define the locoregional recurrence at the surgical field, and distant to describe the presence of metastatic disease outside the urinary tract and out of the locoregional surgical field. Imaging techniques such as computed tomography were applied to assess the presence of metastasis. No patient died during the investigation period nor underwent a cystectomy. Table 1. The baseline clinico-pathologic characteristics in the prospective study of patients with BC.

Variable Patients
All, n (%) 60  Blood samples from 53 BC patients were collected and analyzed between March 2018 and March 2020. Healthy donors (n = 5) were age-matched men without BC who were recruited for normal tissue specimens according to the Ethics Committee of the University of Perugia-approved protocol. Moreover, fresh samples (n = 8) of primary tumors were recovered for comparative analyses.
The present study was approved by the Ethical Committee CEAS UMBRIA (Ethic Committee approval code: URO009-3171/18). The use of patient data and CTCs for research purposes at the University of Perugia were executed pursuant to the Italian legislation and international standards. Anonymous written informed consent was obtained from each participant or their guardians prior to the study enrolment.

Patient Samples
For CTC detection, 7 mL of peripheral blood from each patient was collected into Vacutainer ® Blood Collection Tubes EDTA (BD Biosciences, Franklin Lake, NJ, USA) on the day of surgery. Blood was obtained at the middle of vein puncture after the first 2 mL of blood was discarded to avoid contamination by normal epithelial cells. Peripheral blood was immediately processed within 3 h by means of ScreenCell devices (Sarcelles, France), as previously described [38], to isolate the CTCs using a size-based method already used in rare tumors [39] as well as in more common tumors including BC [40].
Seven blood samples from healthy donors were obtained from the blood bank to establish the normal expression of each marker in the white blood cells as previously described [38].
In addition, fresh BC tissue samples (biopsies), obtained during TURB from eight patients, were used in the study.

RNA Extraction
A Single Shot Cell Lysis Kit (BioRad, Hercules, CA, USA) was employed to extract the total RNA according to the protocol. Messenger RNA from biopsies was extracted by the RNeasy Mini Kit (Qiagen, Milan, Italy). Messenger RNAs from five human normal bladder samples (NHB) were purchased from AMS Biotechnology (Abigdon, UK).

Reverse Transcription
The iScript Advanced cDNA Synthesis Kit (BioRad) was employed to retrotranscribe RNA and the SSOADvanced PreAmp Kit and PrimePCR PreAMP Assays (Bio-Rad) were used to preamplify the cDNA including all of the primers used in the gene expression analysis.

InCTC Assay
The InCTC assay was performed starting from blood diluted with PBS plus 2% fetal bovine serum (FBS), processed with Ficoll-Paque premium gradient medium (GE Healthcare, Milan, Italy). The monolayer containing the CTCs and leukocytes was then collected and washed with RPMI medium without serum. A micro-Boyden transwell chamber (Greiner Bio One, Cassina de Pecchi, Italy) with 8-µm-diameter pore inserts coated with Matrigel was exploited to assess the migration ability of the CTCs. The CTCs were seeded onto the upper chamber in the medium without serum, while 10% FBS was added to the lower chamber. The chambers were then incubated in a humidified incubator with 5% CO 2 . After 24 h, the inserts were washed with PBS, fixed with ice-cold methanol, rinsed with PBS, and swabbed with a cotton swab to remove the non-migrated cells. Toluidine blue staining was performed to visualize the CTCs with light microscopy analysis, with DAPI and antipan CK for confocal analysis. The T24 cell line was used as a positive control. The T24 cell line, cultured in complete RPMI 1640 medium (10% FBS plus 1% penicillin-streptomycin) and grown at 37 • C in a 5% CO 2 incubator, was used as the positive control.

Protein-Protein Interaction (PPI) Network Analysis
The search tool for retrieval of interacting genes (STRING) (https://string-db.org, accessed on 22 October 2018) database, based on known and predicted PPIs, was employed to seek potential interactions between the markers [41]. Text mining, experiments, databases, co-expression, species limited to "Homo sapiens", and an interaction score >0.7 were considered as active interaction sources and applied to construct the PPI networks. The PPI network was visualized by Cytoscape software version 3.6.1. (https://cytoscape.org/, accessed on 22 October 2018) In the networks, proteins are schematized as nodes and interactions as edges.

Prognostic Signature
The univariate Cox proportional hazard regression model was used to assess the association of gene expression in CTC with the recurrence free survival (RFS) of patients included in the cohort. Genes with p ≤ 0.05 were considered statistically significant. Using the regression coefficients from the Cox regression analysis of significant genes, we calculated a risk score for each patient based on their individual expression levels. In order to improve reliability, patients were stratified using the median risk score as a cutoff value, and the BC patients were classified into high-(1) and low-risk (0) groups, and Kaplan-Meier and receiver operating characteristic (ROC) analyses were performed.

Statistical Analysis
The quantitative variables are expressed as the mean values ± standard deviation (SD) or as the median values (interquartile range). Qualitative variables are expressed as absolute and relative frequencies. For the comparisons of the proportions, categorical data were tested by the Fisher's exact test.
The differential expression gene (DEG) analysis was performed by applying the Mann-Whitney U t-test. The correlation matrix was analyzed by the Spearman correlation test.
The Cox proportional hazard regression model was used to perform the univariate and multivariate survival analyses. Only statistically significant variables in the univariate analysis were incorporated into the risk score equation. To assess whether the risk score was independent of other clinical characteristics, multivariable Cox regression analysis and stratification analysis were used.
The differences in the recurrence time of BC patients stratified by low-and high-risk according to gene expression in CTCs was assessed by the Kaplan-Meier method. The RFS was considered as the time from the TURB surgery to the date of the first recurrence or, if the recurrence of tumor was not present, the date of the last follow-up. The statistical significance of the observed differences between groups was determined by the log-rank test.
The grading variable was dichotomized by the time-dependent ROC curve analysis. The cutoff point estimated by the ROC analysis indicates the distinguishing characteristic of the grading variable used to classify participants into the recurrence group. The accuracy of the dichotomous results was calculated through the area under the ROC curve (AUC).
All testing was carried out using GraphPad Prism 9.1 software (GraphPad Software, La Jolla, CA, USA) or SPSS (version 26.0; Chicago, IL, USA) or R (version 4.1.2, Auckland, NZ, USA). All p-values were two-sided, and p ≤ 0.05 was considered to indicate statistical significant.

Selection of Previously Described BC Biomarker Candidates
Analysis of the PubMed search highlighted several bladder-specific transcripts implicated in BC growth and metastasis. Among them, we selected BIRC5, CDH11, EGFR, EPCAM, KRT18, SPP1, TRPA1, TRPC1, TRPC3, TRPC6, TRPV6, TRPM4, TWIST1, TNC, UPK2, VEGFA, VIM, and ZEB1 as potential biomarker candidates in our investigation to develop a prognostic signature model. Figure 1 illustrates the work flow diagram of this study. All testing was carried out using GraphPad Prism 9.1 software (GraphPad Software, La Jolla, CA, USA) or SPSS (version 26.0; Chicago, IL, USA) or R (version 4.1.2, Auckland, NZ, USA). All p-values were two-sided, and p ≤ 0.05 was considered to indicate statistical significant.

Selection of Previously Described BC Biomarker Candidates
Analysis of the PubMed search highlighted several bladder-specific transcripts implicated in BC growth and metastasis. Among them, we selected BIRC5, CDH11, EGFR, EPCAM, KRT18, SPP1, TRPA1, TRPC1, TRPC3, TRPC6, TRPV6, TRPM4, TWIST1, TNC, UPK2, VEGFA, VIM, and ZEB1 as potential biomarker candidates in our investigation to develop a prognostic signature model. Figure 1 illustrates the work flow diagram of this study.  First, through the gene expression profiling interactive analysis (GEPIA2), we analyzed the different expression of these candidates between the normal and tumor bladder tissue ( Figure 2A). Results showed that the BIRC5, EPCAM, KRT18 and SPP1 genes were upregulated, and CDH11, TRPA1, TRPC1, TRPC3, TWIST1, VIM, and ZEB were downregulated whereas EGFR, TNC, TRPC6, TRPV6, TRPM4, UPK2, and VEGFA were not differently expressed in BC with respect to the normal tissue. Then, we constructed the network connecting the selected key genes to the predicted functional associations ( Figure 2B). The PPI network was visualized in the form of a graph network. Most of the proteins were functionally associated, while UPK2 was completely not connected to the network and so it was not included in the successive analyses. Interestingly, the cluster of TRP cation channels was associated with the rest of the group by the gene PLCG1, which codifies for a calcium-dependent enzyme involved in signal transduction and in the regulation of TRP activity [42][43][44]. First, through the gene expression profiling interactive analysis (GEPIA2), we analyzed the different expression of these candidates between the normal and tumor bladder tissue ( Figure 2A). Results showed that the BIRC5, EPCAM, KRT18 and SPP1 genes were  upregulated, and CDH11, TRPA1, TRPC1, TRPC3, TWIST1, VIM, and ZEB were downregulated whereas EGFR, TNC, TRPC6, TRPV6, TRPM4, UPK2, and VEGFA were not differently expressed in BC with respect to the normal tissue. Then, we constructed the network connecting the selected key genes to the predicted functional associations ( Figure  2B). The PPI network was visualized in the form of a graph network. Most of the proteins were functionally associated, while UPK2 was completely not connected to the network and so it was not included in the successive analyses. Interestingly, the cluster of TRP cation channels was associated with the rest of the group by the gene PLCG1, which codifies for a calcium-dependent enzyme involved in signal transduction and in the regulation of TRP activity [42][43][44].

Isolation and Gene Expression Profiling of CTCs from BC Patients
It is well-known that CTCs acquire different features with respect to the primary mor cells by changing gene expression [36]. Thus, we evaluated the gene expression le of the biomarker candidates in CTCs isolated from BC patients. The analyzed cohort cluded 60 newly diagnosed patients with BC, of which the CTCs (53 patients) or biop (8 patients) were evaluated. For one patient, it was possible to collect both the CTCs biopsy. The median age of the enrolled patients was 74 years. The demographic and m clinic-pathological features of the patients enrolled in the study are summarized in Ta 1.
As shown, the BC clinical stage at the time of diagnosis ranged from the Tx to stage. However, the majority of the patients had Ta stage disease. None of the patie had a distal metastasis, nor were positive nodes found at diagnosis. Recurrence was orded in 41.4% patients. The mean follow-up period was 17.6 months, with a media 16.8 months.

Isolation and Gene Expression Profiling of CTCs from BC Patients
It is well-known that CTCs acquire different features with respect to the primary tumor cells by changing gene expression [36]. Thus, we evaluated the gene expression levels of the biomarker candidates in CTCs isolated from BC patients. The analyzed cohort included 60 newly diagnosed patients with BC, of which the CTCs (53 patients) or biopsies (8 patients) were evaluated. For one patient, it was possible to collect both the CTCs and biopsy. The median age of the enrolled patients was 74 years. The demographic and major clinic-pathological features of the patients enrolled in the study are summarized in Table 1.
As shown, the BC clinical stage at the time of diagnosis ranged from the Tx to T4 stage. However, the majority of the patients had Ta stage disease. None of the patients had a distal metastasis, nor were positive nodes found at diagnosis. Recurrence was recorded in 41.4% patients. The mean follow-up period was 17.6 months, with a median of 16.8 months.
The blood samples were processed by ScreenCell devices and the CTCs were detected in all patients with the exception of three. The isolated CTCs were characterized for the expression of pan-CK and for the absence of CD45, as previously described [45] ( Figure 3A). In addition, the cell invasion assay was used to investigate the viability and functionality of the isolated CTCs. For this purpose, CTCs enriched from 30 patients, of which 23 were NMIBC and seven were MIBC, were evaluated by using the in vitro invasion assay. Our results showed that the CTCs from 20 BC patients (66.7%) spread into the Matrigel, indicating an invasive and functional phenotype ( Figure 3B,C). T24 BC cells, known to be high-grade and invasive, were used as the positive control [46]. Interestingly, we found that the group of invading cells was composed of CTCs derived from both the NMIBC (16 patients) and MIBC patients (four patients), suggesting that the dissemination in the blood of cells, potentially able to invade tissue, occurs very early during BC progression ( Figure 3D). ncers 2022, 14, x FOR PEER REVIEW 10 functionality of the isolated CTCs. For this purpose, CTCs enriched from 30 patient which 23 were NMIBC and seven were MIBC, were evaluated by using the in vitro i sion assay. Our results showed that the CTCs from 20 BC patients (66.7%) spread into Matrigel, indicating an invasive and functional phenotype ( Figure 3B,C). T24 BC c known to be high-grade and invasive, were used as the positive control [46]. Interestin we found that the group of invading cells was composed of CTCs derived from both NMIBC (16 patients) and MIBC patients (four patients), suggesting that the dissemina in the blood of cells, potentially able to invade tissue, occurs very early during BC gression ( Figure 3D). Then, through the ddPCR, we analyzed the expression of the selected markers in lated CTCs . BIRC5, CDH11, EGFR, EPCAM, KRT18, SPP1, TRPA1, TRPC1, TRP  TRPC6, TRPV6, TRPM4, TWIST1, TNC, VEGFA, VIM, and ZEB1 genes were heterog ously expressed. Furthermore, we compared the gene expression levels obtained in CTCs, evaluated as fold changes with respect to the normal human bladder RNAs (NH with biopsies. A significantly increased expression of BIRC5, CDH11, EGFR, EPC KRT18, SPP1, TNC, TRPC1, TRPC6, TRPM4, TRPV6, TWIST1, VEGFA, and ZEB1 found in the CTCs with respect to NHB and, even more interestingly, compared with Then, through the ddPCR, we analyzed the expression of the selected markers in isolated CTCs . BIRC5, CDH11, EGFR, EPCAM, KRT18, SPP1, TRPA1, TRPC1, TRPC3,  TRPC6, TRPV6, TRPM4, TWIST1, TNC, VEGFA, VIM, and ZEB1 genes were heterogeneously expressed. Furthermore, we compared the gene expression levels obtained in the CTCs, evaluated as fold changes with respect to the normal human bladder RNAs (NHB), with biopsies. A significantly increased expression of BIRC5, CDH11, EGFR, EPCAM, KRT18, SPP1, TNC, TRPC1, TRPC6, TRPM4, TRPV6, TWIST1, VEGFA, and ZEB1 was found in the CTCs with respect to NHB and, even more interestingly, compared with the biopsies (Figures 4 and S1).
Cancers 2022, 14, x FOR PEER REVIEW 11 of 21 Figure 4. The expression of biomarker candidates in the CTCs isolated from BC patients. The expression was evaluated by ddPCR compared with the NHB and BC biopsies. Genes found to be expressed at significant higher/lower levels in CTCs with respect to the BC biopsies or NHB, are shown. The gene expression levels, evaluated by ddPCR, are expressed as fold changes with respect to NHB, used as calibrators (=1). NS = not statistically significant, p ≤ 0.05 was considered as statistically significant.
Neither TRPA1 nor TRPC3 expression levels differed statistically from NHB or the biopsies (data not shown), so we excluded them from further analysis. No cells with CTClike features were found in the peripheral blood of healthy donors, whose WBCs were similarly analyzed, demonstrating the negligible expression of these genes (data not shown).
In addition, the positive detection of BIRC5, CDH11, EGFR, SPP1, TNC, TRPC6, TRPM4, TRPV6, TWIST1, VEGFA, and ZEB1 was more present in the CTC group compared with the biopsy group (Table 2). Neither TRPA1 nor TRPC3 expression levels differed statistically from NHB or the biopsies (data not shown), so we excluded them from further analysis. No cells with CTC-like features were found in the peripheral blood of healthy donors, whose WBCs were similarly analyzed, demonstrating the negligible expression of these genes (data not shown).

Markers
Biopsy n (%) CTC n (%) p-Value Overall, these data indicate that the increase in the expression levels of several genes involved in the processes of angiogenesis, epithelial-mesenchymal transition (EMT), and cation flux-regulation pathways is associated with the acquisition of the CTC phenotype.

Biomarkers Predictive of Recurrence
Then, to investigate the relationship between the marker gene expression levels in CTCs and tumor stage, we analyzed the mRNA expression levels of BIRC5, CDH11, EGFR, EPCAM, KRT18, SPP1, TRPC1, TRPC6, TRPV6, TRPM4, TWIST1, TNC, VEGFA, VIM, and ZEB1 in the CTCs from the NMIBC and MIBC patients ( Figure S2). Our data showed no statistically significant differences by comparing the NMIBC with MIBC patients. In addition, to investigate the association between the molecular profile and histological classification, we analyzed the expression of the selected markers in the CTCs by stratifying the BC patients according to high (grade HG) and low grade (LG). Our results demonstrated that BIRC5, EGFR, SPP1, TNC, TRPM4, VEGFA, TWIST1, and ZEB1 were upregulated whereas KRT18, TRPV6, and VIM were downregulated in the HG with respect to the LG BC ( Figure 5A). No differences were found in the expression levels of CDH11, EPCAM, TRPC1, and TRPC6, thus they were excluded in the following analysis. This result prompted us to evaluate the correlation index in order to strengthen the relationship between the biomarkers according to their expression levels in the CTCs. The data established the presence of two well-defined clusters of strong association: one group consisting of BIRC5, EGFR, KRT18, SPP1, TNC, TWIST1, TRPM4, VEGFA, and ZEB1 genes and the second of TRPV6 and VIM, suggesting a suitable gene model for CTCs from BC patients ( Figures 5B,C and S3).

Identification of a Four-mRNA Based Signature
Through the univariate Cox analysis, we analyzed whether the biomarker expression levels of the two clusters were associated with the prognosis of BC patients. Our results showed that the RFS was significantly associated with the expression levels of EGFR, TRPM4, TWIST1, and ZEB1 belonging to the same cluster. These data outline a four-mRNA-based model to predict the RFS time in BC patients (Table 3). Furthermore, Kaplan-Meyer analysis was used to validate the association between the expression of four gene candidates and the prognosis of patients. For this purpose, first, we stratified patients into high-and low-expressing, according to the ROC analysis, for the EGFR, TRPM4, TWIST1, and ZEB1 markers. As shown in Figure 6A, patients expressing high levels of EGFR, TWIST1, TRPM4, and ZEB1 had worse RFS compared with those that had low levels. Profiles retrieved from the GTEx showed the insignificant expression of the four selected biomarkers in whole blood with respect to the bladder, making our analysis consistent and reliable ( Figure S4).

The Risk Score Based on Four-mRNA Signature Predicted the RFS of BC Patients: Validation of the Prognostic Model
To facilitate the application of identified RFS-related mRNAs in the clinical management of BC patients, the risk score for each patient was constructed with the regression coefficient from the univariate Cox analysis by using: Risk score = (3.002 × Expression EGFR + 4.514 × Expression TRPM4 + 5.360 × Expression TWIST1 + 3.590 × Expression ZEB1). Then, the patients were stratified into high-and low-risk groups according to the median risk score, which was used as the cutoff point. The Kaplan-Meier analysis showed a shorter RFS in the high-risk group with respect to the low-risk group ( Figure 6B) as also supported by the ROC curve, reaching an AUC value of 0.6863 ( Figure 6C). Overall, these results demonstrated that the risk score displayed a substantially effective performance for RFS prediction and it can be used as a prognostic indicator. Confocal microscopy was also used to assess the expression at the protein level of the four biomarkers. The data showed that the CTCs expressed at the protein levels the four markers EGFR, TRPM4, TWIST1, and ZEB1 ( Figure 6D). Finally, multivariate analysis, by considering the risk score, grade, stage, diabetes type II, hypertension, and RFS, showed that the four-mRNA-based signature was an independent prognostic factor in BC patients (Table 4).

The Risk Score Based on Four-mRNA Signature Predicted the RFS of BC Patients: Validation of the Prognostic Model
To facilitate the application of identified RFS-related mRNAs in the clinical management of BC patients, the risk score for each patient was constructed with the regression coefficient from the univariate Cox analysis by using: Risk score = (3.002 × Expression EGFR + 4.514 × Expression TRPM4 + 5.360 × Expression TWIST1 + 3.590 × Expression ZEB1). Then, the patients were stratified into high-and low-risk groups according to the median risk score, which was used as the cutoff point. The Kaplan-Meier analysis showed a shorter RFS in the high-risk group with respect to the low-risk group ( Figure 6B) as also supported by the ROC curve, reaching an AUC value of 0.6863 ( Figure 6C). Overall, these results demonstrated that the risk score displayed a substantially effective performance for RFS prediction and it can be used as a prognostic indicator. Confocal microscopy was also used to assess the expression at the protein level of the four biomarkers. The data showed that the CTCs expressed at the protein levels the four markers EGFR, TRPM4, TWIST1, and ZEB1 ( Figure 6D). Finally, multivariate analysis, by considering the risk score, grade, stage, diabetes type II, hypertension, and RFS, showed that the four-mRNAbased signature was an independent prognostic factor in BC patients (Table 4).

Discussion
The lack of an accurate classification system predicting the risk of recurrence or the progression of BC patients hinders the search of new biomarkers. Indeed, even though TURBT is the standard of care for the management of BC [47], its main purpose is quite different in NMIBC and MIBC, and NMIBC patients still have a heterogeneous prognosis [48]. The fact that the prognosis for BC patients has not changed for over three decades [49], and the identified molecular markers are rarely discriminatory enough for clinical use [50,51] suggests the need of a deeply understanding if this pathology.
In the last few years, most of the studies have focused their attention on liquid biopsy and the enrichment of CTCs, also in BC [18,52,53]. In particular, the number of CTCs detected in the peripheral blood of BC patients has been considered as a prognostic indicator of cancer prognosis [25,54]. The presence of CTCs predicting decreased time to first recurrence and time to progression identified a subgroup of patients with super-high-risk CTC-positive NMIBC. There was a non-statistical difference in the number of preoperative CTCs between the NMIBC and MIBC patients, suggesting that further investigation is needed [26].
The spread of CTCs in the bloodstream represents a fundamental moment for the spread of cancer, even if it is now accepted that only a small percentage of these cells can survive in the circulation and therefore give rise to metastases [55]. In line with these previous findings, we demonstrated that CTCs, obtained by size-based isolation, are present in the peripheral blood of BC patients. Interestingly, we also showed that the isolated CTCs from the NMIBC and MIBC patients displayed the same in vitro invasive ability. This indicates that their dissemination occurs very early during BC progression and that CTC features are independent from the phenotype of the primary tumor counterpart, underlying the importance to better explore CTCs as a source of potential prognostic biomarkers that are useful immediately after the diagnosis.
Thus, our aim was to identify a specific CTC gene expression profile correlated to the ability of recurrence. First, we analyzed several potential markers such as BIRC5, CDH11,  EPCAM, KRT18, SPP1, TIMP2, UPK2, TRPA1, TRPC1, TRPC3, TRPC6, TRPV6, TRPM4,  TWIST1, TNC, VEGFA, VIM, and ZEB1 disclosed as related to the BC primary tumor in studies published back to 2010. The expression levels of these biomarker candidates were then analyzed in the CTCs from BC patients. No differences were found in the expression by comparing the CTCs from the NMIBC with MIBC patients whereas several markers were expressed at different levels in HG with respect to the LG samples. For this reason, we investigated them in order to evaluate the relation with RFS. The univariate analysis, also supported by the Kaplan-Meyer, showed that in our study including the NMIBC and MIBC patients, four biomarkers-EGFR, TRPM4, TWIST1, and ZEB1-are strongly associated with RFS.
Our data were in line with recent findings that describe the association between high EGFR expression and shorter RFS, not only in the MIBC, but also in stage T1 NMIBC patients, suggesting early stratification is essential in order to select the most appropriate targeted treatment [56]. Encouraging results have been achieved by using small molecules targeting EGFR such as erlotinib [57]. In fact, EGFR signaling regulates numerous cellular pathways associated with tumor invasion and metastasis in preclinical models of BC. Moreover, recently, the targeted toxin EGF-PE40, obtained by combining the ligand EGF with PE40, a truncated version of Pseudomonas Exotoxin A, was generated and the expression of EGFR in BC cells was used to internalize the toxin and promote cell death [58].
The TRPM4 channel is a non-selective monovalent cationic channel expressed in several tissues including smooth muscle and urothelial cells of bladder. The TRPM4 activation leads to sodium ion influx into the cell with consequent plasma membrane depolarization and calcium signaling [59]. It has been suggested that TRPM4 has a key role in normal bladder contraction. However, the absence of information about its involvement in urothelial cancer suggests that further efforts are needed to understand its physiological and pathological functions. On the other hand, it is well-known that TRPM4 is involved in the growth and progression of several cancers. In fact, germline variants of TRPM4, affecting the intestinal mucosal integrity, have been found in colon rectal cancer cells [60]. Moreover, TRPM4 is responsible for stemness mediation in breast cancer. It is overexpressed in breast cancer tumor-spheres and its inhibition displays anti-tumor effects by directly targeting cancer stem cells [61]. Finally, it has also been classified as a cancer driver gene in androgen-independent prostate cancer where it regulates the EMT, which is fundamental for cancer cell migration and invasion. Therefore, given the numerous associations with cancer hallmark functions, TRPM4 has been suggested to be a new potential anticancer drug target [62].
TWIST1 overexpression is frequently detected in several cancers also including BC [63]. Its expression has been associated with the invasive aggressive cancer behavior, recurrence, therapy failure, and poor prognosis for the patient. In fact, it is well-accepted that it is involved in tumor initiation, intravascular migration, EMT, metastasis, and angiogenesis [64]. In addition, it has been demonstrated in BC that TWIST1, by regulating the expression of ABC transporters and the AKT signaling, induces resistance to several chemotherapeutic drugs [63]. Finally, the expression of TWIST1 was significantly enhanced in the metastatic lesions with respect to the BC primary site, consolidating its pivotal role in the metastatic process [65]. Thus, our data, obtained by using CTCs from NMIBC and MIBC, show the inclusion of TWIST1 in the four-based gene signatures, are in agreement with these previous findings, which support the potential diagnostic and prognostic role of TWIST1 in BC [45,65].
Finally, ZEB1, the fourth biomarkers identified in our study, was significantly overexpressed in BC in comparison with the healthy adjacent tissues. Moreover, it is associated with vasculogenic mimicry and with EMT promotion in BC [66]. In fact, its knockdown markedly reduced the formation of vasculogenic mimicry in the BC cell lines. Recently, it has also been demonstrated that the VIM Antisense RNA 1/miR-655/ZEB1 axis regulates metastasis in BC by modulating the EMT. In fact, VIM Antisense RNA 1, by competing with ZEB1 for miR-655 binding, abolishes the miR-655-mediated suppression of ZEB1 that is then able to finally affect EMT [67]. Our results are supported by findings described in breast cancer cells, where ZEB1 expression plays a pivotal role in the maintenance of stem-like features, immune evasion, epigenetic reprogramming, and aberrant cellular polarity [68]. Overall, these findings strongly highlight the implication of ZEB1 in the resistance and survival of disseminated cancer cells.
Identifying the features of primary tumors is insufficient to manage the disease's progression, as single-site tumor biopsies cannot accurately portray the extremely variable cell population in the tumor mass. In this regard, CTCs are found to be genetically different and subjected to a distinct evolutionary process than the primary tumor mass. For this reason, the non-invasive liquid biopsy and the examination of CTCs could represent a good choice to replace the tissue biopsy in order to better evaluate the optimal therapy and predict the prognosis of the patients [69].

Conclusions
The results from this study demonstrated that the expression of the four-gene signature in CTCs was markedly associated with recurrence, suggesting an innovative strategy to manage disease progression without repetitive invasive procedures.
In the near future, the precise stratification of patients based on the EGFR, TRPM4, TWIST1, and ZEB1 gene signature risk score obtained in CTCs could represent, immediately after the diagnosis, a strategy to better plan the clinical evolution of the patient and therefore select the best therapeutic management. Indeed, faster follow-up, early radical surgery, or systemic therapy may be approved in the very high-risk group of patients, leading to the possibility of reducing the number of cystoscopies performed that are onerous for patients and expensive for health care providers. However, we believe that it is urgent to conduct studies with larger cohorts of patients to validate this risk score in the clinical practice.
Supplementary Materials: The following supporting information can be downloaded at: https: //www.mdpi.com/article/10.3390/cancers14133118/s1. Figure S1: The gene expression of selected biomarker candidates in biopsies with respect to the NHB. Figure S2: The gene expression analysis of selected biomarkers; Figure S3: The correlation matrix analysis; Figure S4: The profiles retrieved from the GTEx to support and validate our selection. Informed Consent Statement: Informed consent was obtained from all subjects involved in the study. Written informed consent was obtained from the patients to publish this paper.

Data Availability Statement:
The data that support the findings of this study are available from the corresponding authors upon request.