Prognostic value of CA20, a score based on centrosome amplification-associated genes, in breast tumors

Centrosome amplification (CA) is a hallmark of cancer, observable in ≥75% of breast tumors. CA drives aggressive cellular phenotypes such as chromosomal instability (CIN) and invasiveness. Thus, assessment of CA may offer insights into the prognosis of breast cancer and identify patients who might benefit from centrosome declustering agents. However, it remains unclear whether CA is correlated with clinical outcomes after adjusting for confounding factors. To gain insights, we developed a signature, “CA20”, comprising centrosome structural genes and genes whose dysregulation is implicated in inducing CA. We found that CA20 was a significant independent predictor of worse survival in two large independent datasets after adjusting for potentially confounding factors. In multivariable analyses including both CA20 and CIN25 (a gene expression-based score that correlates with aneuploidy and has prognostic value in many types of cancer), only CA20 was significant, suggesting CA20 captures the risk-predictive information of CIN25 and offers information beyond it. CA20 correlated strongly with CIN25, so a high CA20 score may reflect tumors with high CIN and potentially other aggressive features that may require more aggressive treatment. Finally, we identified processes and pathways differing between CA20-low and high groups that may be valuable therapeutic targets.

the most abundant proteins in the centrosome in several cell lines 11 ) (Supplementary Table 1), along with TUBG1, which encodes the most abundant centrosomal protein and is primarily responsible for microtubule nucleation, key to centrosomal function 11 . Our objective was to test the prognostic value of CA20 after adjusting for potentially confounding factors in multiple breast cancer cohorts and to explore processes, pathways, and oncogenic signatures that are associated with a high CA20 score. Because CA causes CIN, we were also interested in comparing the prognostic value of CA20 with that of the CIN score "CIN25", which correlates with total functional aneuploidy and predicts worse outcomes in a variety of cancers 12 , determining which of these two scores has the most significant impact on outcomes when included together in multivariable models of survival, and comparing processes, pathways, and oncogenic signatures that are enriched in tumors with high CA20 and CIN25 scores.

Results
We tested the ability of CA20 and CIN25 to risk-stratify breast cancer patients in two datasets, the METABRIC and TCGA datasets, comprising n = 1,969 primary breast cancers with breast cancer-specific survival (BCSS) annotation and n = 524 primary invasive breast cancers with OS annotation, respectively. The METABRIC dataset was split into discovery and validation sets. The TCGA dataset was not split because power analysis suggested the subsets would potentially be too small, so bootstrapping was instead used to obtain more reliable estimates of population parameters. METABRIC dataset. Stratification was conducted according to average CA20 and CIN25 scores found in the discovery set as well as optimal cutpoints in CA20 and CIN25 scores found in the discovery set based on the log-rank test. In Kaplan-Meier plots, stratification into high-and low-BCSS groups based on the average and optimal cutpoints in CA20 and CIN25 scores was significant in both the discovery and validation sets (p < 10 −6 for all, Fig. 1; see Tables 1 and 2 for descriptive statistics of study datasets). When both CA20 and CIN25 (both stratified by the average score) were entered as covariates in full multivariable models using discovery set data, only CA20 (stratified by the average score) appeared in the final model, and it was a significant predictor of BCSS (Hazard Ratio [HR] = 2.88, p < 0.001; Table 3). In the validation set too, CA20 (stratified by the average score) remained a significant predictor in the final model (HR = 2.13, p < 0.001). Common significant covariates between discovery and validation set final models included tumor stage and chemotherapy. When both CA20 and CIN25 (both stratified by the optimal cutpoint) were entered as covariates in full multivariable models using discovery set data, both covariates appeared in final models but only CA20 (stratified by the optimal cutpoint) significantly affected BCSS (HR = 2.13, p = 0.006; Table 4). In the validation set, CA20 (stratified by the optimal cutpoint) remained a significant predictor (HR = 1.82, p = 0.028), whereas CIN25 (stratified by the optimal cutpoint) did not significantly impact BCSS. Common significant covariates between discovery and validation final multivariable models of BCSS included tumor stage and chemotherapy, as was found when stratifying by average signature scores. Thus, CA20 (whether stratified by the average score or optimal cutpoint) is a significant predictor of BCSS after adjusting for stage and chemotherapy, whereas CIN25 (whether stratified by the average score or optimal cutpoint) is not an independent predictor in these models.

Figure 1.
Plots of Kaplan-Meier product limit estimates of breast cancer-specific survival of patients in METABRIC discovery and validation sets stratified by (A,E) CA20 (average value), (B,F) CA20 (optimal threshold), (C,G) CIN25 (average value), and (D,H) CIN25 (optimal threshold), respectively. Average values and optimal thresholds were determined using the discovery set.   Table 3. Final multivariable Cox proportional-hazards models of breast cancer-specific survival including CA20 and CIN25 (stratified by average scores) in full models using METABRIC data. HR = Hazard Ratio; CI = Confidence Interval. CA20 score was highly correlated with CIN25 score (ρ = 0.93, p < 10 −6 ), which may reveal that breast tumors with high CA20 scores have high levels of CIN. Although breast cancer subtype was not a common independent predictor of outcomes, we were interested to test whether CA20 and CIN25 scores differed grade-wise between TNBCs and non-TNBCs, which differ in aggressiveness. No grade 1 TNBCs were present in the dataset for comparison, but we found that average CA20 and CIN25 scores were higher in TNBCs than non-TNBCs in both grade 2 and 3 tumors per two-tailed independent samples t-tests, equal variances assumed (p < 0.001 for all) (Fig. 2), consistent with the more aggressive behavior of TNBCs compared with non-TNBCs and mirroring what we previously found for CA as assessed by microscopy 3 .

TCGA dataset.
To confirm the prognostic value of CA20 in a separate cohort, we analyzed the TCGA breast dataset. Stratification was conducted according to the average CA20 and CIN25 scores found in the entire dataset as well as optimal cutpoints in CA20 and CIN25 found in the entire dataset based on the log-rank test. In Kaplan-Meier plots, stratification into high-and low-OS groups based on CA20 (average and optimal cutpoints) was significant (p = 0.025 and p = 0.024, respectively), with high CA20 conferring a worse prognosis. For comparison, stratification by CIN25 (optimal cutpoint) was significant (p = 0.029), whereas stratification by CIN25 (average cutpoint) was not (Fig. 3). In stage-adjusted models, high CA20 scores (based on both average and optimal cutpoints) were associated with 2.72-and 2.79-fold worse OS (bootstrap-p = 0.016 and 0.008, respectively). For comparison, in stage-adjusted models, high CIN25 scores (based on both average and optimal cutpoints) were also associated with worse OS, HR = 2.31 and 4.65 (bootstrap-p = 0.026 and 0.035, respectively). However, as in the METABRIC dataset, when both CA20 and CIN25 (stratified by average or optimal cutpoints) were entered along with stage in full models, following backward variable selection (based on an α = 0.10 removal criterion), only CA20 and stage remained as predictors in full models (CA20 [stratified by average cutpoint]: HR = 2.80, p = 0.008; CA20 [stratified by optimal cutpoint]: HR = 2.55, p = 0.009), and they remained significant following bootstrapping (Table 5).  Table 4. Final multivariable Cox proportional-hazards models of breast cancer-specific survival including CA20 and CIN25 (based on optimal thresholds) in full models using METABRIC data. HR = Hazard Ratio; CI = Confidence Interval. Although age at diagnosis was not a significant predictor of BCSS in the METABRIC dataset in final models, we recognized the possibility that it could confound analyses of OS in this independent dataset. We thus refit multivariable models entering CA20 and CIN25 (stratified by average or optimal cutpoints), AJCC stage, and age at diagnosis. In final models, CA20 remained a significant predictor, along with stage and age but not CIN25, and the hazard associated with high CA20 was even greater than in models not adjusted for age (CA20 [stratified by average cutpoint]: HR = 3.82, p = 0.001; CA20 [stratified by optimal cutpoint]: HR = 3.67, p = 0.002); furthermore, significance was retained after bootstrapping (Table 6). Because all the cases in the multivariable models were annotated with age at diagnosis, our sample size and, thus, statistical power were not diminished. Therefore, the prognostic value of CA20 adjusting for stage was upheld in this separate dataset adjusted for confounding variables, suggesting broad clinical utility for this score to predict outcomes in female breast cancer patients. Similar to our findings in the METABRIC dataset, in the TCGA dataset CA20 was very strongly correlated with CIN25 (ρ = 0.95, p < 10 −6 ), suggesting that breast tumors with high CA20 scores also have high levels of CIN.
Finally, we were interested in exploring differences in biological processes, molecular pathways, and oncogenic signatures between CA20-high and low groups (defined by the average CA20 value), which may reveal potentially actionable biology. To this end, we performed Gene Set Enrichment Analysis (GSEA) 13 using the    Table 6. Final multivariable Cox proportional-hazards models of overall survival including CA20 and CIN25 (based on average or optimal thresholds), AJCC stage, and age at diagnosis in full models using TCGA data along with p-values and 95% confidence intervals (CIs) for final model covariate hazard ratios obtained via simple bootstrapping. References for hazard ratios (HRs) are CA20 (low) and stage I/II. TCGA dataset and explored differentially enriched biological processes, Reactome pathways, and oncogenic signatures. For the CA20-high group, 262 biological process gene sets were enriched at false discovery rate (FDR) q < 0.05 (Supplementary Table 2), but no such gene sets were significantly enriched in the CA20-low group (Supplementary Table 3). Among the most significant results, the CA20-high group was enriched in DNA repair processes, the DNA integrity checkpoint, many cell cycle processes (e.g., mitotic nuclear division, cell cycle phase transition, cell division, spindle assembly, regulation of sister chromatid segregation, mitotic spindle organization), and regulation of microtubule polymerization/depolymerization. Regarding Reactome pathways, the CA20-high group was enriched in 96 gene sets at FDR q < 0.05 (Supplementary Table 4), but the CA20-low group was not enriched in any such gene sets (Supplementary Table 5). Top enriched Reactome pathways in the CA20-high group exhibited much overlap with biological processes, including DNA repair and cell cycle pathways. For the purposes of comparison, we also compared biological processes and Reactome pathways between CIN25-high and low groups (stratified by the average CIN25 score) (Supplementary Tables 6-9), and it was found that results overlapped greatly with those from the CA20 analyses. Regarding enriched biological processes, only one gene set found in the CA20-high group (<1% of gene sets) was not found in the CIN25-high group, and only 14 gene sets found in the CIN25-high group (~5% of gene sets) were not found in the CA20-high group (Supplementary Table 10). Similarly, only four Reactome pathways enriched in the CA20-high group and four in the CIN25-high group differed (~4% of gene sets) (Supplementary Table 12). We also explored differences in oncogenic signature gene sets. We found that the CA20-high group was enriched in 13 such gene sets, including genes upregulated upon overexpression of E2F1, stimulation with sonic hedgehog (SHH) protein, and loss of retinoblastoma protein (pRb) (Supplementary Table 13), whereas the CA20-low group was not significantly enriched in any of such gene sets (Supplementary Table 13). The CIN25-high group was enriched in 12 oncogenic signature gene sets (Supplementary Table 14), all of which were found in the CA20-high group, whereas no such gene sets were significantly enriched in the CIN25-low group (Supplementary Table 15). These data suggest the CA20-and CIN25-high groupings may capture rather similar molecular tumor profiles. Finally, we tested whether the CIN25-high group was enriched in the centrosome gene ontology cellular component, and we found that it was at FDR q < 0.001 (Normalized Enrichment Score = 2.29), suggesting this group is enriched in centrosomal genes, consistent with the strong correlation we found between CA20 and CIN25 scores.

Discussion
CA is a well-characterized hallmark of cancer 14 , especially breast cancer. Indeed, ≥75% of breast tumors (ductal carcinomas in situ, adenocarcinomas, invasive ductal carcinomas, or breast tumors not otherwise specified) exhibit CA 1 . Because CA promotes CIN and other aggressive phenotypes, it may be a driving force in tumorigenesis and tumor evolution that can offer insights into the clinical course of breast tumors, but only a few studies have investigated the potential prognostic value of CA. Our lab previously developed a four-gene signature, the CAI, which we demonstrated could stratify n = 162 breast cancer patients into two groups with significantly different clinical outcomes in Kaplan-Meier analysis, with high CAI based on an optimal cutpoint correlating with worse OS 3 . In the same study, we also found a non-significant trend among n = 120 breast cancer patients towards worse progression-free survival (PFS) in Kaplan-Meier analysis for tumors with high levels of CA (defined as the sum of the percentage of cells with >2 centrosomes and the percentage of cells with abnormally voluminous centrosomes based on microscopy, using an optimal cutpoint). Another study of n = 362 breast tumors found that large centrosomal size was not associated with OS or recurrence-free survival (RFS) after adjusting for tumor stage and subtype in multivariable Cox models; however, it is not known whether 2D (i.e., cross-sectional) measurements reliably estimate centrosome size, given that centrosomes are 3D structures, and numerical CA was not considered in multivariable analyses 15 . In the same study, however, Kaplan-Meier analyses revealed that high numerical CA (defined as >2 centrosomes per cell on average) was associated with worse BCSS, OS, and RFS. Thus, there is limited evidence that CA may be associated with worse outcomes in breast cancer, but it is unclear what impact CA has on survival after adjusting for potential confounders and what biological processes and pathways could be targeted therapeutically in tumors with high levels of CA.
To shed light on these questions, we developed the CA20 score based on genes encoding centrosome structural proteins and genes that have been demonstrated to induce CA following experimental perturbations in their expression. As we found for CA previously 3 , CA20 (and CIN25) were higher in the aggressive TNBC subtype than non-TNBCs in grade-matched comparisons. In analyses of two large and well-annotated breast cancer datasets (the METABRIC and TCGA breast datasets), we found that high CA20 score was associated with worse BCSS and OS after adjusting for potentially confounding factors, suggesting that CA20 could be a useful clinical tool to identify breast cancer patients at greater risk of poor outcomes. When both CA20 and CIN25 were factored into multivariable models, only CA20 was significantly associated with outcomes. This finding suggests that when CA20 is accounted for CIN25 no longer holds prognostic value. Given that we found a very strong correlation between CA20 and CIN25 in breast tumors and it has been shown by others that CA and CIN are correlated in breast tumors 15 , it is tempting to speculate that CA20 captures CIN, thus rendering CIN25 redundant, and perhaps also captures other aggressive phenotypes not encompassed by CIN25 that are consequences of CA. Given that CIN engenders karyotypic diversity within tumors, we assert that CA20 may perhaps even serve as an indirect measure of intratumor heterogeneity in breast tumors. The overlap in biological processes, Reactome pathways, and oncogenic signatures that are enriched in CA20-and CIN25-high groups is striking given that the two signatures only share one gene in common (CDK1) and suggests that they reflect relatively similar molecular tumor biology (namely, potential activation of DNA repair pathways, perhaps to cope with DNA damage occurring due to chromosome missegregation, enhanced cell cycle kinetics and microtubule dynamics, and activated E2F1 signaling), although perhaps with subtle but prognostically important qualitative and quantitative distinctions.
An exciting avenue for future research would be to test whether breast tumors with high CA20 are more susceptible to E2F1 or SHH inhibitors, drugs targeting DNA repair mechanisms (e.g., PARP inhibitors), chemotherapeutics that target the cell cycle (e.g., taxanes), or centrosome declustering drugs (such as griseofulvin, noscapinoids, PJ34, and KifC1/HSET inhibitors), which preferentially eliminate cells with CA by forcing them to construct a multipolar spindle during mitosis [16][17][18][19][20] . Because most normal cells do not have amplified centrosomes, declustering drugs exhibit low to no apparent toxicity to them. It will also be important to validate (through careful microscopy and rigorous quantitation) that CA20 scores indeed correlate with CA in breast tumors in future studies.

Methods
Dataset details and power analyses. Microarray datasets were chosen based on their availability in Oncomine 21 and the presence of annotation regarding survival time (measured in days) and statuses and signature gene expression levels. Three microarray datasets met these criteria, including the METABRIC 22 , TCGA 23 , and Esserman 24 breast datasets; however, power analysis suggested the Esserman dataset was too small, so it was excluded from analyses (see "Esserman dataset" below). The clinical data and log 2 median-centered signature gene expression levels of the METABRIC and TCGA datasets were thus downloaded from Oncomine. METABRIC dataset: Normal breast, benign breast neoplasms, and cases without BCSS annotation were excluded from analyses, resulting in a sample size of n = 1,969 primary breast cancers. A majority of the cases were annotated for AJCC stage and whether adjuvant chemotherapy was given. The dataset was then split randomly (via random number assignment) and approximately equally into discovery and validation sets (n = 985 and n = 984, respectively; Table 1, descriptive statistics). Neither significant differences nor non-significant trends (i.e., 0.05 < p < 0.10) were found between these two sets for continuous variables (age, CA20 score, and CIN25 score; 2-tailed t-tests), ordinal variables (Nottingham grade, tumor stage, CA20 group [optimal], CA20 group [average], CIN25 group [optimal], and CIN25 group [average]; Mann-Whitney tests), or nominal variables (breast cancer subtype, chemotherapy, radiotherapy, and hormone therapy; Chi-square tests) (data not shown). TCGA dataset: Normal breast specimens, metastases, and male breast cancers were excluded from analyses, resulting in a sample size of n = 524 primary invasive breast cancers ( Table 2, descriptive statistics). OS annotation was incomplete, so we supplemented it with clinical data downloaded from the TCGA data portal, after which all cases had OS time and status. 395 cases had AJCC stage annotation, but information about adjuvant chemotherapy was not available. We analyzed the METABRIC data to estimate whether this sample size would potentially achieve statistical power ≥0.80 in a study of the effect of CA20 on OS in stage-adjusted models with an average follow-up time of approximately 3 years, as in the TCGA study. Among METABRIC patients with invasive breast cancers (n = 1,030), the overall probability of an event (death) within 3 years was p E = 0.12, the probability of belonging to the CA20 (optimal)-high group was p H = 0.70, and the relative risk of death was HR = 2.34. Based on these data, it was estimated that a sample size of n = 339 would be needed to detect a HR = 2.34 with a Type I error rate of α = 0.05 and Type II error rate of β = 0.20, based on the formula to calculate the one-sided sample size in Cox proportional-hazards models 25 45 . In addition, the gene encoding the primary centrosome structural protein, TUBG1, was included, resulting in a set of 20 genes. The CIN25 gene signature is described in Carter et al. 12 . For both datasets, many genes were represented by multiple probes. To select the probe most likely to represent the gene, probes were filtered by rational selection processes that differed by dataset since they are based on different platforms (the Illumina HumanHT-12 V3.0 R2 Array for the METABRIC dataset and the Agilent custom 244 K for the TCGA dataset). The CA20 and CIN25 scores were calculated as the sum of the normalized (log 2 median-centered) expression levels of the signature genes. METABRIC dataset: Probes were filtered by preferentially selecting those targeting all isoforms (A designation). For some genes, only probes targeting some isoforms (S designation) were available. In addition, probes mapping only to the gene of interest according to a BLAST-like Alignment Tool (BLAT) search against the reference genome GRCh38 using Ensembl 46 were preferentially selected. When multiple probes mapped to the gene of interest, the average expression level was calculated to represent that gene. TCGA dataset: In the absence of A and S designations, probes were filtered by performing a BLAT search as for the METABRIC data, also with averaging of normalized expression levels when multiple probes mapped to the gene of interest. Scores exhibited negative values, so for ease of interpretation, scores were converted to non-negative values by adding the minimum score value to all scores, which did not alter the results of statistical analyses.
Scientific RepoRts | 7: 262 | DOI:10.1038/s41598-017-00363-w Survival analyses. Stratification of cases into high-and low-survival groups both by the average (as performed in the CIN25 analyses previously 12 ) and by an optimal cutpoint (based on the most significant log-rank test statistic found using Cutoff Finder 47 ) per the Kaplan-Meier method. Prior to fitting Cox proportional-hazards models, the proportional-hazards assumption was tested by defining each covariate as a function of time and entering this time-dependent term into a simple Cox model and determining whether there was a significant hazard in the discovery and validation sets. For no covariate was the assumption violated (data not shown). Spearman correlation (2-tailed) was performed to determine the correlation between CA20 and CIN25 scores. IBM SPSS Statistics version 21 was used for all analyses, and p < 0.05 was considered statistically significant. METABRIC dataset: Multivariable Cox models were fit using both discovery set data via backward-stepwise elimination of covariates (subject to an α = 0.10 removal criterion) and validation set data by entering final discovery model covariates. Full multivariable discovery model covariates included age at diagnosis (years), Nottingham grade (1, 2, or 3), AJCC stage (0, I, II, or III/IV, the latter two categories combined due to the relatively small number of stage IV cases), breast cancer subtype (luminal: ER and/or PR+, HER2−; HER2-enriched: ER/PR+/−, HER2+; triple-negative: ER/PR/HER2−), chemotherapy (yes/no), hormone therapy (yes/no), and radiotherapy (yes/no). Also, depending on the analysis, the full model either contained CA20 and CIN25 categorized based on the average score as found in the discovery set or CA20 and CIN25 categorized based on the optimal cutpoint as found in the discovery set. TCGA dataset: To confirm the prognostic ability of CA20 in an independent dataset, multivariable Cox models were fit using TCGA data by entering CA20 and CIN25 (average or optimal, depending on the model) and AJCC stage (categorized as I/II vs. III/IV due to relatively low numbers of stage I and IV patients) into full models (chemotherapy information was not available). Covariates were then subjected to backward-stepwise elimination (α = 0.10 removal criterion). Multivariable models were also fit including age at diagnosis (years).
To more robustly estimate population parameters, the final model covariates were entered into Cox models with simple bootstrapping (1,000 iterations).
Grade-wise comparison of average CA20 and CIN25 between TNBCs and non-TNBCs. Using the METABRIC dataset, as grade information was not available for the TCGA dataset, we compared average CA20 and average CIN25 between TNBCs and non-TNBCs grade-wise using two-tailed independent samples t-tests, guided by F-test results, and p < 0.05 was considered statistically significant.