Circulating Tumor DNA Analyses Predict Disease Recurrence in Non-Muscle-Invasive Bladder Cancer

Circulating tumor DNA (ctDNA) can be a prognostic biomarker for non-muscle-invasive bladder cancer (NMIBC); however, targeted sequencing has not been performed to detect ctDNA in NMIBC. We applied targeted sequencing based on an 861-gene panel to determine mutations in tumor tissue DNA and plasma ctDNA in 82 NMIBC patients receiving transurethral resection (TUR) of bladder followed by immunotherapy. We detected 476 and 165 somatic variants in tumor DNA from 82 NMIBC patients (100%) and ctDNA from 54 patients (65.85%), respectively. Patients with high heterogeneity in tumor DNA had a significantly shorter disease-free survival than those with low heterogeneity. Tumor-derived alterations were detectable in plasma of 43 patients (52.44%). The concordance of somatic variants between tumor DNA and plasma ctDNA were higher in patients with T1 stage (p < 0.0001) and tumor size ≥3 cm (p = 0.0002). Molecular tumor burden index (mTBI) in ctDNA positively correlated with larger tumor size (p = 0.0020). A higher mTBI was an independent predictor of recurrence after TUR of bladder followed by immunotherapy. Analysis of ctDNA based on targeted sequencing is a promising approach to predict disease recurrence for NMIBC patients receiving TUR of bladder followed by immunotherapy.


INTRODUCTION
Bladder cancer is a common malignant tumor, with an estimated 549,000 new cases and 200,000 deaths worldwide every year (1). The most common histological type of bladder cancer is urothelial carcinoma (91.4%), which can be divided into non-muscle-invasive bladder cancer (NMIBC) and muscle-invasive bladder cancer (MIBC), accounting for 74 and 25.2%, respectively (2). NMIBC is characterized by a high recurrence rate of 50-70% and a 5-year survival rate of 90% (3). Progression to MIBC is observed in 20-30% of patients with NMIBC (4). Ultimately, as many as 10-15% of NMIBC patients will die of bladder cancer (5). Thus, it is urgent to develop a biomarker for discriminating high-risk from low-risk NMIBC patients, and for prediction of disease recurrence.
In the past 5 years, considerable advancements have been made with liquid biopsy, which is an important method in precision medicine, especially in oncology. The advantage of liquid biopsy is that it only requires the patient's blood for noninvasive detection without invasive harvesting of the patient's tissue and can monitor gene mutations in patients in real time (6,7). Circulating tumor DNA (ctDNA) is a DNA fragment secreted into peripheral blood by apoptotic tumor cells (8). ctDNA only constitutes a small fraction of circulating cell-free DNA (cfDNA) in peripheral blood; the majority of cfDNA stems from leukocytes. Therefore, it is difficult to detect somatic variants in ctDNA, especially for early-stage cancers.
Currently, plasma ctDNA has exhibited great potential as a biomarker in several cancers, particularly in advanced metastatic cancer. For example, a high concordance rate of the mutation status and changes in copy number between liquid and solid biopsy have been reported for metastatic castration-resistant prostate cancer (9). The European Medicines Agency approved liquid biopsy for epidermal growth factor receptor for lung cancer patients (10,11). However, there are only two reports on ctDNA in NMIBC. Birkenkamp-Demtröder et al. established methods for non-invasive disease surveillance in patients with NMIBC, using highly sensitive and tumor-specific personalized plasma-and urine-based assays, and demonstrated that somatic variants of ctDNA were detectable in plasma and urine samples from patients with NMIBC (12). Christensen et al. developed droplet digital polymerase chain reaction (ddPCR) analyses and screened ctDNA for FGFR3 and PIK3CA hotspot mutations in urine and plasma from NMIBC and MIBC patients undergoing radical cystectomy (13). They demonstrate that increased levels of FGFR3 and PIK3CA mutated DNA in urine and plasma are indicative of later progression and metastasis in bladder cancer. Nevertheless, these two studies both utilized ddPCR assays, and whether next-generation sequencing (NGS) approaches can be applied in screening ctDNA for somatic variations in plasma from NMIBC patients has not been explored.
In this study, we performed targeted sequencing based on an 861-gene panel to screen matched plasma ctDNA and tumor DNA for somatic mutations in patients with NMIBC, evaluated the effects of clinicopathological features on the concordance of somatic variants between ctDNA and tumor DNA, and analyzed whether ctDNA could be a biomarker for prediction of disease recurrence after transurethral resection (TUR) of bladder followed by immunotherapy.

Patients and Sample Collection
This study included 82 NMIBC patients who underwent TUR of bladder followed by bladder transfusion immunotherapy [agent: Bacille Calmette-Gueŕin (BCG)]; mode of administration: bladder transfusion immunotherapy; number of doses/cycles: the period of perfusion is usually once a week, a total of six to eight times, and then once a month, a total of 10 times) to avoid recurrence between May 2015 and June 2016 at Lanzhou University Second Hospital. Archival slides of patients were evaluated by two pathologists. We obtained the following clinicopathological information: gender, age upon diagnosis, t u m o r g r a d e , t u m o r s i z e , a n d t u m o r s t a g e . T h e clinicopathological information was supplied in Table 1. Approval for this study was obtained from Medical Ethics Committee of Shenzhen Second People's Hospital. Patients signed an informed consent. The study was conducted in compliance with the principles of the Declaration of Helsinki and local ethical and legal requirements.

DNA Extraction
Tumor tissue samples were obtained by TUR of bladder, and at least 10 ml of Peripheral blood was collected in an EDTA tube prior to surgery, centrifuged twice at 4°C for 10 min (3,000 and 15,000 revolutions per minute). Tissues, plasma, and peripheral blood lymphocytes (PBLs) were stored at −80°C prior to DNA extraction. cfDNA was extracted from plasma samples with the QIAamp Circulating Nucleic Acid Kit (Qiagen, Hilden, Germany). DNA of PBLs was extracted from peripheral blood cells using the QIAamp DNA Blood Mini Kit (Qiagen). Tissue DNA was extracted using the QIAamp DNA Mini Kit (Qiagen). All DNA was extracted according to the manufacturer's instructions. DNA of PBLs was sequenced as the normal control sample.
Targeted Capture, and Next-Generation Sequencing DNA samples were used to construct libraries and perform targeted capture sequencing. To design a capture region for bladder cancer, we selected 861 genes covering recurrent mutations from the Catalogue of Somatic Mutations in Cancer;

Sequencing Data Analysis
Terminal adaptor sequences and low-quality reads were removed from raw data. BWA (version 0.7.12-r1039) was employed to align the clean reads to the reference human genome (hg19). Picard (version 1.98) was used to mark PCR duplicates. Realignment and recalibration were performed using GATK (version 3.4-46-gbc02625). Single nucleotide variants (SNVs) were called using MuTect2 (version 3.5.0). Small insertions and deletions (Indels) were also called using MuTect2. Variants detected in matched control DNA samples with one or more reads indicative of indels at the same location or in the 40-bp flanking regions of experimental samples or residing near regions with low complexity or short tandem repeats were filtered. Somatic copy-number alterations were identified with CONTRA (version 2.0.8) (14). The final candidate variants were all manually verified in the Integrative Genomics Viewer (IGV) (15).

Statistical Analysis
All statistical analyses were performed using GraphPad Prism (version 8.0; GraphPad Software, La Jolla, CA, USA) software. Statistical significance was defined as a two-sided P-value of <0.05. Fisher's exact tests were used to analyze categorical variables. The non-parametric Mann-Whitney test was used to calculate the difference in mTBI between groups. Time to recurrence and follow-up were calculated from the start of immunotherapy. Survival fractions were estimated using the Kaplan-Meier method and differences between groups were identified using the log-rank test. Multivariate Cox proportional hazards analysis was performed using SPSS software (version 22.0; STATA, College Station, TX, USA). The mutation heat map of tumor DNA was plotted using "Complex Heatmap" package in R software (version 3.4.1) (16).

Demographic Information of NMIBC Patients
NMIBC patient characteristics are summarized in Table 1
We analyzed clonality of mutations using PyClone software, which is based on a Bayesian clustering algorithm (17), using variants from 82 patients, which revealed 153 clonal mutations and 323 subclonal mutations. We defined intratumoral heterogeneity as the proportion of mutations that were inferred to be subclonal in a tumor sample. The proportion of subclones ranged from 0 to 94%, with a median of 60%. By using the median proportion of subclones as the cut-off value, patients were divided into high heterogeneity group (the proportion of subclones ≥60%, n = 42) and low heterogeneity group (the proportion of subclones <60%, n = 40). Patients with high heterogeneity had a significantly shorter disease-free survival (DFS) than those with low heterogeneity [median DFS 36 vs. 56 months, hazard ratio (HR) 2.082, 95% confidence interval (CI): 1.137-3.813, p = 0.0204] (Figure 2). Multivariable Cox analysis showed that heterogeneity was an independent prognostic indicator for DFS [HR 1.982 (95% CI: 1.024-3.835), p = 0.0420] ( Table 2).

The Relationship of Clinicopathological Features With ctDNA Detection and Detectability of Tumor-Derived Variations in Plasma
We detected 165 somatic mutations with a median of two per sample (range, 1 to 22) in plasma ctDNA from 54 patients (65.85%) ( Table S2). Tumor DNA could be traced in matched ctDNA from 43 patients (52.44%). We evaluated the association between clinicopathological features and ctDNA detection. ctDNA was defined as positive if somatic mutations could be detected in plasma. The positive rate of ctDNA in T1 patients exceeded that in Ta patients (85.29 vs. 52.08%, p = 0.0021) ( Figure 3A). ctDNA was more frequently detected in patients with tumor size ≥3 cm than in patients with tumor size <3 cm (96.00 vs. 52.63%, p < 0.0001) ( Figure 3C).
The effects of clinicopathological features on detectability of tumor-derived variations in plasma were analyzed. Tumorderived mutations were more often detected in plasma of patients with T1 stage in comparison to that of Ta patients (67.65 vs. 41.67%, p = 0.0257) ( Figure 3D). Patients with tumor size ≥3 cm presented with more tumor-derived mutations in blood than patients with tumor size <3 cm (84.00 vs. 38.60%, p = 0.0002) ( Figure 3F). However, we observed no significant association between tumor grade and ctDNA detection or detectability of tumor-derived variations in plasma ( Figures 3B, E).

Clinicopathological Associations of Detectability of Tumor-Derived Clonal or Subclonal Variations in Plasma
To determine the relationship between clinicopathological features and detectability of tumor-derived clonal or subclonal variations in plasma, we analyzed clonality of mutations using variants from 43 patients with tumor-derived variations in plasma. We found 74 clonal and 202 subclonal mutations in tumor DNA from 43 patients (100%) and 39 patients (90.70%), respectively, while 77 clonal and 68 subclonal mutations were found in ctDNA from 43 patients (100%) and 20 patients (46.51%), respectively. Of 41 clonal mutations from tumor DNA in T1 patients, 24 (58.54%) was detected in matched ctDNA; however, the detection rate for patients with Ta stage was only 27.27% (9/33) (p = 0.0098) ( Figures 4A-C). The detection rates for patients with tumor size ≥3 cm were significantly higher than that in patients with tumor size <3 cm (62.16 vs. 27.03%, p = 0.0047) (Figures 4G-I). The detection rates for HG patients exceeded that for LG patients  Association Between Clinicopathological Features and mTBI mTBI, a measure of molecular tumor burden, was calculated, and the relationship between clinicopathological features and mTBI was analyzed. We found a significant association between mTBI and tumor size. mTBI was significantly higher in patients with tumor size ≥3 cm than in patients with tumor size <3 cm (1.38 vs. 0.42%, p = 0.0020) ( Figure 5).

DISCUSSION
Targeted sequencing has been applied for ctDNA detection in advanced bladder cancer patients (18)(19)(20). For instance, Christensen et al. reported that ctDNA assessment for early risk stratification, therapy monitoring, and early relapse detection in patients with localized advanced bladder cancer was feasible (18). However, targeted sequencing has not been performed in NMIBC patients, partly owing to the low detection of tumor DNA mutations in blood in NMIBC patients. This study is the first to perform targeted sequencing of ctDNA for NMIBC patients. Our cohort included NMIBC patients at stage Ta and T1. Consistent with previous reports (12), we found that ctDNA from plasma samples could be detected in NMIBC patients with Ta stage, with a detection rate of 52.08%. These results showed that DNA may be released into circulation by some transporter, such as exosomes, even though Ta tumors had an intact basal membrane. The detection rate of ctDNA and the concordance rate of somatic mutations between ctDNA and tumor DNA were higher in T1 patients than in Ta patients, suggesting that invasive state of NMIBC could increase release of ctDNA into blood, and ctDNA levels are reflecting invasiveness of disease. We also found that the larger tumor size was significantly associated with the higher ctDNA detection rate and the higher concordance rate of somatic mutations between ctDNA and tumor DNA. These results showed that larger tumor size represented higher tumor burden, leading to increase in ctDNA release into circulation. Ma et al. reported that mTBI had advantage in assessing therapeutic response compared with single gene mutations (21). Wang et al. found that mTBI could serve as a predictive biomarker of treatment outcome for advanced gastric cancer patients receiving chemotherapy (22). A higher mTBI in baseline ctDNA was a risk factor for prognosis. Thus, we evaluated whether mTBI could predict treatment outcome in NMIBC patients receiving TUR of bladder followed by immunotherapy. We found that higher mTBI was associated with shorter DFS in such patients, and could be an independent prognostic indicator for DFS.
Consistent with previous studies (23,24), the pathway analysis of somatic variants for tumor DNA showed that alterations in either FGFR3 or the receptor tyrosine kinase/phosphatidylinositol 3-kinase pathway were identified more frequently in LG-NIMBC in comparison to HG-NIMBC, whereas variations in either TP53 or cell cycle regulation genes were more common in HG disease than in LG disease. These results suggest that the development of LG-NMIBC and HG-NMIBC may be controlled by the receptor tyrosine kinase/phosphatidylinositol 3-kinase pathway and TP53/ cell cycle regulation genes, respectively. Thus, LG-NMIBC and HG-NMIBC patients should be treated with different targeted therapies. For instance, FGFR3 inhibitor may not be effective for HG-NMIBC.
The study by Christensen et al. utilized ddPCR analyses and screened ctDNA for FGFR3 and PIK3CA hotspot mutations in NMIBC (13). This study had an obvious limitation: variants in FGFR3 gene tend to be detected in LG-NMIBC. In our study, FGFR3 and PIK3CA alterations were seen in 28.04 and 26.83% of NMIBC patients, respectively. Thus, the analyses targeting FGFR3 and PIK3CA hotspot mutations can only cover a small proportion of patients with NMIBC.
Liu et al. defined intratumoral heterogeneity as the proportion of mutations that were subclonal in a tumor sample, and found that high post-treatment heterogeneity was associated with A B FIGURE 6 | The relationship between mTBI in ctDNA and disease recurrence after TUR of bladder followed by immunotherapy. Kaplan-Meier curves depict diseasefree survival of 54 patients receiving TUR of bladder followed by immunotherapy, stratified on the basis of mTBI threshold of (A) 0.85% (the median mTBI) or (B) 1%. HR, hazard ratio. Differences between groups were identified using the log-rank test.  (25). We found that high pre-treatment heterogeneity was associated with worse DFS for NMIBC patients receiving TUR of bladder followed by immunotherapy.
In summary, we found that the detection rate of ctDNA and the concordance of somatic variants between tumor DNA and plasma ctDNA were both significantly associated with tumor stage and tumor size. Furthermore, mTBI in ctDNA could serve as an independent prognostic indicator of DFS for NMIBC patients receiving TUR of bladder followed by immunotherapy, and intratumoral heterogeneity in tumor DNA could also independently predict treatment outcome for such patients.

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.

ETHICS STATEMENT
The studies involving human participants were reviewed and approved by Shenzhen Second People's Hospital. The patients/ participants provided their written informed consent to participate in this study.