Key Points

• Seventeen upregulated genes, ESR1, CCND1, SCUBE2, PGR, ERBB4, THBS1, GATA3, BCL2, TBC1D9, THSD4, STC2, CCDC170, STK32B, NBEA, PLAT, IL6ST, and NAT1 were identified as the genes associated with resistance and connected with other nodes in the network analysis.

• ESR1, CCND1, and SCUBE2 emerged as the top three key genes associated with resistance to taxane-based neoadjuvant therapy.

• SCUBE2 displayed a high predictive power comparable to ESR1, and better than CCND1, the two commonly accepted markers in breast cancer.

• The predictive ability of SCUBE2 was significantly high in estrogen receptor-positive and human epidermal growth factor receptor 2-positive breast cancers.

Introduction

Breast cancer is the most common cancer in the world and the leading cause of death from cancer in women (1). The most common histopathologic subtypes of breast cancer are invasive ductal carcinoma, invasive lobular carcinoma, and mixed ductal/lobular carcinomas. Increased knowledge on the expression status of estrogen receptor (ER), progesterone receptor (PR), and human epidermal growth factor receptor 2 (HER2) allowed molecular classification of breast cancer into molecular subtypes: hormone receptor-positive, HER2 positive, and triple-negative breast cancer (TNBC) (2, 3).

Chemotherapy, radiotherapy, and surgery are the mainstays of treatment for breast cancer. Neoadjuvant chemotherapy became the standard of care in early breast cancer, increasing the chances of breast-conserving surgery and allowing total tumor resection in locally advanced breast cancer. Neoadjuvant chemotherapy also introduces the possibility to tailor the adjuvant chemotherapy regimens based on the response to the regimes used in the neoadjuvant setting (2, 3).

The molecular classification of breast cancer enabled tailoring of the neoadjuvant chemotherapy, based on the expression status of receptors and improving treatment outcomes. Today, endocrine therapy is an important contributor to therapeutic success in hormone receptor-positive breast cancer, and the incorporation of anti-HER2 antibody trastuzumab improved treatment outcomes substantially in HER2-positive breast cancers. Despite that, taxane-based chemotherapy is still critical in the initial management of breast cancer patients (2, 3).

Taxanes, such as paclitaxel or docetaxel, are the main components of taxane-based neoadjuvant chemotherapy in breast cancer. They act by blocking the depolimerization of microtubules, thus inhibiting cell proliferation (4). Different chemotherapeutics were adopted into the taxane-based neoadjuvant regimens over the years. Usually, three months of taxane chemotherapy is followed or preceded by anthracycline-based therapy. Formerly, the general practice was to combine 5-fluorouracil and cyclophosphamide with an anthracycline (adriamycin or epirubicin) in anthracycline-based chemotherapy (3, 5). However, current guidelines mostly recommend adriamycin/cyclophosphamide or epirubicin/cyclophosphamide, since 5-fluorouracil does not increase the therapeutic efficacy significantly (6).

The success of taxane-based neoadjuvant chemotherapy is high in HER2-positive and TNBCs. However, the number of ER+/HER2- breast cancer patients that benefit from taxane-based neoadjuvant chemotherapy are limited. Moreover, the annualized recurrence rate in breast cancer for the first 5 years was calculated as 10.4%, and the risk of recurrence continues beyond 20 years, despite the combination of high-efficacy chemotherapeutics with distinct mechanisms of action in chemotherapy (7). Therefore, further investigation of the molecular markers responsible for chemoresistance is critical in breast cancer. Such investigations may offer new drug targets for overcoming resistance to taxane-based chemotherapy.

A great interest to identify predictive markers for chemotherapy and endocrine therapy led to the development of several risk score tests, such as the Oncotype DX Breast Recurrence Score test of 21 genes, MammaPrint 70-gene assay, and EndoPredict 12-gene Molecular Score. These risk score tests were validated as predictors of response to neoadjuvant therapy in ER+/HER2- negative breast cancers (8, 9). They were proven to be successful in selecting patients who will benefit from neoadjuvant chemotherapy and decreased the cost of breast cancer management in countries including the United Kingdom and Germany, where these tests are reimbursed by public health insurance (10, 11). Despite that, these tests are not incorporated into the routine management and health insurance systems in many developing countries, which limits their use (12). Moreover, the therapeutic potential of targeting the genes included in these scoring tests has not been completely addressed.

The aim of this study was to identify markers of resistance to taxane-based neoadjuvant therapy, with predictive power and therapeutic potential comparable to ER, the most reliable marker in breast cancer. Such markers may guide therapeutic decision making, especially in countries where risk score tests are not incorporated into routine care. To this end, we analyzed three comprehensive breast cancer cohorts with gene expression profiling data by robust bioinformatics tools. We included studies which utilized pathological complete response (pCR) as the surrogate of responsiveness to taxane-based neoadjuvant chemotherapy, since the number of patients in studies which utilized relapse free survival as the surrogate of responsiveness to taxane-based neoadjuvant chemotherapy was significantly lower.

Materials and Methods

Patient Data and Identification of Differentially Expressed Genes

Three GEO expression profiling datasets (GSE20194, GSE25066, GSE32646) (13, 14, 15), which include gene profiling data for breast cancer patients that have undergone taxane-based neoadjuvant chemotherapy were analyzed in the study (https://www.ncbi.nlm.nih.gov/geo/). Affymetrix Human Genome U133A Array was utilized in GSE20194, and GSE25066 datasets, and Affymetrix Human Genome U133 Plus 2.0 Array was utilized in GSE32646. To determine genes associated with chemoresistance to taxane-based chemotherapy, the patients were stratified into pCR and residual disease (RD) groups. Pathological complete response denoted patients without residual cancer in the breast and lymph nodes. Only patients for which information on the pathological response was available were included in the analysis.

In the GSE20194 dataset, the patients received paclitaxel for three months followed by 5-fluorouracil, cyclophosphamide, and adriamycin before surgery (13). The GSE25066 dataset included samples from patients who received a taxane for three months followed by 5-fluorouracil, adriamycin, and cyclophosphamide (FAC), or 5-fluorouracil, epirubicin, and cyclophosphamide (FEC); or received four cycles of adriamycin and cyclophosphamide followed by four cycles of taxane (14). In the GSE32646 dataset, the samples were from primary breast cancer patients who had undergone neoadjuvant chemotherapy with weekly paclitaxel for three months followed by FEC every three weeks for three months (15). The neoadjuvant treatment protocols and distribution of the patients based on the ER/HER2 and pCR/RD status in each dataset are listed in Table 1.

From the GSE20194 dataset, we included 261 samples from patients who received taxane-based neoadjuvant chemotherapy followed or preceded by anthracycline-based therapy (FEC or FAC). In this dataset, eight patients received anti-HER2 therapy, and three patients received endocrine therapy in the neoadjuvant setting in addition to taxane-based neoadjuvant chemotherapy. These patients were not included in our analysis. Additionally, four patients who received only FEC or FAC, one patient who received only taxol and one patient for which information was not available about therapy, were excluded from the analysis. The patients in the GSE25066 dataset all received taxane and anthracycline based neoadjuvant chemotherapy. ER positive patients received endocrine therapy in the adjuvant setting but not in the neoadjuvant setting. GSE32646 dataset included patients who all received taxane and anthracycline based neoadjuvant chemotherapy. The authors did not mention any use of anti-HER2 therapy or endocrine therapy in these patients. The three datasets did not include patients who received any platinum chemotherapeutics.

The differentially expressed genes (DEGs) in the RD group were identified using the GEO2R web tool (https://www.ncbi.nlm.nih.gov/geo/geo2r/). Log transformation and force normalization were applied to the data. The Benjamini & Hochberg (False discovery rate) method was used to adjust the p-values (adjusted p-value significance cut-off = 0.05). The genes were filtered based on their log2-fold change (logFC) values. The genes with the log FC value >0.5 were accepted as the upregulated genes and with the log FC value <-0.5 were accepted as the downregulated genes. Then Venn Analysis was performed on jvenn (an interactive Venn diagram viewer) (http://jvenn.toulouse.inra.fr/app/index.html) to detect DEGs common to all three datasets.

Functional Annotation and Enrichment Analysis

To identify the gene ontologies and pathways that the DEGs were enriched, the upregulated or downregulated gene lists in non-responsive patients were analyzed on The Database for Annotation, Visualization, and Integrated Discovery (DAVID) (Version 6.8) (16) (https://david.ncifcrf.gov/). To understand the gene ontologies (GO-CC: Cellular compartments, GO-MF: Molecular functions, and GO-BP: Biological processes) and pathways at which the DEGs enriched, gene ontology (GO) analysis and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analysis were performed on DAVID (p-value significance cut-off = 0.05).

Protein-protein Interaction Network Analysis

To identify the key markers associated with chemoresistance to taxane-based chemotherapy, the protein-protein interaction (PPI) network of the protein products of the DEGs was constructed on The Search Tool for the Retrieval of Interacting Genes (STRING Version 11.0) (17) (https://string-db.org/). The minimum required interaction score was set as medium confidence (0.4). The resulting PPI network was analyzed in Cytoscape (Version 3.8.2) and the topological parameters of each protein were extracted (18) (https://cytoscape.org/).

Gene Expression Profiling and Receiver Operating Characteristic Curve Analysis

To validate the upregulation of key genes in breast cancer patients who did not respond to taxane-based therapy, the gene expression profiles of breast cancer patients were analyzed on the ROC-plotter (19). Mann–Whitney test was used for statistical analysis. To validate the value of the key markers in predicting resistance to taxane-based chemotherapy, the receiver operating characteristic (ROC) curves of the key genes were analyzed in ROC-plotter (19). A total of 872 breast cancer patients (228 responders and 644 non-responders) who received taxane-based chemotherapy were included in both analysis types. To exclude the confounding effects that may be caused by other therapies, we did not include patients who received endocrine therapy, anti-HER2 therapy or carboplatin but only included patients who received taxane-based chemotherapy in ROC curve analysis. Pathological complete response was considered as the criteria for responsiveness or non-responsiveness in both analysis types.

Results

Identification of Genes Associated with Resistance to Taxane-Based Neoadjuvant Chemotherapy

To identify the genes associated with resistance to taxane-based neoadjuvant chemotherapy in breast cancer patients, we analyzed GSE20194, GSE25066, and GSE32646 datasets in GEO2R. In total, 182 samples from patients with the pCR and 699 samples from patients with the RD were analyzed. Table 1 lists the number of samples with pCR or RD in each dataset.

First, we investigated the genes differentially expressed in primary breast tumors from patients unresponsive to taxane-based neoadjuvant chemotherapy compared with the patients that responded to therapy (Figure 1a-c). We identified 60 common genes differentially expressed in unresponsive patients in all three datasets (Figure 1d). Among these, 39 DEGs were upregulated and 21 DEGs were downregulated in patients with RD (Figure 1e-f).

Functionally Enriched Pathways and Gene Ontologies Associated with Resistance to Taxane-Based Neoadjuvant Chemotherapy

To identify the pathways and gene ontologies at which the DEGs were enriched, we analyzed the gene lists for commonly upregulated genes and downregulated genes in DAVID. The 39 genes upregulated in all three datasets generated eight statistically significant clusters. ADCY1, APBB2, BCL2, CCND1, ESR1, GATA3, IL6ST, NBEA, PGR, and TSPAN1 were the constituents of the top cluster among these eight clusters (Table 2). These genes were enriched in four KEGG pathways (chemical carcinogenesis-receptor activation, endocrine resistance, estrogen signaling pathway, and pathways in cancer) and two biological processes (response to xenobiotic stimulus and response to drug). The 21 genes downregulated in all three datasets did not generate clusters in functional enrichment analysis. Therefore, we continued further investigation with the upregulated genes list.

Determining The Key Hub Genes Associated with Chemoresistance

To determine the key players in resistance to taxane-based neoadjuvant chemotherapy in breast cancer, we analyzed the PPI network of the 39 common genes upregulated in resistant breast tumors at all datasets on Cytoscape (Figure 2). Among the protein products of the 39 upregulated genes, 17 proteins exhibited connectivity with at least one other protein in the PPI network. The topological parameters for these 17 connected proteins in the network are listed in Table 3. Analysis of the network with Cytoscape revealed ESR1 (estrogen receptor 1), CCND1 (cyclin D1), and SCUBE2 (signal peptide-CUB-epidermal growth factor–like domain-containing protein 2), as the top three central genes associated with resistance.

Validating The Predictive Power of Key Hub Genes as Markers of Resistance to Taxane-Based Neoadjuvant Chemotherapy

Among the three markers identified in our analysis, ESR1 and CCND1 are already known to be associated with resistance to anthracycline-taxane-based neoadjuvant chemotherapy, endocrine therapy, and immune checkpoint inhibitors (20, 21). However, there are contradictory findings on the role of SCUBE2 in breast cancer, and its predictive capacity for chemoresistance in breast cancer is not clear (22, 23, 24).

To compare the predictive power of SCUBE2 with that of ESR1, and CCND1 in breast cancer patients who underwent taxane-based neoadjuvant chemotherapy, we analyzed the differential expression of these genes in non-responders versus responders and the ROC plots in a validation set of 228 responders and 644 non-responders (Figure 3). All three genes were upregulated in breast cancer patients, who did not respond to taxane-based chemotherapy. SCUBE2 displayed the highest fold increase in non-responders compared to ESR1 and CCND1. The ROC analysis of these genes indicated that the predictive power of SCUBE2 can be as high as ESR1 and better than CCND1, as a marker of resistance to taxane-based chemotherapy.

A recent study detected SCUBE2 as one of the four drug resistance markers in ER-positive breast cancers (25). However, the authors have not limited their test cohort to patients who received taxane-based neoadjuvant therapy but included patients treated with any neoadjuvant modality. To assess whether SCUBE2 has different predictive power for taxane-based therapy resistance in distinct subtypes of breast cancer, we investigated the differential expression and ROC plots of SCUBE2 in non-responders versus responders who had different molecular subtypes of breast cancer. In our analysis, SCUBE2 displayed the highest fold change in HER2-positive breast cancers, compared to ER-positive/HER2-negative and TNBC subtypes, being insignificant in TNBC (Figure 4). However, the sensitivity of SCUBE2 as a marker for resistance to taxane-based therapy was highest in ER-positive/HER2-negative breast cancer. These findings suggested that SCUBE2 can be used as a novel marker with predictive strength comparable to ESR1 in ER-positive/HER2-negative and HER2-positive breast cancers.

Discussion and Conclusion

Advancements in the molecular dissection of cancers pave the way for personalized medicine in cancer therapy. Evidence on the key role of ER and HER2 in breast cancer progression enabled the incorporation of endocrine and anti-HER2 therapies with chemotherapy regimens in breast cancer management (2, 3). Despite these approaches improving treatment outcomes significantly, nearly 30% of breast cancer patients experience recurrence due to resistance to chemotherapy (26). Therefore, further dissection of the molecular markers responsible for chemoresistance is of critical importance in breast cancer.

In this study, we analyzed three breast cancer cohorts with gene expression profiling data, using up-to-date bioinformatics tools to identify key markers of resistance to taxane-based neoadjuvant therapy, which is the standard of care in early breast cancer patients. We identified 39 genes upregulated in breast cancer patients who did not respond to taxane-based neoadjuvant therapy. These genes were highly enriched in gene ontologies and KEGG pathways such as “response to xenobiotic stimulus”, “chemical-carcinogenesis-receptor activation”, and “endocrine resistance, confirming that the genes we identified are associated with resistance to therapy. Among these genes, CCND1, BCL2, ADCY1, ESR1, and PGR were also enriched in “endocrine resistance” and “estrogen receptor signaling” suggesting them as markers of resistance to both chemotherapy and endocrine therapy.

In network analysis, we detected that the protein products of 17 upregulated genes, namely ESR1, CCND1, SCUBE2, PGR, ERBB4, THBS1, GATA3, BCL2, TBC1D9, THSD4, STC2, CCDC170, STK32B, NBEA, PLAT, IL6ST, and NAT1 displayed connectivity with others. ESR1, PGR, BCL2, and SCUBE2 are being tested as a part of the 21-gene OncotypeDx Risk score test and 41-gene Biomark Assay (27). Despite that, the remaining 13 genes also have a high potential as markers of resistance to taxane-based neoadjuvant therapy in breast cancer. This should be addressed in future studies. Since we aimed to identify markers with a predictive power comparable to ER, in this study we focused on the top three genes with the highest degree and centrality in the network analysis. The most central gene was ESR1, the gene coding for ER, as expected. The other two central genes were CCND1 and SCUBE2.

ESR1 is a key marker for prognosis and responsiveness to therapy in breast cancer. ER-positive and ER-negative breast cancers display distinct gene expression profiles. ER-positive breast cancer is known to be more resistant to chemotherapy compared to ER-negative breast cancers. This knowledge of the impact of ESR1 on poor response to chemotherapy lead to the incorporation of endocrine therapy into the neoadjuvant setting in ER-positive early breast cancer, which improved treatment outcomes substantially (20).

CCND1 (cyclin D1) is a cell cycle protein that activates cyclin-dependent kinases (CDK4 and CDK6) and promotes the transition from G1 to S. CCND1 is overexpressed in almost half of breast cancer cases and is associated with poor prognosis, especially in ER-positive breast cancers (28). High CCND1 expression was also correlated with a poor response to anthracycline-taxane-based neoadjuvant chemotherapy, endocrine therapy, and immune checkpoint inhibitors (21). Since CCND1 is also a positive regulator of ER, both direct and ER-mediated actions of CCND1 can contribute to chemoresistance. In parallel to the growing evidence on the significance of CCND1 mediated activation of CDK4 and CDK6, several CDK inhibitors such as Ribociclib, and Abemaciclib are approved for breast cancer treatment (29).

Our study confirmed the association of ESR1 and CCND1 with resistance to taxane-based therapy in breast cancer. However, the most important contribution of our study is the identification of SCUBE2 as a relatively novel marker for taxane-based neoadjuvant therapy and a potential drug target. Our analysis suggested that SCUBE2 may be used as a predictive marker, especially in ER-positive/HER2- negative and HER2-positive breast cancers with a sensitivity and specificity similar to ESR1.

SCUBE2 encodes a secreted glycoprotein with epidermal growth factor-like repeats and a CUB domain (CUB: complement C1r/C1s, Uegf, Bmp1), which interacts with the cell surface (30). It is involved in the regulation of different molecules altered in cancer, such as sonic Hedgehog, and GRB7. SCUBE2 is used as a biomarker in various cancers, namely endometrium cancer, non-small cell lung carcinoma, colorectal cancer, glioma, and breast cancer. Genetic alterations in SCUBE2 were observed in uterine carcinomas, gastric cancer, melanoma, glioma, colorectal cancer, and many other cancers. SCUBE2 exhibited tumor suppressor function in glioma, non-small cell lung cancer and colorectal carcinoma. However, several reports suggest that SCUBE2 may display tumor suppressor or oncogenic effects in breast cancer in a context-dependent manner (31).

SCUBE2 was suggested to suppress the proliferation of breast cancer cells via inhibition of BMP, an inducer of cell proliferation in the MCF-7 metastatic breast cancer cell model and a mouse xenograft model (22). SCUBE2 induced an epithelial phenotype, suppressing epithelial-mesenchymal transition, invasion, and migration of MDA-MB-231 invasive ductal breast carcinoma cell line (23). Additionally, SCUBE2 positivity was associated with better disease-free survival in breast cancer patients with primary invasive ductal carcinoma (22). Despite these findings, a more recent study in MDA-MB-231 cells and in vivo models reported that increased expression of SCUBE2 in breast cancer stem cells induces epithelial-mesenchymal transition, increased tumor growth, and metastasis via activation of Notch signaling. Additionally, ectopic overexpression of SCUBE2 led to resistance to paclitaxel in TNBC cells (24).

Although there is a discrepancy about the exact molecular role of SCUBE2 in breast cancer, our analysis demonstrated that SCUBE2 may be a key marker for chemoresistance to taxane-based neoadjuvant therapy in breast cancer. Its predictive specificity and sensitivity were as high as ESR1, a well-established marker for chemoresistance in breast cancer. Ruey-Bing Yang’s group demonstrated that SCUBE2 acts as a coreceptor for vascular endothelial growth factor receptor 2 (VEGFR2), potentiating angiogenesis (32). The group demonstrated that knock out of SCUBE2 suppressed angiogenesis and tumor growth in melanoma and Lewis Lung carcinoma xenograft models, and an anti-SCUBE2 antibody displayed synergistic action with the anti-VEGF antibody in an orthotopic pancreas cancer model (33). Further research in breast cancer may reveal the molecular mechanisms by which SCUBE2 contributes to resistance to taxane-based chemotherapy and provide further insight into its molecular functions. Such insight may open the door for the development of novel molecular targeted agents against SCUBE2.

In our study, we analyzed breast cancer samples from all receptor subtypes as a pool to identify markers of resistance to taxane-based therapy, that can be utilized as a predictor and a new therapeutic target in large groups of patients. Despite that, a more detailed analysis of breast cancer subtype-specific cohorts, and a comparison of the markers for different subtypes would improve the efforts to predict responsiveness and personalize therapy in distinct breast cancer subtypes. This will be addressed in our future studies.

Another limitation of the study is the use of pCR as the surrogate of sensitivity to taxane-based neoadjuvant chemotherapy, since the number of patients in studies that utilized relapse free survival as the surrogate was much lower. pCR is mostly preferred in clinical trials to speed up the drug registry process. However, there are controversies about the efficiency of pCR as a surrogate of survival, and there seems to be differences in its surrogacy in different breast cancer subtypes (34, 35, 36). Therefore, the value of SCUBE2 to predict resistance to taxane-based therapy should also be validated in large cohorts using overall survival as the surrogate.

In conclusion, our study identified SCUBE2 as a novel marker for resistance to taxane-based therapy with a predictive power comparable to ESR1 and even better than CCND1 in breast cancer. Further investigations into the molecular functions of SCUBE2 in specific breast cancer subtypes may provide the opportunity to develop new, targeted therapies that can overcome resistance to taxane-based therapy in breast cancer.

Availability of Data and Materials

The GSE20194, GSE25066, and GSE32646 datasets we analyzed during the current study are available in Gene Expression Omnibus (GEO) repository (https://www.ncbi.nlm.nih.gov/geo/geo2r/).

Acknowledgments

The author gratefully acknowledges the use of the services and infrastructure of the Koc University Research Centre for Translational Medicine (KUTTAM).

Ethics Committee Approval: The author declares that ethics approval and consent to participate are not applicable, since the study is conducted using publicly available gene expression datasets.

Informed Consent: Not necessary.

Peer-review: Externally and internally peer-reviewed.

Conflict of Interest: No conflict of interest was declared by the author.

Financial Disclosure: The author declared that this study received no financial support.