Clinical utility of the S3-score for molecular prediction of outcome in non-metastatic and metastatic clear cell renal cell carcinoma

Background Stratification of cancer patients to identify those with worse prognosis is increasingly important. Through in silico analyses, we recently developed a gene expression-based prognostic score (S3-score) for clear cell renal cell carcinoma (ccRCC), using the cell type-specific expression of 97 genes within the human nephron. Herein, we verified the score using whole-transcriptome data of independent cohorts and extend its application for patients with metastatic disease receiving tyrosine kinase inhibitor treatment. Finally, we sought to improve the signature for clinical application using qRT-PCR. Methods A 97 gene-based S3-score (S397) was evaluated in a set of 52 primary non-metastatic and metastatic ccRCC patients as well as in 53 primary metastatic tumors of sunitinib-treated patients. Gene expression data of The Cancer Genome Atlas (n = 463) was used for platform transfer and development of a simplified qRT-PCR-based 15-gene S3-score (S315). This S315-score was validated in 108 metastatic and non-metastatic ccRCC patients and ccRCC-derived metastases including in part several regions from one metastasis. Univariate and multivariate Cox regression stratified by T, N, M, and G were performed with cancer-specific and progression-free survival as primary endpoints. Results The S397-score was significantly associated with cancer-specific survival (CSS) in 52 ccRCC patients (HR 2.9, 95% Cl 1.0–8.0, PLog-rank = 3.3 × 10–2) as well as progression-free survival in sunitinib-treated patients (2.1, 1.1–4.2, PLog-rank = 2.2 × 10–2). The qRT-PCR based S315-score performed similarly to the S397-score, and was significantly associated with CSS in our extended cohort of 108 patients (5.0, 2.1–11.7, PLog-rank = 5.1 × 10–5) including metastatic (9.3, 1.8–50.0, PLog-rank = 2.3 × 10–3) and non-metastatic patients (4.4, 1.2–16.3, PLog-rank = 1.6 × 10–2), even in multivariate Cox regression, including clinicopathological parameters (7.3, 2.5–21.5, PWald = 3.3 × 10–4). Matched primary tumors and metastases revealed similar S315-scores, thus allowing prediction of outcome from metastatic tissue. The molecular-based qRT-PCR S315-score significantly improved prediction of CSS by the established clinicopathological-based SSIGN score (P = 1.6 × 10–3). Conclusion The S3-score offers a new clinical avenue for ccRCC risk stratification in the non-metastatic, metastatic, and sunitinib-treated setting. Electronic supplementary material The online version of this article (10.1186/s12916-018-1088-5) contains supplementary material, which is available to authorized users.


Background
Clear cell renal cell carcinoma (ccRCC) is the most common subtype of renal cell carcinoma, with a currently increasing incidence [1][2][3]. Approximately 30% of patients develop metastases and, despite the implementation of targeted therapies, the 5 year survival rate of patients with metastatic disease remains below 20%. Thus, stratification of patients with ccRCC into different molecularly defined groups to identify patients at risk of worse outcome is increasingly important in the perspective of personalized medicine. With this in mind, several prognostic scores have been developed based on, for example, pathological features, gene expression, or DNA methylation status [4][5][6]. One of the most widely applied score established on clinicopathological data is the SSIGN (stage, size, grade, and necrosis) score [7,8], whereas the ClearCode34 score, which predicts two ccRCC subtypes (ccA/ccB), has been suggested for prediction of survival using gene expression data [9,10]. Moreover, Rini et al. [11] proposed a 16-gene score to predict recurrence in ccRCC patients. In general, prognostic signatures using RNA-seq data hold great promise for precision oncology, as previously demonstrated for lung adenocarcinoma [12]. We recently developed an in silico prediction score (named S3-score) for ccRCC, based on the gene expression of 97 signature genes and the similarity of gene expression between tumor cells and their proposed normal cell of origin in the nephron [13]. The S3-score outperforms several other scores [13], including the ClearCode34 model, and significantly improves the predictive value of the SSIGN score and the original ccA/ccB assignment based on clustering [14]. Moreover, compared with the ccA/ccB signature, the S3-score is slightly less dependent on the tumor section investigated [13] and, in consequence displays little intra-tumor heterogeneity. This is of importance because, in a recent study investigating the ccA/ccB signature [10], approximately one-quarter of metastatic tumors (two of nine patients) displayed intra-tumor heterogeneity and, in 43% of the cases, patient-matched primary and metastatic tumors displayed different molecular ccA/ccB subtypes. In this context, a recent multiregion sampling process using a protein-based prognostic model was described, enabling the study of the impact of intra-tumor heterogeneity on risk stratification of sunitinib-treated metastatic patients [15].
As our S3-score was evaluated only in silico using data from The Cancer Genome Atlas (TCGA), we now intended to verify the performance of the score using newly generated whole transcriptome data of an independent cohort of ccRCC patients, including metastases derived from ccRCC, to determine the concordance of the score prediction in primary tumors and ccRCC-derived metastases. Moreover, we evaluated whether the score predicts outcome in sunitinib-treated ccRCC patients. Finally, our objective was to improve the clinical applicability of the S3-score by reducing the number of genes necessary for calculation of the score and by using the more cost-effective real-time PCR technology.

Study cohorts
The study investigated different ccRCC cohorts listed in Table 1 and Additional file 1: Figure S1.
First, our 97 gene-based S3-score (S3 97 ), which was developed using publically available gene expression data of a ccRCC cohort from TCGA (n = 463) ( Table 1) [16] was evaluated in a set of 52 primary ccRCC patients (ccRCC cohort 1). These 52 primary tumor samples were collected from non-metastatic and metastatic patients with ccRCC histology, treated at the Department of Urology, University Hospital Tuebingen, Germany. Patient characteristics are provided in Table 1. The use of the tissue was approved by the ethics committee of the University of Tuebingen and informed written consent was provided by each subject prior to surgical resection. Cancer-specific survival (CSS) was used as the endpoint in the survival analysis of these ccRCC patients.
In addition, publicly available gene expression data from an independent cohort of primary tumors obtained from sunitinib-treated ccRCC patients (n = 53, sunitinib-treated cohort) (Table 1) [17] were used in the analysis. This cohort consisted of ccRCC patients with synchronous or metachronous metastases, who received first-line sunitinib treatment (dosing schedule: 50 mg/ day, 4 weeks on/2 weeks off; at least one 28-day cycle of sunitinib treatment completed; prior cytokine therapy allowed) [17]. Primary ccRCC tissue samples were collected from patients undergoing nephrectomy prior to sunitinib treatment [17]. Further details of these study patients are outlined in Beuselinck et al. [17]. Progression-free survival was used as the endpoint in the survival analysis of sunitinib-treated ccRCC patients.
Next, publicly available gene expression data of TCGA [16] from the cohort of ccRCC patients (n = 463) ( Table 1) were used as a development cohort to define a modified S3-score, which requires a reduced number of genes for clinical application. This S3 15 -score was validated in an extended cohort of 108 metastatic and non-metastatic ccRCC patients treated at the Department of Urology, University Hospital Tuebingen, Germany (extended ccRCC cohort 2, n = 108) ( Table 1). CSS was used as endpoint in the survival analysis. Kaplan-Meier curves of CSS for ccRCC cohorts 1 and 2, as well as for the TCGA cohort are shown in Additional file 1: Figure S2.
In addition, metastases samples (n = 22) derived from 15 patients treated at the Department of Urology, University Hospital Tuebingen, Germany, were collected, including matched primary tumor and metastases samples from five patients of our ccRCC cohorts 1 and 2 (Additional file 1: Table S1 and S5). In part, several regions from one metastasis were collected. Further details about metastases are given in Additional file 1: Table S1 and S5. Use of the tissue was approved by the ethics committee of the University of Tuebingen and informed written consent was provided by each subject prior to surgical resection.
Additional file 1: Figure S1 shows an overview about the workflow of data analyses including the different cohorts and technologies used in the present study.

Gene expression analyses and S3-score calculation
Total RNA was isolated from fresh-frozen ccRCC and metastasis tissue using the mirVana™ miRNA Isolation T primary tumor, N regional lymph node, M distant metastasis present at diagnosis, G grading, NA not available a ccRCC cohort published by Beuselinck et al. [17]; descriptive data were not available b This extended ccRCC cohort 2 includes the 52 patients from ccRCC cohort 1 Kit (Life Technologies) as previously described [18,19]. Genome-wide transcriptome analyses were performed using Human Transcriptome Array HTA 2.0 (Affymetrix) according to the manufacturer's protocol. Further processing of microarray data were performed as previously described [18] (Additional file 1: Supplementary methods). Gene expression data (generated using HuGene 1.0ST Affymetrix array) from 53 sunitinib-treated ccRCC patients were downloaded from ArrayExpress (E-MTAB-3267). Quantitative real-time PCR (qRT-PCR) was performed using TaqMan technology on a BioMARK System (Fluidigm) as described previously [18,19]. TaqMan gene expression assays for 97 genes of the S3-score, as well as five genes used for normalization were purchased from Life Technologies. Further details about calculation of the S3-score based on interprofile correlations and development of a S3-score calculation model for use of qRT-PCR data are provided in the Additional file 1: Supplementary methods.

ClearCode34 and SSIGN calculation
The SSIGN score was calculated as denoted in Zigeuner et al. [8]. The ClearCode34 classifier, as introduced by Brooks et al. [9], was applied on the set of matched primary tumors and metastases of our present cohort for which genome-wide expression data measured by HTA 2.0 microarrays were available (Additional file 1: Supplementary methods).

Statistical analyses
All statistical analyses were performed with R-3.3.3, including additional packages (Additional file 1: Supplementary methods) [20]. Survival analyses for endpoints CSS or progression-free survival were conducted by Kaplan-Meier curves and corresponding log-rank tests as well as uni-and multivariate Cox models. Comparisons of Cox models were performed by analysis of deviance. All statistical tests were two sided. Statistical significance was defined as P < 0.05 (Additional file 1: Supplementary methods).

Evaluation of the S3-score in ccRCC primary tumors and primary tumors of patients treated with sunitinib
We previously developed the S3-score in silico using RNA-seq data from the TCGA (Table 1) [13]. The S3-score was calculated based on 97 genes by correlating tumor expression to the expression in the eight nephron regions. In the present work, we first evaluated this 97 gene-based S3-score (S3 97 ) in our own cohort, consisting of 52 ccRCC samples (ccRCC cohort1) ( Table 1) for which genome-wide expression data using transcriptome arrays were generated. Partitioning of the ccRCC samples by means of the cut-off value that was established in our previous work [13] resulted in two groups with significantly varying CSS (Fig. 1a); i.e. patients with a high S3 97 -score had an decreased risk for cancer-related death compared to patients with low S3 97 -scores. Furthermore, univariate Cox regression including evaluation of the predictive ability according to Harrell's c-index indicated a significant association of the S3 97 -score with patient survival (HR 2.9, 95% Cl 1.0-8.0, P Log-rank = 3.3 × 10 -2 ) (Additional file 1: Table S2). a b Fig. 1 a Cancer-specific survival (CSS) of clear cell renal cell carcinoma (ccRCC) tumors predicted by the S3 97 -score. Kaplan-Meier curves showing CSS of ccRCC cohort 1 (n = 52). Groups are defined by the cut-off of the S3 97 -score, as determined by conditional inference tree models in Büttner et al. [13]. b Validation of the S3 97 -score in sunitinib-treated ccRCC patients. Kaplan-Meier curves showing progression-free survival predicted by the S3 97 -score in an independent cohort of primary tumors of sunitinib-treated patients (n = 53). Groups are defined by the same cut-off value of the S3 97 -score as in a. HR hazard ratio, CI confidence interval. Detailed information on statistical methods is provided in Additional file 1: Supplementary methods Since survival prediction might be majorly influenced by treatment with, for example, sunitinib, we next investigated whether prediction of survival is possible in sunitinib-treated patients. Using a cohort of 53 sunitinib-treated metastatic ccRCC patients with publically available microarray data [17], we calculated the S3-score based on the 97 signature genes. Partitioning of the sunitinib-treated patients by means of our established cut-off value resulted in two groups with significantly varying progression-free survival (Fig. 1b); i.e., patients with a high S3 97 -score had increased progression-free survival probability after sunitinib treatment compared with patients with a low S3 97 -score. Furthermore, univariate Cox regression including evaluation of the predictive ability according to Harrell's c-index, indicated a significant association of the S3 97 -score with patient survival after treatment with sunitinib (HR 2.1, 95% CI 1.1-4.2, P Log−rank = 2.2 × 10 -2 ) (Additional file 1: Table S2).

Refinement of the 97 gene-based S3 97 -score for clinical application
Based on our results, the S3 97 -score has the ability to significantly predict not only CSS in ccRCC patients, but also the progression-free survival in sunitinib-treated individuals. However, calculation of the score was currently based on gene expression data of 97 marker genes, generated through genome-wide transcriptome analyses (RNA-seq or microarray). Thus, for clinical application and utility of the S3-score, expression analyses using quantitative real-time PCR (qRT-PCR), as well as a reduced number of genes, would be more appropriate. Therefore, we aimed to develop a new calculation model of the S3-score. First, the expression of the 97 signature genes, which constituted the basis for the development of the new prediction approach, and the expression of five normalization genes, was quantified by qRT-PCR in our extended ccRCC cohort of 108 non-metastatic and metastatic samples (ccRCC cohort 2) ( Table 1). In order to ensure minimum failure rates in future applications, all assays that failed at least once were excluded. Moreover, we considered only genes that were (after normalization) comparable with respect to mean expression and variation of expression between the RNA-seq data from the TCGA cohort and the RT-PCR values (Additional file 1: Supplementary methods and Additional file 1: Figure S3). In total, the resulting set of variables used for model selection included 41 genes. Subsequently, a linear model using RNA-seq data from the TCGA cohort was created that reconstructs the correlation-based S3-scores.
As expected, a higher incidence of advanced stage tumors as well as metastatic tumors occurred in the S3 15 -low group with poor survival (Additional file 1: Figure S5 and Table S4). Next, we compared the S3 15 -score with clinicopathological prediction factors (T, N, M, G). Multivariate Cox regression indicated that the S3 15 -score is able to significantly improve the predictive ability of the clinicopathological parameters (Table 3). Additionally, the multivariate Cox model outperformed the univariate model (TNMG vs. TNMG+S3: P χ2 = 3.98 × 10 -4 ). Moreover, the S3 15 -score significantly improved CSS prediction when added to the Cox model initially including only the clinicopathologic-based SSIGN score (P χ2 = 1.6 × 10 -3 ) (Fig. 2c).
Evaluation of the S3 97 -score and S3 15 -score in metastases derived from ccRCC patients Tumor heterogeneity of the original S3-score has been previously evaluated to assess whether a single tumor sample is sufficient for prediction of survival [13]. We now aimed to investigate the S3 97 -score and S3 15 -score in metastases in order to evaluate the concordance between primary tumors and metastases. First, we analyzed the S3 97 -score in metastases using microarray data. For a total of 15 ccRCC patients, genome-wide expression data were generated from metastases samples, including five metastatic patients from our ccRCC cohorts with matched primary tumor and metastases samples, as well as three patients for whom several metastases were available. Calculation of the S3 97 -score individually for tumor and metastases resulted in similar risk prediction (Fig. 3a). One patient (P4) was assigned to the high risk group (low S3 97 -score) with worse prognosis using either metastases or tumor tissue (Additional file 1: Table S5). Three patients (P1, P3, P5) showed a high S3 97 -score in tumor as well as metastasis tissue (Additional file 1: Table S5). S3-score was discordant between the primary tumor and its metastasis in only one sample (P2). For different metastases derived from the same patient (P7, P8), as well as for four regions of one metastasis (P6), the S3 97 -score values were also comparable (Fig. 3a). Using the ClearCode34 signature, recently used also for metastatic tissue [10], revealed similar results and classification of tumor/metastases pairs into different molecular subtypes as the S3 97 -score (Additional file 1: Table  S5). We additionally performed qRT-PCR quantification and calculation of the S3 15 -score in a subset of metastases samples. Regarding the five metastatic patients from our ccRCC cohorts for whom primary tumors as well as metastases were available, we found that, except for one case, the S3 97 -score tendency was preserved using the improved S3 15 -score (Fig. 3b).

Discussion
Several risk scores based on gene expression data have been developed for prediction of patient survival in ccRCC [4]. We recently developed a novel prediction score, named the S3-score, based on the similarity of gene expression in the tumor to its cell of origin in the nephron region [13,21]. Thus, in contrast to other scores, risk prediction using the S3-score is related to biologic alterations of the cell of origin of ccRCC. The S3-score outperformed other scores or signatures based on gene expression data or clinicopathological variables [13,21] and was even able to improve the predictive value of the clinically validated SSIGN score [7,8]. Moreover, evaluation of tumor heterogeneity of our S3-score showed that only a few samples displayed heterogeneity [13], which indicates that risk prediction with our score is largely independent from the tumor region investigated.
Generally, most of the scores developed using gene expression data are thus far not introduced into clinical practice because they have not been generated to evaluate individual patients. Thus, for clinical application, the prediction scores need to be validated in several studies defining optimal cut-off values for classification of individual patients into subtypes. Moreover, prediction scores, typically developed using genome-wide gene expression data, need to be evaluated using different technologies and gene expression platforms. Since our S3-score, which is based on the expression of 97 signature genes, was originally developed using RNA-seq data from the TCGA, we first evaluated its predictive ability in the present work using gene expression data generated through microarray technology in our own ccRCC cohort. Here, we showed not only that a platform transfer to microarray data is possible, but also that the S3 97 -score significantly predicts CSS in our cohort.
In contrast to other prediction scores such as the 16-gene signatures [11], the 97 marker genes were not selected based on pathway analyses (e.g., including genes related to inflammation or immune response) and subsequent optimization for prediction of prognosis, but were originally selected to show that tumor aggressiveness in RCC correlates with the level of divergence from its cell of origin within the nephron region. Noticeably, we observe an overlap of one vascular pathway gene (PPAP2B) in the 97 marker genes and those genes from the 16-gene signature described by Rini et al. [11]. Further studies are warranted to compare the predictive ability of both scores.
Because metastases might represent the most aggressive phenotypes of a heterogeneous tumor, herein, we were interested in inter-tumor or metastases heterogeneity, using gene expression data generated by microarray technology once again. Interestingly, the predictive S3 97 -score was comparable between matched tumor and metastases, or matched metastases pairs. In only one case (Additional file 1: Table S5) classification differed between metastases and tumors.
Since data on treatment outcome are limited in the TCGA cohort originally used to develop the S3-score, we were not previously able to evaluate the effect of tyrosine kinase inhibitor (TKI) treatment on outcome prediction. Therefore, herein, we investigated the S3 97 -score using microarray data from a cohort of sunitinib-treated ccRCC patients. In this cohort, the S3 97 -score was significantly associated with progression-free survival of patients, indicating that our score enables even prediction of sunitinib outcome. Whether the same holds true for immunotherapy in the form of T cell immune checkpoint inhibitors like nivolumab needs to be investigated in future studies. Preliminary investigation of the S3 97 -score in metastatic RCC patients treated with nivolumab [22] shows that the S3 97 -score did not differ significantly in pre-and post-treatment biopsies (Additional file 1: Figure S6), indicating that there was no influence of treatment with nivolumab on the S3-score.
Taken together, we provide evidence that the S3 97 -score is more widely applicable than originally intended. To provide a more cost-effective approach for clinical application of the S3-score in individual patient samples, such as even formalin-fixed paraffin-embedded samples, we improved the S3 97 -score by reducing the number of signature genes from 97 to 15 especially for expression analyses through RT-PCR. Our improved a b c Fig. 2 Cancer-specific survival (CSS) of clear cell renal cell carcinoma (ccRCC) tumors of our validation cohort predicted by the simplified RT-PCRbased S3 15 -score. Kaplan-Meier curves showing CSS of (a) the extended ccRCC cohort 2 (n = 108) and (b) the non-metastatic (n = 92) and metastatic subsets (n = 15) of this cohort. Groups are defined by the cut-off of the S3 15 -score, as determined by conditional inference tree models. c S3 15 -score significantly improves the established SSIGN prediction score. χ 2 statistic values depict the improvement of the model likelihood when risk classification based on the S3 15 -score (red) was added to the Cox model initially including only the SSIGN score (blue; left) or vice versa (right). χ 2 test P values are shown in the bars. HR hazard ratio, CI confidence interval. Detailed information on statistical methods is provided in Additional file 1: Supplementary methods S3 15 -score was validated using RT-PCR technology in a cohort of 108 ccRCC cases, clearly indicating that the S3 15 -score was associated with CSS in the complete cohort, as well as non-metastatic and metastatic subsets. Moreover, the S3 15 -score improves prediction of CSS by the currently clinically applied SSIGN score, which is based on clinical parameters and pathologic features. Finally, the S3 15 -score allows risk prediction in tumor and metastases tissue.
In summary, we found that our score enables valid prediction of patient outcome even if applied to different sample types (e.g., primary and metastatic tissue) and independent cohorts (e.g., patients treated with TKIs). Moreover, different platforms (RNA-seq, microarray) and technologies more appropriate for clinical utility (qRT-PCR) can be used for prediction of patient risk by the S3-score. Further prospective studies are warranted to assess the implementation  15 -scores were determined based on gene expression data measured by RT-PCR; T, primary tumor; N, regional lymph node; M, distant metastasis present at diagnosis; G, grading; tumors with grade "G2-3" and "G4" were added to "G3". Tumors with no grading information or metastasis status "MX" were disregarded a b