Cell-Free DNA Variant Sequencing Using Plasma and AR-V7 Testing of Circulating Tumor Cells in Prostate Cancer Patients

Prostate cancer (PCa) is the second most common malignant cancer and is a major cause of morbidity and mortality among men worldwide. There is still an urgent need for biomarkers applicable for diagnosis, prognosis, therapy prediction, or therapy monitoring in PCa. Liquid biopsies, including cell-free DNA (cfDNA) and circulating tumor cells (CTCs), are a valuable source for studying such biomarkers and are minimally invasive. In our study, we investigated the cfDNA of 34 progressive PCa patients, via targeted sequencing, for sequence variants and for the occurrence of CTCs, with a focus on androgen receptor splice variant 7 (AR-V7)-positive CTCs. The cfDNA content was associated with overall survival (OS; p = 0.014), disease-specific survival (DSS; p = 0.004), and time to treatment change (TTC; p = 0.001). Moreover, when considering all sequence variants grouped by their functional impact and allele frequency, a significant association with TTC (p = 0.017) was observed. When investigating only pathogenic or likely pathogenic gene variants, variants of the BRCA1 gene (p = 0.029) and the AR ligand-binding domain (p = 0.050) were associated with a shorter TTC. Likewise, the presence of CTCs was associated with a shorter TTC (p = 0.031). The presence of AR-V7-positive CTCs was associated with TTC (p < 0.001) in Kaplan–Meier analysis. Interestingly, all patients with AR-V7-positive CTCs also carried TP53 point mutations. Altogether, analysis of cfDNA and CTCs can provide complementary information that may support temporal and targeted treatment decisions and may elucidate the optimal choice within the variety of therapy options for advanced PCa patients.


Introduction
Prostate cancer (PCa) is a major cause of disease and mortality among men worldwide each year, and 1.4 million men are diagnosed with PCa. Approximately 375,000 men died in 2020 of PCa [1]. PCa is recognized as a genetically heterogeneous disease, comprising a large scope of malignancies, from indolent localized cancers that may never progress to rapidly progressing castration-resistant PCa (CRPC) [2]. In accordance with the recent knowledge of clinicopathologies and genetics, the World Health Organization's (WHO) classification of prostatic cancers has been revised, and five grade groups have been Abbreviations: None-no pretreatment; Depriv-hormone deprivation; AE-antiandrogens abiraterone/enzalutamide; Tax-taxanes; Mult-multiple (both AE and Tax). Chemotherapy/hormone following None; AE following Depriv; AE following Tax; AE following Mult; Tax following Depriv; Tax following AE; Tax following Mult.

Sampling of Blood, Isolation of CTCs, and Processing of Plasma
Blood samples were collected via venipuncture in K 3 EDTA anti-coagulation tubes (Sarstedt AG & Co KG, Nümbrecht, Germany), stored at 4 • C, and processed within 4 h after the blood was drawn. Plasma was prepared by centrifugation at 2000× g for 10 min and stored as aliquots at −80 • C. Cell-free DNA was isolated from 4-8 mL plasma via affinity-based binding to magnetic beads (QIAamp MinElute ccfDNA Kit, QIAGEN, Hilden, Germany), according to the manufacturer's instructions, and was quantified via fluorometry, using a Quantus fluorometer (Promega, Maddison, WI, USA). Plasma cfDNA concentrations were normalized to the plasma input volume and are given in ng/mL.
Library preparation was performed using a QIAseq Targeted DNA Panel Library Prep Kit. An amount of 10 ng, or the maximum when 10 ng of starting material was not available, was enzymatically fragmented. Ends were repaired and 3 adenylated. Adapters were ligated to the overhang, and the reactions were cleaned up. Target regions were enriched via target-specific PCR (22 cycles). After bead-based reaction cleanup, molecules were enriched via 22 cycles of universal PCR, and a final cleanup was performed. Library preparation was quality controlled using capillary electrophoresis (Agilent DNA 7500 Chip). High-quality libraries were pooled, based on equimolar concentrations using a bioanalyzer automated electrophoresis system (Agilent Technologies). The library pool(s) were quantified using qPCR, and the optimal concentration of the library pool was used to generate the clusters on the surface of a flow cell before sequencing on a NextSeq (Illumina Inc., San Diego, CA, USA) instrument (2 × 151 bp, 2 × 8 bp), according to the manufacturer's instructions (Illumina Inc., San Diego, CA, USA).

Data Analysis and Bioinformatics Analysis
Primary bioinformatic analysis was carried out with the CLC Genomic Workbench. Read quality control was facilitated with the "QC for Sequencing Reads" tool. Adapters were trimmed by the "Trim Reads" tool, and alignment, UMI grouping, and variant calling were facilitated by the "Identify QIAseq DNA Somatic Variants" workflow.
Finalized variant files were introduced to the Qiagen Gene Globe data analysis application. By using the QCI translational application (QIAGEN Clinical Insight Interpret 8.0.20210827, Hilden, Germany) with precurated databases for somatic variations, the results were further refined. Briefly, all common SNPs with a known allele frequency above 1% were omitted. Furthermore, any sequence variations located within intronic sequences more than 20 bp from the exon-intron border were also omitted. Additionally, sequence variants with a known pathogenic function were explicitly included in the analysis. According to the QCI translational precured databases, all sequence variants were classified by their biological function (loss of function, gain of function, or normal function), as well as their pathological impact (pathogenic, likely pathogenic, benign, likely benign, or uncertain impact). The final variation set included a total of 2505 sequence variations distributed over 98 genes.

Isolation of mRNA and qRT-PCR
CTCs were isolated from 5 mL of whole blood via positive immunomagnetic selection by applying the AdnaGen prostate cancer select system (QIAGEN GmbH, Hilden, Germany), targeting EpCAM and PSMA. Cellular mRNAs were isolated via immunomagnetic enrichment (oligo dT(25)-coated magnetic beads) and reverse transcribed. The mRNAs of CD45, PSA, PSMA, AR, AR-V7, and GAPDH were quantified using the AdnaTest Prostate-CancerPanel AR-V7 (Qiagen, Hilden, Germany). PSMA and GAPDH detection were considered for verification of CTCs, and additional detection of AR-V7 was considered for verification of AR-V7-positive CTCs. Exclusive detection of CD45 or GAPDH (or both) was considered as verification of lymphocytes.

Statistical Analysis
Correlations between the detection of cfDNA, CTCs, and sequencing results and clinicopathological data were calculated using Spearman's correlation test, a chi-squared test, or a Mann-Whitney test. The associations of cfDNA, CTCs, AR-V7-positive CTCs, and sequencing results with overall survival (OS), disease-specific survival (DSS), and time to treatment change (TTC) were determined by univariate analyses (Kaplan-Meier analysis with log-rank test and Cox's regression hazard models). A p value less than 0.05 was considered statistically significant. Statistical analyses were performed using the SPSS 21.0 software package (SPSS Inc., Chicago, IL, USA). For hierarchical clustering, the Euclidean distance and average linkage were used, and the analysis was performed using R statistical framework Ver. 3.2.1 (R Foundation for Statistical Computing, Vienna, Austria. http://www.R-project.org/ (accessed on 20 October 2021)).

Results
We studied 39 peripheral blood samples from 34 PCa patients with progressive PCa for plasma cfDNA sequence variants. Out of these for 37 samples (33 patients), the presence of CTCs was analyzed with a focus on CTCs expressing the androgen splice variant AR-V7 (AR-V7-positive CTCs). An overview of the clinicopathological data, treatment schemes, and study data is given in Table 1. The patients were grouped according to their pretreatment into five groups, as follows: no pretreatment; basic hormone deprivation (e.g., leuprorelin); antiandrogens abiraterone/enzalutamide (AE); taxanes; or multiple (both AE and taxanes). In addition, the treatment (performed subsequent to blood sampling and prior to therapy change) consisted of three schemes: AE, taxanes, or chemohormone therapy. According to pretreatment and treatment, we separated the patients into seven groups, as follows: chemohormone therapy without pretreatment; AE following hormone deprivation; AE following taxanes; AE following multiple pretreatments; taxanes following hormone deprivation; taxanes following AE; and taxanes following multiple pretreatments. The first clinical end point was the time to treatment change (TTC), which was defined as the interval between blood sampling (including the scheduled therapy scheme) and the initiation of the following therapy scheme. Other clinical endpoints, including the interval between diagnosis and the last date of patient information, were obtained to define overall survival (OS) and disease-specific survival (DSS). In the Kaplan-Meier analysis, we studied the association between these pretreatment-treatment groups and prognosis, i.e., OS, DSS, and TTC ( Figure 1).
Concerning the underlying prognostic characteristics, we discovered in a pairwise comparison that patients treated with AE following multiple pretreatments had poorer survival (both OS and DSS) than patients receiving taxanes following multiple pretreatments (both p = 0.043) or patients receiving taxanes following AE (both p = 0.017; Figure 1).
Kaplan-Meier analysis of the association between the pretreatment-treatment groups and TTC was performed with the chemohormone therapy group as the reference. Here, all groups, with the exception of patients receiving AE following hormone deprivation, showed a significantly shorter TTC, i.e., AE following taxanes (p = 0.034), taxanes following hormone deprivation (p = 0.004), taxanes following AE (p = 0.019) and taxanes following multiple pretreatments (p = 0.028; Figure 1).
Next, we were interested in whether a liquid biopsy-i.e., analysis of cfDNA, CTCs, and AR-V7-positive CTCs-could provide further information concerning molecular changes and the prognosis of PCa patients. Concerning the underlying prognostic characteristics, we discovered in a pairwise comparison that patients treated with AE following multiple pretreatments had poorer survival (both OS and DSS) than patients receiving taxanes following multiple pretreatments (both p = 0.043) or patients receiving taxanes following AE (both p = 0.017; Figure  1).
Kaplan-Meier analysis of the association between the pretreatment-treatment groups and TTC was performed with the chemohormone therapy group as the reference. Here, all groups, with the exception of patients receiving AE following hormone deprivation, showed a significantly shorter TTC, i.e., AE following taxanes (p = 0.034), taxanes following hormone deprivation (p = 0.004), taxanes following AE (p = 0.019) and taxanes following multiple pretreatments (p = 0.028; Figure. 1).
Next, we were interested in whether a liquid biopsy-i.e., analysis of cfDNA, CTCs, and AR-V7-positive CTCs-could provide further information concerning molecular changes and the prognosis of PCa patients.   Table 2). The level of cfDNA was correlated with the follow-up time after blood sampling (r s = −0.333; p = 0.047) and the presence of AR-V7positive CTCs (r s = 0.346; p = 0.036) tended to be correlated with the TTC (r s = −0.304; p = 0.071; Table 3). However, the concentration of cfDNA was not correlated with age, PSA level, or the pretreatment scheme (data not shown).  Next, we stratified patients, according to the mean of their cfDNA concentration, into two groups (≤ 17.93 ng/mL vs. > 17.93 ng/mL). Prognosis was analyzed via Kaplan-Meier analysis ( Figure 1). A higher level of cfDNA was associated with a shorter OS (p = 0.014; 42.2 vs. 134.0 months) and a shorter DSS (p = 0.004; 42.2 vs. 138.7 months; Table 4). Next, we were interested in whether the level of cfDNA was associated with the TTC. A higher level of cfDNA was associated with a significantly shorter TTC (p = 0.001; 6.3 vs. 16.2 months) than that observed in patients with a lower level of cfDNA ( Figure 2; Table 4). In univariate Cox's regression analysis, the patients with higher cfDNA levels had a 5.53-fold increased risk of death (p = 0.027) and a 7.87-fold risk of disease-specific death (p = 0.013; Table 5), compared with patients with lower levels of cfDNA. Likewise, the risk for an earlier treatment change in the group with higher cfDNA levels was 4.11-fold increased (p = 0.004) compared with the group with lower levels of cfDNA (Table 5).   Significant values are in bold. Estimated mean survival is given in months.   Targeted sequencing of cfDNA was performed for a total of 99 genes. After applying the data filtering schemes described in the methods section, a total of 2505 sequence variations were identified in our patient cohort. Of these, the classification was as follows: 1835 uncertain variants, 164 likely benign variants, 31 benign variants, 387 likely pathogenic variants, and 61 pathogenic variants (Table S1; Figure S1). The occurrence of individual sequence variations in the 39 samples differed from a single detection in one sample for 2049 variants, two independent detections for 129 variants, and three independent detections for 54 variants. A total of 273 variants were detected in four or more of the samples. Considering the pathological impact-i.e., only the likely pathogenic and pathogenic variants-448 variants distributed on 68 genes were identified. At a frequency of ≥2; 375 variants were found distributed on 66 genes. By combining the pathological impact (likely pathogenic or pathogenic) and the frequency of occurrence, 64 variants distributed on 22 genes were detected (Table S1). The variants comprised SNVs, insertions or deletions of variable length, and substitutions of ≥2 base pairs. Interestingly, most of the variants resulted in a loss of function and were related to genes considered tumor suppressors or DNA repair genes, or both. Taken together, all samples analyzed likely harbored pathogenic and pathogenic variants.
Considering the patient with three sequential samples analyzed, it was of interest that the number of detected variants decreased with tumor progression. However, pathogenic variants that occurred in the first sample were retained in the second and third samples, e.g., in the PTEN gene (c.916G > A and c.926G > A), in the CHEK2 gene (c.1556C > T), and in the TP53 gene (c.857A > C).

Frequency of Variants
Of the pathogenic and likely pathogenic variants (n = 424), a total of 56 (13.2%) were located in the MUC16 gene and 33 variants (7.8%) were located in SYNE1 (spectrin repeatcontaining nuclear envelope; synonymous: Nesprin1), and with 32 variants (7.5%), the AR gene contained the third highest variant frequency. However, when we normalized the frequency of pathogenic and likely pathogenic variants per covered sequence (in kb), the TP53 gene, with a ratio of 10.85 (15 variants/1.383 kb), showed the highest ratio; the AR gene, with a ratio of 10.65 (32 variants/3.005 kb), exhibited the second highest ratio; and the CHEK2 gene, with a ratio of 7.85 (15 variants/1.911 kb), exhibited the third highest ratio ( Figure S2); whereas, MUC16, with a ratio of 1.26, and SYNE1, with a ratio of 1.18, presented rather moderate ratios (Table S1).

Association between Sequence Variants and Prognosis
Considering only pathogenic and likely pathogenic variants for single genes, several genes, such as TP53, AR, RAD51, and RB, did not show an association with prognosis (OS, DSS, and TTC; data not shown). However, patients with BRCA1 pathogenic and likely pathogenic variants, but not those with BRCA2 variants, showed a shorter TTC (p = 0.029; Figure 3); whereas, OS or DSS were not associated with pathogenic or likely pathogenic variants in either gene. Keup et al. reported that three variants in the MUC16 gene in metastatic breast cancer were associated with a longer survival time after diagnosis of metastasis [37]. Interestingly, in our study, patients with likely pathogenic variants for the MUC16 gene showed a longer OS and a longer survival time after diagnosis of metastasis (p = 0.028 and p = 0.025; Figure 4). We also tested the large gene SYNE1, but pathogenic and likely pathogenic variants were not associated with prognosis.
Cells 2021, 10, x FOR PEER REVIEW 10 of 20 (p = 0.028 and p = 0.025; Figure 4). We also tested the large gene SYNE1, but pathogenic and likely pathogenic variants were not associated with prognosis.  Through a more detailed analysis, Conteduca et al. found that two AR gene mutations in the ligand binding domain (LBD; 2105T > A (p. L702H) and 2632A > G (p. T878A)) were associated with a shorter OS in PCa [33]. We detected these two AR point mutations in four patients, and two of these patients possessed both mutations simultaneously. In addition, two patients showed other AR-LBD gene variants (Table S1). In our study, for all patients with AR-LBD mutations, an association with a shorter survival time after (p = 0.028 and p = 0.025; Figure 4). We also tested the large gene SYNE1, but pathogenic and likely pathogenic variants were not associated with prognosis.  Through a more detailed analysis, Conteduca et al. found that two AR gene mutations in the ligand binding domain (LBD; 2105T > A (p. L702H) and 2632A > G (p. T878A)) were associated with a shorter OS in PCa [33]. We detected these two AR point mutations in four patients, and two of these patients possessed both mutations simultaneously. In addition, two patients showed other AR-LBD gene variants (Table S1). In our study, for all patients with AR-LBD mutations, an association with a shorter survival time after Through a more detailed analysis, Conteduca et al. found that two AR gene mutations in the ligand binding domain (LBD; 2105T > A (p. L702H) and 2632A > G (p. T878A)) were associated with a shorter OS in PCa [33]. We detected these two AR point mutations in four patients, and two of these patients possessed both mutations simultaneously. In addition, two patients showed other AR-LBD gene variants (Table S1). In our study, for all patients with AR-LBD mutations, an association with a shorter survival time after blood sampling (p = 0.036) and a trend toward a shorter TTC were found (p = 0.050; Figure 5), when compared with patients without mutations in the AR-LBD. blood sampling (p = 0.036) and a trend toward a shorter TTC were found (p = 0.050; Figure  5), when compared with patients without mutations in the AR-LBD. To assess both the functional impact of detected variants and the frequency of variant alleles, we calculated a gene wise and sample variation score. Hereby, the inferred functional impact of the variant was rated as 0 (variant not present), 1 (normal function), or 2 (loss of function or gain of function). These ratings were multiplied by the allele frequency of the variant allele. To calculate a gene wise, size-normalized score, the variation values were summarized for the observed gene and divided by the sequence length assessed by targeted sequencing. Hierarchical clustering revealed three distinct clusters ( Figure 6). The three clusters only tended to be correlated with OS (rs = −0.300; p = 0.095) or DSS (rs = −0.345; p = 0.053; Table 3). However, the three clusters showed a significant association with TTC (p = 0.017; Figure 7). To assess both the functional impact of detected variants and the frequency of variant alleles, we calculated a gene wise and sample variation score. Hereby, the inferred functional impact of the variant was rated as 0 (variant not present), 1 (normal function), or 2 (loss of function or gain of function). These ratings were multiplied by the allele frequency of the variant allele. To calculate a gene wise, size-normalized score, the variation values were summarized for the observed gene and divided by the sequence length assessed by targeted sequencing. Hierarchical clustering revealed three distinct clusters ( Figure 6). The three clusters only tended to be correlated with OS (r s = −0.300; p = 0.095) or DSS (r s = −0.345; p = 0.053; Table 3). However, the three clusters showed a significant association with TTC (p = 0.017; Figure 7).

Isolation of CTCs and Correlation with Clinicopathological and Molecular Data
Out of the 39 samples, 37 samples were also available to assess the presence of CTCs and the expression of AR-V7. In 11 samples, we detected CTCs, and in 4 of these samples, AR-V7-positive CTCs were found (Table S2). The presence of CTCs was correlated with the variant groups (r s = 0.472; p = 0.004; Table 3). The presence of AR-V7-positive CTCs was correlated with the time of blood sampling to follow-up (r s = −0.377; p = 0.028) and TTC (r s = −0.378; p = 0.028). However, neither the presence of CTCs nor that of AR-V7-positive CTCs was correlated with AR-LBD mutations (data not shown).
TP53 gene alterations are the most commonly occurring gene alterations in PCa [38]. TP53 alterations can be differently associated with prognosis in cancer [39]. Accordingly, we separated our patients into 3 groups according to their pathogenic and mostly pathogenic TP53 gene variants: (Group 1) without gene variants, (Group 2) with frameshift mutations, and (Group 3) with point mutations. The 3 TP53 variant groups were associated with AR-V7-positive CTCs (r s = 0.523; p = 0.001; Table 3) and trended to be associated with the occurrence of CTCs (r s = 0.309; p = 0.063). Interestingly, CTCs occurred predominantly in the 2 groups with TP53 sequence variants (7/11) and in 5 cases in Group 3 (with TP53 point mutations; Table S3). All 4 samples with AR-V7-positive CTCs were found in Group 3 (p < 0.001; Table S3).

Associations between CTCs and Prognosis
Based on Kaplan-Meier analysis, the presence of CTCs was not associated with OS or DSS but was associated with a shorter TTC (p = 0.031; 7.6 vs. 15.1 months; Table 4; Figure 8) compared with patients with no CTCs. In univariate Cox's regression analysis, patients with CTCs possessed a 2.32-fold increased risk for a shorter TTC (p = 0.046; Table 5) compared with patients without CTCs.

Associations between CTCs and Prognosis
Based on Kaplan-Meier analysis, the presence of CTCs was not associated with OS or DSS but was associated with a shorter TTC (p = 0.031; 7.6 vs. 15.1 months; Table 4; Figure 8) compared with patients with no CTCs. In univariate Cox's regression analysis, patients with CTCs possessed a 2.32-fold increased risk for a shorter TTC (p = 0.046; Table  5) compared with patients without CTCs. The presence of AR-V7-positive CTCs was also not associated with OS but showed an association trend with DSS (p = 0.097; 57.7 vs. 134.8 months) and was significantly associated with TTC (p < 0.001; 4.0 vs. 14.2 months; Table 4; Figure 8); however, it must be considered that only 4 samples from 3 patients showed AR-V7-positive CTCs. In the univariate Cox's regression analysis, the presence of AR-V7-positive CTCs was associated The presence of AR-V7-positive CTCs was also not associated with OS but showed an association trend with DSS (p = 0.097; 57.7 vs. 134.8 months) and was significantly associated with TTC (p < 0.001; 4.0 vs. 14.2 months; Table 4; Figure 8); however, it must be considered that only 4 samples from 3 patients showed AR-V7-positive CTCs. In the univariate Cox's regression analysis, the presence of AR-V7-positive CTCs was associated with a 14.75-fold increased risk of treatment change (p = 0.001) but not with an increased risk of disease-specific death (Table 5).

Discussion
We studied a group of PCa patients with progressive disease in different pretreatment and treatment stages. As expected, the chemohormone therapy group, being in the beginning of therapy, showed the longest OS, DSS, and TTC. In contrast, but equally expected, the AE group following multiple pretreatments presented with the shortest OS, DSS, and TTC. However, we were interested in whether molecular characterization of cfDNA and CTCs could provide insight into associations with the prognosis of these PCa patients. We performed targeted cfDNA sequencing for sequence variant identification and AR-V7 testing of CTCs from the patients.
Quantitative assessment of cfDNA can be used as a diagnostic biomarker and to assess tumor burden and response to therapy [40]. In our study, the amount of cfDNA found (approximately 18 ng/mL plasma) was in the range found previously in many other cancer types, including PCa (average 18 ng/mL to 31.9 ng/mL) [41,42]. However, the level was somewhat below what is described on average for CRPC [41]; however, it is well known that the amount of cfDNA is prone to change during treatment [43]. The concentration of cfDNA was associated with the prognosis of PCa patients. A higher cfDNA concentration-i.e., above the mean value (17.93 ng/mL)-was associated with a shorter OS and DSS, which is in accordance with previous findings of other authors [27,29,31,34,[43][44][45]. In addition, higher levels of cfDNA were associated with a shorter TTC. However, our cfDNA concentrations are not baseline values and therefore mainly reflect therapy response. Several studies report a reduction in cfDNA concentrations between baseline and after therapy for responders, e.g., after chemotherapy [44] or PARP inhibition [34].
However, in addition to the amount of cfDNA, the mutational spectrum has relevance for prognosis and therapy decisions [46][47][48]. Considering only the pathogenic and likely pathogenic variants, all patients showed such variants in their cfDNA, and as expected, the TP53 gene and AR gene exhibited these variants most frequently in our study. This finding is not unexpected; Robinson et al. showed that mutations in the TP53 and AR genes are the most frequent in mCRPC [38], and cf/ctDNA has been demonstrated to be representative of multiple metastases [30,32]. Recently, TP53 status in cfDNA of CRPC patients has been reported as an independent prognostic factor for progression-free survival [49]. However, TP53 mutations can be an early event in some PCa, while in others, an enrichment in metastases can be found; thus, their role in malignant progression needs further clarification [50]. In our study, patients with BRCA1 pathogenic and likely pathogenic variants but not those with BRCA2 variants showed a shorter TTC. Our finding is comparable to the results of Kohli et al., who found that mutations in multiple DNA repair genes (ATM, BRCA1, BRCA2, and CHEK2) were associated with time to ADT failure and survival in metastasized hormone-sensitive PCa [45]. Two AR gene mutations in the ligand binding domain (2105T > A (p. L702H) and 2632A > G (p. T878A)) have been associated with resistance to abiraterone therapy and with a shorter OS in PCa [33]. Mutations in the AR-LBD region are also associated with shorter progression-free survival (PFS) [48]. In our study, the presence of four mutations in the LBD of the AR gene, including the two previously mentioned mutations, tended to be associated with a shorter TTC than observed in patients without mutations in the AR-LBD.
Interestingly, patients with likely pathogenic variants in the MUC16 gene showed a longer OS and a longer time from metastasis to last follow-up. This finding is comparable to the results of Keup et al., who found that three MUC16 variants were associated with a longer time from metastasis to last follow-up in metastatic breast cancer patients [37]. To test whether this is an effect of the large gene size of MUC16, we also analyzed the large-sized gene SYNE1. However, pathogenic and likely pathogenic variants in the SYNE1 gene were not associated with the prognosis of our PCa patients. In gastric cancer, the occurrence of MUC16 gene mutations is correlated with tumor mutational burden and associated with better OS [51]. Recently, Yu et al. showed an association between immune checkpoint inhibitor therapy and survival (OS, PFS) in patients with non-small cell lung cancer with MUC16 variants [52].
When considering all variants, i.e., their functional impact and allele frequency, they were associated with TTC. Group 3, with the highest rate, showed the shortest time for treatment change. A short TTC was correlated with poor OS and DSS in our study ( Table 3). Considering that our patient cohort consisted of primarily late-stage progressive disease patients, a short TTC was mainly related to acute therapy failure. However, more studies with larger cohorts are necessary to investigate the sequence of therapies and their effect on gene variants and resistance mechanisms.
In addition, we identified several pathogenic and likely pathogenic gene variants in genes that are of potential relevance for the following existing therapies: (i) DNA repair genes, such as CHEK2, BRCA1/2, ATM, MSH2, MSH6, PALB2, RAD51, ERCC3, MRE11, and NBN; (ii) genes in the PIK3-AKT pathway, such as PIK3CA and PTEN; and (iii) genes in the CDK4/6-RB pathway, such as RB1 and CDKN2A (p16). The occurrence of mutations in DNA repair genes is related to PARP inhibitor therapy, and PTEN deficiency is implicated in AKT-PIK3 pathway blockade, e.g., with ipatasertib in combination with abiraterone [32,34,45,[53][54][55][56]. Aberrations in the CDK4/-6-RB1 pathway may suggest application of CDK4/-6 inhibitors [46]. In addition, a clinical activity of pembrolizumab has been shown in metastatic PCa patients with a high microsatellite instability detected in ctDNA [57]. Altogether, our results and previous findings in the literature indicate that characterization of gene variants in cfDNA is relevant for prognosis, therapy prediction, and therapy monitoring.
We also analyzed the occurrence of CTCs and AR-V7-positive CTCs in our PCa patient cohort. We detected CTCs in 29.7% (11/37) of samples and AR-V7-positive CTCs in 10.8% (4/37) of samples. The number of samples with CTCs is expected to be higher in reality, because several PCa patients have undetectable CTCs, despite progressive disease [58], and CTCs with epithelial to mesenchymal transition usually escape CTC detection [13]. In our study cohort, the occurrence of CTCs and AR-V7-positive CTCs was associated with a shorter TTC (p = 0.031 and p < 0.001). The prognostic role of CTCs in PCa has been clearly shown, as recently reviewed [9]. It has been reported that the prognosis of PCa patients worsens from patients without CTCs to patients with AR-V7-negative CTCs to patients with AR-V7-positive CTCs [59]. Detection of the androgen receptor splice variant AR-V7 has been shown to confer resistance to abiraterone and enzalutamide in PCa [21]. Detection of AR-V7 in CTCs at the RNA and protein levels was associated with shorter PFS and OS after abiraterone or enzalutamide treatment [16]. However, the detection of AR-V7 does not preclude a response to next-generation androgen deprivation therapy [23]. Recently, a study showed that early CTC declines were good predictors of improved outcomes in mCRPC patients treated with docetaxel [60]. Altogether, CTC count has been described as a reliable biomarker for treatment response and prognosis in patients receiving chemotherapy or AR-targeting therapies [61]. Interestingly, in our study, the occurrence of AR-V7-positive CTCs was associated with TP53 point mutations. AR-V7 splicing is regulated by the histone demethylase JMJD1A. JMJD1A has been found to interact with and to promote the recruitment of a member of the heterogeneous nuclear ribonucleoprotein (HNRNP) family, HNRNPF, to the cryptic exon 3b on androgen receptor pre-mRNA for generation of AR-V7 [62]. Another member of the HNRNP family, HNRNPK, is known to function as a cofactor for TP53, e.g., in the DNA damage response [63]. HNRNPK is also well known as a regulator of androgen receptor translation [64]. Recently, mutant TP53 has been shown to have an effect on RNA splicing via the RNA binding protein HNRNPK in pancreatic cancer [65]. Both HNRNPK and HNRNPF are part of the splicing machinery [66]. It is tempting to speculate that a mutated TP53 protein may play a role in the splicing of AR-V7 via HNRNPs.
We found that the level of cfDNA and the occurrence of AR-V7-positive CTCs were correlated and that both were associated with TTC. In addition, the molecular analyses of cfDNA and CTCs also provided complementary information, such as genetic variants and the presence of AR-V7-positive CTCs, which can support therapy decisions. This is in accordance with findings in the literature. As recently reviewed, on the one hand, the assessment of genomic aberrations in cf/ctDNA can potentially predict therapy response and detect mechanisms of resistance, and on the other hand, the presence of CTCs can be a surrogate for disease stage or progression, and can measure therapy response [40,61].
Our study has several limitations. Although the targeted resequencing panel covers 99 genes, more genes are expected to show PCa-relevant gene variants. With the applied CTC detection technique, we were unable to enumerate CTCs or to follow changes in CTC numbers during therapy. In addition to gene variants, splicing variants and epigenetic mechanisms-such as DNA methylation, histone modifications, and expression of noncoding RNAs-can also affect gene or protein expression and may reflect the therapy response, which remained unconsidered in our study. In addition, standardization and the clinical validation of liquid biopsy-based assays to detect predictive biomarkers of resistance is essential.
However, the finding that there are associations of gene variants in cfDNA or CTC occurrence with TTC is of clinical interest. This finding can aid in understanding therapy resistance mechanisms and drive therapeutic decisions in the future, ultimately improving outcomes for PCa patients. Considering that we found variants that are associated with a longer OS and a longer survival time after diagnosis of metastasis (MUC16), it is also conceivable that knowledge of specific variants could lead to less overtreatment or fewer later treatment changes, avoiding stress and new side effects in patients.
Supplementary Materials: The following data are available online at https://www.mdpi.com/ article/10.3390/cells10113223/s1. Table S1: Overview of sequence variants, Table S2: Detection of CTCs and AR-V7+ CTCs with the Adna-Test ProstateCancerPanel AR-V7, Table S3: Cross table: Relationship between TP53 gene variants and CTCs or AR-V7-positive CTCs, Figure S1A: Classification of called variants in cfDNA according to their known functional impact, Figure S1B: Classification of called variants in cfDNA according to variation type, Figure S2: Frequency of pathogenic and likely pathogenic variants per sequence range for the top 15 genes.