A Six-Epithelial–Mesenchymal Transition Gene Signature May Predict Metastasis of Triple-Negative Breast Cancer

Introduction

Triple-negative breast cancer (TNBC) is characterized immunohistochemically by the absence of estrogen receptor (ER), progesterone receptor (PR) and human epidermal growth factor receptor-2 (HER2), though HER2 deficiency is also defined as the lack of genomic amplification of the gene encoding HER2 as measured using fluorescence in situ hybridization (FISH).¹ It constitutes approximately 15–20% of all breast cancers with a higher incidence in young individuals and is commonly of higher grade.² TNBC is also more aggressive and has greater potential to metastasize to distant sites. Patients with TNBC often suffer from worse clinical outcomes as evidenced by the higher relapse and low survival rates.^3,4 Moreover, conventional prognostic markers, such as patient age, tumour size and lymph node involvement are seldomly associated with distant metastasis in TNBC patients.⁵

Despite the advance in the treatment of other types of breast cancers using endocrine therapy and targeted therapy,^6,7 non-surgical management of TNBC remains largely limited to chemotherapy, alone or in combination with poly (ADP-ribose) polymerase (PARP) inhibitors in patients with germline breast cancer gene (BRCA) mutations.⁸ In particular, neoadjuvant chemotherapy (NACT) is often applied to TNBC patients prior to surgery with an aim to reduce local tumour masses and thus increase the operability and achieve a better cosmetic outcome.⁹ A pathological complete response (pCR), which is defined as the absence of residual invasive lesions in both breasts and axilla after NACT,¹⁰ is associated with favourable disease-free survival (DFS) and overall survival (OS),^11,12 and is thus a prognostic marker of TNBC patients.¹¹ Moreover, pCR is commonly used as a surrogate endpoint in clinical trials for developing new NACT drugs or approaches in the treatment of TNBC.¹³ Nevertheless, fewer than half of TNBC patients achieve pCR.¹⁴ Of note, a proportion TNBC patients with residual diseases do not relapse and can achieve long-term survival.¹⁵ This calls for further development of tools to predict clinical outcomes of TNBC patients who do not achieve pCR to NACT.

The diverse responses and outcomes in TNBC patients treated with NACT are closely associated with the biologically heterogeneous nature of the disease.^2,16 Indeed, according to gene expression profile analyses, TNBC has been dissected into either four or seven genetically distinct subtypes.^17,18 For example, the seven-subtype classification independently predicted pCR to NACT in TNBC patients. However, it was not associated with distant metastasis-free survival (DMFS) and OS.¹⁹ Similarly, the PAM50 subtype predictor that characterizes intrinsic subtypes of breast cancer using a 50-gene assay predicted pCR to NACT across all breast cancer subtypes, but none of the PAM50 signatures at diagnosis is associated with pCR when the analysis was confined to TNBC.²⁰ Moreover, a number of studies have demonstrated the potential of epigenetic signatures and tumour infiltrating immunophenotypes as predictive and prognostic tools in TNBC patients treated with NACT.^21,22

In this study, we have compared gene expression profiles in paired TNBC tissues from patients pre- and post-NACT. We report here that a group of six-epithelial–mesenchymal transition (EMT) genes were expressed at higher levels in TNBC tumours before and after NACT in patients who developed distant metastasis than those who did not, and that these EMT genes were further upregulated in post-NACT residue tumours compared with prior-NACT samples. Interrogating Gene Expression Omnibus (GEO) TNBC datasets validated the 6-EMT gene signature to be predictive of metastasis after NACT and moreover demonstrated that this signature may also serve as a predictive tool of metastasis in TNBC patients treated with surgery followed by adjuvant chemotherapy.

Patients and Methods

Patients and Samples

The discovery cohort consisted of ten retrospectively recruited TNBC patients who received NACT and displayed locally invasive residual disease in the breast and/or regional lymph nodes. All patients received the EC-T regimen (4 cycles of doxorubicin combined with cyclophosphamide followed by 4 cycles of docetaxel or paclitaxel according to the National Comprehensive Cancer Network Guideline²³) before surgery at the Department of Breast Surgery of Shanxi Bethune Hospital (Taiyuan, Shanxi, China) between 2014 and 2018. Paired formalin-fixed paraffin-embedded (FFPE) pre-NACT core needle biopsy tissues and post-NACT samples were retrieved from the Department of Pathology of Shanxi Bethune Hospital. All patients gave their written informed consent with the study independently approved by the ethics committee/institutional review board of the Shanxi Bethune Hospital. Validation cohorts used were the MD Anderson Cancer Center (MDACC)-based (Texas, Houston, USA) study that included 178 TNBC patients who received NACT before surgery (GEO dataset GSE25066) and an ICO-UMGC-based (Saint Herblain, France) cohort that encompassed 107 TNBC patients treated with adjuvant chemotherapy following surgery (GSE103091).

Profiling of RNA Expression

Total RNA of FFPE TNBC tissues were subjected to the Affymetrix Human Transcriptome Array to profile the expression of mRNAs after verification of the diagnosis using H&E staining and immunohistochemistry (IHC) with or without FISH tests. Five 10μm FFPE sections were collected for each sample and total RNA extracted using the RecoverAllTM Total Nucleic Acid Isolation Kit (Thermo Scientific, AM1975) with purification using the QIAGEN RNeasy^® Mini Kit (Qiagen, 74,104). RNA samples were quantified by NanoDrop ND-2000 (Thermo Scientific) and all considered adequate for microarray analysis (the optical density 260/280 >1.8; total RNA yield >100ng). Samples were then transcribed to double stranded cDNAs, labeled with biotin, and hybridized onto the microarray (GeneChip^® Hybridization, Wash, and Stain Kit, Affymetrix, 900,720) followed by scanning with the Affymetrix Scanner 3000 (Affymetrix).

Differential Expression and Gene Signature Analysis

Raw microarray data were extracted using the Affymetrix GeneChip Command Console software (v4.0, Affymetrix) and normalized with RMA Method by Expression Console software (v1.3.1, Affymetrix). Differentially expressed transcripts were identified using the Limma package in R platform (v3.6.2)²⁴ comparing patients with and without development of metastasis after NACT. Fold changes were determined for each case pre- and post-NACT and differences with a |logFC| (log2 fold change)≥1 and a p-value <0.05 were considered significant. The differentially expressed genes were illustrated in volcano plots and heatmaps using ggplot2 and pheatmap packages in R platform, respectively. Gene set enrichment analysis (GSEA) was performed to identify the biological pathways and gene signatures associated with metastasis after NACT by the ClusterProfiler package in R software,²⁵ with the Hallmark gene set, downloaded from MSigDB (Molecular Signatures Database; https://www.gsea-msigdb.org/gsea/msigdb/index.jsp). Gene sets with a normalized enrichment score (NES) ≥1 or ≤-1 and a false discovery rate (FDR) <0.05 after 1000-time permutations were considered enriched significantly. The distribution of common genes was depicted in Venn diagrams (http://bioinformatics.psb.ugent.be/webtools/Venn). The correlation among the selected genes was assessed and visualized by the Hmisc package with Pearson analysis in the R platform.

Statistical Analyses

Prognostic analyses were conducted using R software and survival curves generated with GraphPad Prism 8 software. Mann–Whitney U-test was used to determine the difference between groups. The Cox proportional hazard model was employed to assess the prognostic value of the identified signature. The weighted linear prognostic model combined with gene expression in validation cohorts were used to calculate the prognostic scores of the selected gene signature, and Uno’s concordance index (C-index)²⁶ was used to quantify the accuracy of predicted models in different validated datasets. A patient was classified into the high-risk (with a prognostic score of the signature higher than the cut-off value) or the low-risk (with a prognostic score lower than the cut-off value) group. A cut-off value calculated by Youden index from the time-dependent receiver operating characteristic (ROC) curve.²⁷ The log-rank method (Mantel-Haenszel test) was employed to compare the DMFS between the high- and low-risk group, which was depicted with Kaplan–Meier estimators of the survival function.²⁸ A two-tailed P value <0.05 with 95% CI (Confidence intervals) was considered significant.

Unless otherwise stated, statistical analyses were performed by SPSS software, version 23.0 (IBM Corp., Armonk, NY).

Results

Patient and Disease Characteristics

The demographic information of TNBC patients included in the discovery cohort and the clinicopathological characteristics of their tumours are summarized in Table 1. All patients were initially diagnosed through H&E staining of core needle biopsy tissues in conjunction with immunohistochemistry (IHC) to detect ER, PR and HER2. FISH analysis was conducted to further clarify HER2 mRNA positivity if IHC staining of the HER2 protein in a case was scored ++ according to the guidelines of the American Society of Clinical Oncology (ASCO).²⁹ All patients underwent NACT with doxorubicin/cyclophosphamide followed by docetaxel or paclitaxel (EC-T) before surgery. The patients were reassessed using ultrasonography and mammography after completion of NACT. According to the response evaluation criteria in solid tumours (RECIST),³⁰ five patients reached partial response (PR) and 4 exhibited stable disease (SD), whereas the other 1 developed progressive disease (PD). The median follow-up time of these patients was 25 months, ranging from 14 to 67 months. One patient developed multiple metastatic lesions in the lung, brain, and bone in the 14th month during follow-up, and another was detected with multiple metastases in the lung and bone in the 17th month. The remainders were distant metastasis-free at the latest follow-up.

Table 1 Characteristics of TNBC Patients and Tumors

Post-NACT Metastasis of TNBC Is Associated with High Expression of EMT Genes

We carried out microarray-based transcriptomic profiling of paired TNBC tissues before and after NACT (Figure 1). To identify genes potentially involved in metastasis after NACT followed by surgery, we divided patients into metastasis and non-metastasis groups based on whether metastatic lesions were detected at the latest follow-up. Comparative analyses identified 699 significantly upregulated genes and 129 downregulated genes in tumours from the metastasis group at diagnosis (logFC ≥ 1 or =< −1, P<0.05). Similarly, 687 genes were found to be expressed at higher levels whereas 62 genes at lower levels in tumours of the metastasis group after NACT (logFC ≥ 1 or =< −1, P<0.05) (Figure 1A and B; Supplementary Table 1). Hierarchical clustering analysis showed that the tumours in the metastasis group displayed strong similarity in their transcriptomic profiles, whereas the tumours in the non-metastasis groups exhibited wide variations in the composition of their transcriptomes (Figure 1C and D). Consistent with the critical role of EMT in metastasis,^31,32 gene set enrichment analysis (GSEA) based on hallmark gene sets revealed that the EMT pathway was highly augmented among the 10 most significantly enriched pathways in metastasis tumours before and after NACT (Figure 2 and Supplementary Table 2).

Figure 1 Comparison of transcriptomic profiles of TNBCs between the metastasis and non-metastasis group. (A and B) Volcano plots of transcripts of TNBCs of metastasis group relative to non-metastasis group at diagnosis (A) or post-NACT (B). Red and blue dots represent up-regulated and down-regulated genes with |log2(fold change)|≥1.0 and p values <0.05, respectively. (C and D) Heat maps with clusters of differentially expressed genes with |log2(fold change)|≥1.0 and p values <0.05 at diagnosis (C) and post-NACT (D). Abbreviations: M, metastasis; N, non-metastasis; TNBC, triple-negative breast cancer; NACT, neoadjuvant chemotherapy.

Figure 2 GSEA plots depicting the enrichment of signal pathways in tumors of the metastasis relative to non-metastasis group. (A and B) Ridgeline plots of the 10 most significantly enriched pathways based on the hallmark gene sets at diagnosis (A) and post-NACT (B). (C) The EMT pathway is positively enriched in tumors of the metastasis group relative to non-metastasis group at diagnosis (left) and post-NACT (right). Abbreviations: M, metastasis; N, non-metastasis; NES, normalized enrichment score; NACT, neoadjuvant chemotherapy; EMT, epithelial–mesenchymal transition.

Differential Responses of EMT Genes to EC-T NACT

Having found that the EMT pathway is positively enriched in TNBC tumours that developed metastatic lesions, we examined the potential effect of EC-T NACT on the expression of EMT genes. Strikingly, the EMT pathway responded differently to EC-T NACT in the two metastasis cases. In metastasis case 1, the EMT pathway was further positively enriched post-NACT whereas in metastasis case 2 this pathway was negatively enriched in the post-NACT tumour compared with the pre-NACT tumour (Figure 3A). There were subsets of EMT genes in cases of the metastasis group remained unchanged or were downregulated after NACT (Figure 4A; Supplementary Table 3). Nevertheless, there were 15 EMT genes that were commonly increased in the two metastasis cases, including 2 genes that were upregulated in both cases before NACT, 6 genes that were upregulated in both cases after NACT and 7 genes, before and after NACT (Figure 4A; Supplementary Table 3). Among these 15 EMT genes were 7 genes including lumican (LUM), secreted frizzled-related protein 4 (SFRP4), collagen type VI alpha 3 chain (COL6A3), matrix metalloproteinase-2 (MMP2), C-X-C Motif Chemokine Ligand 12 (CXCL12), decorin (DCN), and high-temperature requirement A serine peptidase 1 (HTRA1) that were commonly expressed at higher levels in each of the tumour in the metastasis group before and after NACT (Figure 4A), suggestive of the potential importance of these genes in metastasis of TNBC after EC-T NACT. As depicted in the correlation matrices (Figure 4B and C), although the expression levels of these 7 EMT genes were largely correlated with one another before NACT (Figure 4B), the expression of COL6A3 was disassociated with the others after NACT (Figure 4C). This finding implied that the other 6 genes (LUM, SFRP4, MMP2, CXCL12, DCN, and HTRA1) were most likely associated with TNBC metastasis after NACT. Similar to the differential responses of EMT genes in the metastasis cases, the EMT genes the non-metastasis group were either increased or decreased after NACT in a patient-dependent manner (Figure 3B).

Figure 3 GESA plots based on the hallmark gene sets depicting the differential enrichment of the EMT pathway caused by NACT in the individual patients from the metastasis group (A) and the non-metastasis group (B). Abbreviations: B, biopsy samples before the NACT treatment; S, surgery specimens after finishing the NACT treatment; NES, normalized enrichment score; EMT, epithelial–mesenchymal transition; NACT, neoadjuvant chemotherapy.

Figure 4 (A) The Venn diagram displaying the overlaps among the upregulated EMT genes in pre-NACT (purple) and post-NACT (pink) NACT tumors in the metastasis group. The EMT genes that were upregulated in each of the post-NACT compared with pre-NACT tumors in the metastasis group were also depicted (green and yellow). (B and C) Pearson correlation matrices describing the relationships among the indicated 7 genes in paired tumors pre- (B) and post-NACT (C). The Pearson ρ coefficients were exhibited in various size of circles as well as the shade of colors. Numbers represent the exact values of coefficients. *P<0.05; **P <0.01; ***P <0.001. Abbreviations: M, metastasis; EMT, epithelial–mesenchymal transition; NACT, neoadjuvant chemotherapy.

The Six-EMT Gene Signature Predicts Post-NACT Metastasis of TNBC

To test the strength of the 6-EMT gene signature in predicting metastasis of TNBC after NACT, we analysed transcriptomic data at diagnosis from the MD Anderson Cancer Center (MDACC)-based cohort of 178 TNBC patients who received sequential taxane and anthracycline neoadjuvant chemotherapy, inclusive of 113 patients who did not achieve pCR (GSE25066).³³ Among these 113 patients, 55 developed distant metastases, whereas the remaining 58 patients were metastasis-free at the last follow-up (Table 2).

Table 2 Patients’ Characteristics in the MDACC Validation Sets

A time-dependent ROC curve with an AUC (area under curve) at 0.63 (p=0.014) demonstrated the ability of the 6-EMT signature to predict metastasis at 18 months after initial diagnosis in the 113 non-pCR patients (Figure 5A). When the Cox proportional hazard model was applied, an optimal cut-off value of 0.1261 was derived through the Youden index for stratifying patients with and without development of metastasis after NACT (sensitivity: 47.27%; specificity: 79.31%).²⁷ By use of this cut-off value, the 113 patients were divided into high- and low-risk groups. A patient with a predicted score higher than or equal to 0.1261 was assigned to the high-risk group, whereas a patient with a predicted score lower than 0.1261, the low-risk group. The distribution of the prognostic score and the status of metastasis of each patient are shown in Figure 5B. Notably, patients with higher prognostic scores (n=38) appeared prone to metastasis (p=0.0008, HR=0.36, 95% CI: 0.200–0.658) (Figure 5C). Together these results suggest that the 6-EMT gene signature at diagnosis may be useful for prediction of metastasis of TNBC patients who do not achieve pCR to NACT.

Figure 5 Validation of the 6-EMT gene signature in the MDACC-based cohort of the GEO dataset (GSE25066). (A) A time-dependent ROC curve for 18 months among the TNBC patients with residual disease after NACT. The cut-off value was generated with the highest specificity and sensitivity (Youden Index). (B) The distribution of the prognostic values of the TNBC patients with residual disease (top panel) and the corresponding distant relapse-free survival for each patient (bottom panel). Patients were divided into high-risk (right part) and low-risk (left part) groups by a line in the diagram, according to the cut-off value (0.126). (C) Kaplan–Meier analysis of the probability of metastasis-free survival of the patients with residual disease post-NACT in the MADACC cohort. Abbreviations: EMT, epithelial–mesenchymal transition; ROC, receiver operating characteristic curve; AUC, area under curve; RD, residual disease; HR, hazard ratio; TNBC, triple-negative breast cancer; NACT, neoadjuvant chemotherapy.

The 6-EMT Gene Signature Predicts Metastasis in TNBC Patients Without NACT

To further validate the significance of the 6-EMT gene signature in predicting TNBC metastasis, we applied the Cox proportional hazard model using the cut-off value of 0.1261 to the ICO-UMGC-based cohort comprised 107 TNBC patients treated with adjuvant chemotherapy following surgery, including 31 patients who developed distant metastases during follow-up (GSE103091, Table 3).³⁴ Similar to the observations made with the MDACC-based cohort, the 6-EMT signature classifier divided the 107 patients into high- and low-risk groups based on the optimal cut-off value (0.1261; Figure 6A) with a C-index at 0.66, and the patients with higher prognostic scores (n=48) appeared prone to metastasis (p=0.032, HR=0.46, 95% CI: 0.225–0.937) (Figure 6B). When time-dependent ROC curves created at 1.5, 3, and 5 years were used to evaluate the power of the 6-gene signature to predict metastasis, the AUC increased progressively from 0.52 at 1.5 years to 0.70 at 5 years (Figure 6C). These results suggest that the 6-EMT gene signature may also be useful as a predictive marker in TNBC patients who undergo surgery followed by adjuvant chemotherapy.

Table 3 Patients’ Characteristics in the Extended Validation Sets

Figure 6 The 6-EMT gene signature may predict metastasis of TNBC in patients treated with surgery followed by adjuvant chemotherapy as shown by analysis of the ICO-UMGC-based cohort of the GEO dataset (GSE103091). (A) The distribution of the prognostic values of the 107 TNBC patients who were classified into high-risk and low-risk groups (top panel), and the corresponding distant metastasis-free survival status for each patient (bottom panel). (B) Kaplan–Meier analysis of the probability of metastasis-free survival of the patients in the ICO-UMGC-based cohort. (C) Time-dependent ROC curves for patients in the ICO-UMGC-based cohort. AUCs at 1.5, 3, and 5 years were used to evaluate the accuracy of the prognostic model. Abbreviations: EMT, epithelial–mesenchymal transition; TNBC, triple-negative breast cancer; ROC, receiver operating characteristic; AUC, area under curve; HR, hazard ratio; NACT, neoadjuvant chemotherapy.

Discussion

Curative treatment of TNBC remains an unmet medical need worldwide and solutions to this problem are hampered by the disease’s inherent heterogeneity.^2,16,35 Moreover, genomic evolution resulting from selection pressure afforded by chemotherapy also contributes to the often dismal outcome of patients.³⁵ Although considerable efforts have been made to identify biomarkers and/or molecular targets in order to advance treatment of TNBC,¹⁴ longitudinal information concerning the dynamic changes of molecular determinants of diverse TNBC biological behaviours remains scant. By use of paired pre- and post-NACT TNBC tissues, we have identified a 6-EMT gene signature that was not only predictive of metastasis of TNBC after NACT, but also potentially a useful tool for identification of TNBC patients with high metastatic risk who are treated with surgery followed by adjuvant chemotherapy.

As EMT is an enabling feature of cancer cell invasion and metastasis,^31,32 it is not surprising that the EMT pathway was positively enriched in TNBCs of the metastasis group compared with those of non-metastasis group before and after NACT. However, it was intriguing that the EMT pathway was either further enhanced or attenuated in TNBCs in both the metastasis and non-metastasis groups, suggesting that the pathway responded to the chemotherapeutic drugs in a highly tumour-dependent manner, consistent with diverse responses of TNBCs to chemotherapy commonly observed in clinical practice.³⁶ Nonetheless, we identified a subset of interrelated EMT genes that were increased in post-NACT relative to pre-NACT tumours. Moreover, these genes were expressed at higher levels in the tumours of the metastasis group compared with those in the non-metastasis group before and after NACT, pointing to the potential importance of these genes in regulating metastasis. Indeed, interrogation of a MDACC-based TNBC cohort that received NACT revealed that the 6-EMT gene signature was predictive of DMFS in those who did not achieve pCR. Unlike those patients in our discovery cohort who received the same NACT regimen, the MDACC patients underwent varying NACT protocols, suggesting that the predictive value of this signature is independent of chemotherapeutic drugs applied.

Further application of 6-EMT gene signature to an ICO-UMGC-based cohort suggested it may also be useful as a predictive marker of metastasis in TNBC patients who do not receive NACT. Although many genetic and epigenetic features have been demonstrated to be potential predictive tools in TNBC,^15,37 a predictive biomarker that can be applied broadly to TNBC regardless of the treatment approaches is uncommon. Nevertheless, given the small size of our discovery cohort and the non-coherent nature of the validation cohorts, it is apparent that the usefulness the 6-EMT gene signature needs to be further confirmed in large-scale prospective cohorts of patients. Similarly, whether this signature identifies tumours that overlap with known molecular subtypes of TNBC needs to be clarified. For example, the recently defined claudin-low subtype of TNBC express high levels of EMT genes and is associated with low pCR rates after NACT.³⁸ Should the 6-EMT gene signature withstand further large-scale in-depth validation, it may serve as a useful tool for identifying TNCB patients with a high risk of distant metastasis after NACT who may benefit more from other therapeutic modalities.

Of the proteins encoded by the 6 EMT genes, MMP2 is a member of type IV collagenases that degrades extracellular matrix,³⁹ whereas CXCL12 is a chemokine (C-X-C motif) ligand⁴⁰ and decorin (DCN), a member of the extracellular matrix small leucine-rich proteoglycan family in the tumor microenvironment that influences the biological behaviour of various types of tumours.⁴¹ While the roles of these proteins in breast cancer metastasis have been well established,^42–44 the role of HTRA1, a member of the trypsin family of serine proteases, in regulating metastasis of breast cancer remains unclear. Intriguingly. HTRA1 loss was associated with aggressive behaviour and poor survival in breast cancer suggesting tumour suppressive functions.⁴⁵ In contrast, the SFRP family member SFRP4 promotes cell growth and resistance to chemotherapy in prostate cancer,⁴⁶ providing a plausible explanation for its inclusion in the 6-EMT gene signature identified for breast cancer. Finally, LUM that encodes an extracellular matrix protein has been shown to inhibit or even reverse several metastatic features that EMT endows in breast cancer.⁴⁷ Regardless, our results suggest that these six proteins may all contribute to breast cancer metastasis. Nevertheless, further studies are required to clarify how they contribute, particularly HTRA1, SFRP4, and LUM, and whether they functionally cooperate in breast cancer metastasis. Similarly, a caveat of this study is that the transcriptomic data analysed involved bulk tumours which are not homogenous. It is thus possible that data could reflect gene expression by non-tumour cells such as stromal cells.⁴⁸ Notwithstanding the potential use of the 6-EMT gene signature in diagnostic specimens, single-cell sequencing analysis could be instrumental in delineating the function of these genes in breast cancer metastasis, including their role in response to NACT.

Conclusion

The high expression of a signature of six EMT genes (LUM, SFRP4, MMP2, CXCL12, DCN and HTRA1) at diagnosis may be of predictive value of metastasis in TNBC patients who did not achieve pCR to NACT and in those who are treated with surgery in combination with adjuvant therapy.