Overexpression of EZH2/NSD2 Histone Methyltransferase Axis Predicts Poor Prognosis and Accelerates Tumor Progression in Triple-Negative Breast Cancer

Two histone methyltransferases, enhancer of zeste homolog 2 (EZH2) and nuclear SET domain-containing 2 (NSD2), are aberrantly expressed in several types of human cancers. However, the regulatory relationship between EZH2 and NSD2 and their prognostic values in breast cancer (BC) have not been fully elucidated. In this study, we demonstrated that EZH2 and NSD2 were overexpressed in BC compared with benign lesions and normal tissues using tissue microarray, immunohistochemistry, and bioinformatic databases. Both EZH2 and NSD2 expression were associated with pathological grade of tumor and lymph node metastasis. A comprehensive survival analysis using Kaplan-Meier Plotter database indicated that EZH2 expression was negatively correlated with relapse-free survival (RFS), overall survival (OS), distant metastasis-free survival (DMFS), and postprogression survival (PPS) in 3951 BC patients, and NSD2 expression was negatively correlated with RFS and DMFS. Notably, EZH2 and NSD2 expression were coordinately higher in triple-negative breast cancer (TNBC) than that in other subtypes. Stable knockdown of EZH2 using lentiviral shRNA vector significantly reduced the proliferation, migration and invasion abilities of TNBC cell line MDA-MB-231 and MDA-MB-468, and downregulated NSD2 expression as well as the levels of H3K27me3 and H3K36me2, two histone methylation markers catalyzed by EZH2 and NSD2, respectively. By contrast, overexpression of EZH2 using adenovirus vector displayed an inverse phenotype. Furthermore, knockdown of NSD2 in EZH2-overexpressing cells could dramatically attenuate EZH2-mediated oncogenic effects. Bioinformatic analysis further revealed the function and pathway enrichments of co-expressed genes and interactive genes of EZH2/NSD2 axis, suggesting that EZH2/NSD2 axis was associated with cell division, mitotic nuclear division and transition of mitotic cell cycle in TNBC. Taken together, EZH2/NSD2 axis may act as a predictive marker for poor prognosis and accelerate the progression of TNBC.


INTRODUCTION
Breast cancer (BC) is the most common malignancy among females and the incidence of BC has increased steadily in recent decades (1). About 12.2% of all newly diagnosed BC cases and 9.6% of all BC-related death cases worldwide are from China (1). The survival and prognosis of BC patients have been improved recently by surgery, radiotherapy, chemotherapy, endocrinotherapy, and molecular targeted therapy. However, metastasis and recurrence remain the clinical challenges and the leading causes of high mortality rates for BC patients, especially for those with triplenegative breast cancer (TNBC), since TNBC is insensitive to most available therapies (2). Currently, the molecular pathogenesis of TNBC is poorly understood and novel therapeutic targets for TNBC need to be explored. Hence, it is essential to identify biomarkers and therapeutic targets involved in unfavorable progression of TNBC.
Enhancer of zeste homolog 2 (EZH2) is a SET-domain containing histone methyltransferase and the catalytic subunit of the polycomb repressive complex 2 (PRC2) which inhibits gene transcription through trimethylation of Lysine 27 on histone H3 (H3K27me3). Previous studies have demonstrated the dual nature of EZH2, as it can function as a gene activator or repressor involving in cell proliferation, differentiation and cell cycle (3). Dysregulated EZH2 in tumorigenesis affects cell proliferation, apoptosis, epithelial-to-mesenchymal transition (EMT), invasion, immune response, radiotherapy, and chemotherapy resistance (3,4). EZH2 is overexpressed and contributes to aggressive progression and poor prognosis in diverse types of human malignancies, including breast, lung, prostate, gastric, colorectal cancer, melanoma, leukemia, and lymphoma (4)(5)(6)(7).
The enzymatic activity and potential therapeutic tractability of these histone methyltransferases attract broad attention in oncology, but the relationship between them is rarely understood. A previous study has shown that EZH2 may function upstream of NSD2 in prostate cancer (18). However, the relationship between EZH2 and NSD2 in BC remains unclear, and the clinicopathologic significances and prognostic values of them in BC have not been fully illuminated.
This study aimed to explore the roles and mechanisms of EZH2 and NSD2 in the pathogenesis of BC. The protein expression patterns of EZH2 and NSD2 in BC were investigated by immunohistochemistry (IHC) on tissue microarray (TMA) series and analyzing Human Protein Atlas database. The transcriptional levels of EZH2 and NSD2 in BC were evaluated using Oncomine and UALCAN databases. The correlations between their expression and clinicopathological characteristics were analyzed. Kaplan-Meier Plotter database was used to assess the prognostic values of EZH2 and NSD2 in BC. The interaction of EZH2 and NSD2 was illustrated in vitro in the proliferation, migration, and invasion of TNBC cells. Furthermore, bioinformatics analysis was conducted to envisage the possible functions and related pathways of EZH2/NSD2 axis in TNBC.

Tissue Specimen Collection
Paraffin-embedded surgical tissue specimens of 146 BC cases were collected from the Department of Pathology at First Affiliated Hospital of Dali University (Dali, China), including 126 cases of invasive ductal carcinoma (IDC), 12 cases of ductal carcinoma in situ (DCIS), and 8 cases of invasive lobular carcinoma (ILC). None of the patients received any preoperative adjuvant therapy. Two pathologists confirmed the pathological diagnoses of all the specimens. For controls, 24 cases of benign breast lesions were included.

Construction of TMA
The construction of TMA followed the guidelines reported previously (19). The hematoxylin and eosin (HE) stained sections of source paraffin-embedded tissue blocks were used to identify the representative areas of the specimens by the pathologists. A TMA recipient block was constructed by a TMA instrument. A tiny core with 5 millimeters in diameter was extracted from each source tissue block and then inserted into the recipient block. Sections cut from the TMA block were used for HE staining and immunohistochemical analysis.

HE Staining
For histological analysis under the microscope, the TMA blocks were sectioned into 4 mm and stained with HE according to a routine protocol (20).

IHC and Protein Expression Analysis
IHC was performed on the TMA sections to assess the protein levels of EZH2 and NSD2 in BCs compared with benign breast lesions. Deparaffinization and rehydration, antigen retrieval, and endogenous peroxidase repression were performed on 4 mmthick TMA slices. Then the slices were incubated with anti-EZH2 (1:500, Cell Signaling Technology, #5246) or anti-NSD2 (1:400, Abcam, #ab75359) antibody at 4°C overnight and then labeled using an EnVision HRP Kit (Dako) at room temperature for 30 min. The signal was visualized with diaminobenzidinechromogen, followed by counterstaining with hematoxylin. Human Protein Atlas database (https://www.proteinatlas.org/) was further used to analyze the differential protein expression of EZH2 and NSD2 in BCs and normal breast tissues. The expression was defined using the criteria set in the database. Protein expression score is based on immunohistochemical data manually scored with regard to staining intensity (negative, weak, moderate, or strong) and fraction of stained cells (<25, 25-75, or >75%). Each combination of intensity and fractions is automatically converted into a protein expression level score as follows: negative-not detected; weak <25%-not detected; weak combined with either 25-75% or >75%-low; moderate <25%low; moderate combined with either 25-75% or >75%-medium; strong <25%-medium, strong combined with either 25-75% or >75%-high.

Gene mRNA Expression Analysis
Oncomine database (http://www.oncomine.org) was used to analyze the differential mRNA expression of EZH2 and NSD2 in BCs and normal breast tissues. UALCAN database (http:// ualcan.path.uab.edu/) was used to analyze the correlation between EZH2 and NSD2 and the differential expression of them among different molecular subtypes of BC.

Survival Analysis
The data of genes expression and survival information of 6,234 BC patients contained in Kaplan-Meier Plotter database were downloaded from the Gene Expression Omnibus (GEO), TCGA, and European Genome-Phenome Archive. A total of 35 datasets containing 3,951 BC patients were utilized to assess the associations of EZH2 (probe ID, 203358_s_at) or NSD2 (probe ID, 209053_s_at) mRNA expression with overall survival (OS), relapse-free survival (RFS), distant metastasis-free survival (DMFS), and postprogression survival (PPS). Hazard ratio (and 95% confidence intervals) and logrank P-values were calculated and displayed. The cases were categorized into highor low-expression group using the median expression as the cutoff value.

Cell Culture
Human TNBC cell line MDA-MB-231 was purchased from the cell bank of Institute of Zoology of the Chinese Academy of Sciences (Kunming, China) and MDA-MB-468 was purchased from Cobioer Biosciences Company (Nanjing, China). Cells were cultured in DMEM medium (Gibco) supplemented with 10% FBS (Gibco), 100 U/ml penicillin (Beyotime) and 0.1 mg/ml streptomycin (Beyotime) at 37°C with 5% CO 2 .

Adenovirus Overexpression Vector and Transient Transfection
The EZH2 adenovirus expression vector was constructed by inserting EZH2 cDNA fragment into an adenovirus shuttle vector (Adv/CMV/IRES), named Adv-EZH2. The packing and purification of the recombinant adenovirus vectors were performed by GenePharma Company (Shanghai, China). Cells were infected with recombinant adenovirus vector at a multiplicity of infection of 10.

CCK8 Assay
The proliferation ability of cells was measured using a CCK8 Assay Kit (Beyotime) according to the manufacturer's protocols. Briefly, 4×10 3 stably transfected or transient transfected cells were seeded into 96-well plates and CCK8 assay was performed at different time points of 0, 24, 48, 72, and 96 h. Ten microliters of CCK8 reagent was added into each well and incubated for 1 h. Subsequently, the absorbance at 450 nm wavelength was measured on a microplate spectrophotometer (Gene Company).

Wound Healing Assay
Wound healing assay was performed to observe the migration ability of cells. 7x10 5 cells were plated in six−well plates and a 10 ml pipette tip was used to scratch a vertical wound. The detached cells were removed and the adherent cells were incubated with serum-free medium. Images of the scratched fields (100×) were captured at 0 and 24 h using a light microscope (Olympus BX53), and the area of the wound was measured with ImageJ software. At least three scratched fields were recorded in each group.

Transwell Assay
Transwell chambers (Corning) with 8 µm pores were used to evaluate the invasion ability of cells. 2×10 5 cells in serum-free medium were seeded into the upper chamber with Matrigel (BD Biosciences, #356234). Medium containing 10% FBS was added to the lower chambers. After 8 h, the cells remaining on the upper membrane were removed, while the cells invading through the membrane were fixed with paraformaldehyde, stained with crystal violet (Beyotime), and observed under a light microscope (200×). 33% acetic acid solution was used to elute the staining and the absorbance at 570 nm wavelength was measured on a microplate spectrophotometer (Gene Company).

Protein-Protein Interaction Network Construction and Functional Enrichment Analysis
The co-expressed genes of EZH2 and NSD2 with a Pearson correlation (≥0.5) in BC were analyzed using UALCAN database. STRING online system (https://string-db.org/) was used to construct a protein-protein interaction (PPI) interactions network of EZH2, NSD2 and their co-expressed genes. Subsequently, these genes were subjected to gene ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analysis in Database for Annotation, Visualization and Integrated Discovery (DAVID) (https://david.ncifcrf.gov/ home.jsp). P < 0.05 was set as the cutoff value.

Statistical Analysis
Statistical analysis was performed using SPSS 21.0 software and datas were expressed as the mean ± standard deviation (SD). Differences between two groups were assessed by Student's t test. The chi-square test was used to determine differences of EZH2 or NSD2 expression among groups with different clinicopathological characteristics. The spearman rank correlation test was used to analyze the expression correlation between EZH2 and NSD2. Kaplan-Meier analysis followed by the logrank test was used to analyze the association between EZH2 or NSD2 expression and the survival rate. P < 0.05 indicated a statistically significant difference.

Construction of TMA Containing BCs and Benign Lesions
Firstly, the TMA blocks were constructed using paraffinembedded surgical tissue specimens composed of 24 cases of benign breast lesions and 146 cases of BC, including 126 cases of IDC, 12 cases of DCIS, and 8 cases of ILC. HE staining showed the structural integrity and histological features of core lesion regions in TMA sections. The microscopic appearance of IDC was heterogeneous, ranging from tumors with low-grade nuclei and well-developed tubule formation to tumors consisting of sheets of anaplastic cells ( Figure 1A). IDC cells invaded into the stroma with a surrounding desmoplastic response ( Figure 1A). DCIS cells were confined to the duct and the nuclear appearance ranged from bland or monotonous to pleomorphic, whereas the basement membrane all remains intact ( Figure 1A). ILC cells usually aligned in chains or strands and invaded into stroma individually ( Figure 1A).

EZH2 and NSD2 Are Coordinately Overexpressed in BC Tissues
To investigate the protein expression patterns of EZH2 and NSD2 in BC, IHC assay was performed in TMA sections. Both EZH2 and NSD2 showed a nuclear staining pattern ( Figure 1B). The protein expression levels of EZH2 and NSD2 in benign lesions were significantly lower than that in cancers ( Figure 1B and Table 1). EZH2 expression displayed strong correlation with NSD2 expression in all BC samples (Spearman correlation coefficient r s = 0.694, P < 0.001, Table 2). In different histological types of samples, EZH2 expression was strongly correlated with NSD2 expression in IDC (Spearman correlation coefficient r s = 0.745, P < 0.001) and DCIS (Spearman correlation coefficient r s = 0.683, P < 0.05) ( Table 2). However, there was no correlation between their expression in ILC probably due to the limited samples (Spearman correlation coefficient r s = −0.143, P > 0.05, Table 2). Besides, analysis using Human Protein Atlas database also revealed that EZH2 and NSD2 were predominantly located in the cell nucleus, and the medium/high expression rates of EZH2 and NSD2 in BC were 83.33% (10/12 cases) and 100% (10/10 cases), respectively ( Figure 1C).
To examine the mRNA expression of EZH2 and NSD2 in breast cancer, Oncomine database was used to analyze their mRNA expression. According to four different studies, the mRNA expression level of EZH2 was upregulated in BC by 2.241 to 8.874 folds relative to that in normal breast tissues (Figures 2A-D). Meanwhile, the mRNA expression level of NSD2 was increased in BC by 2.121 to 4.797 folds ( Figures 2E-H). A positive correlation between EZH2 and NSD2 mRNA expression in BC was also demonstrated using UALCAN database (Spearman correlation coefficient r s = 0.45, P < 0.05, Figure 2I).

High Expression of EZH2 and NSD2 Are Associated with Poor Prognosis of BC
To evaluate the potential roles of EZH2 and NSD2 in BC, the relationship between EZH2 or NSD2 protein expression and the clinicopathological characteristics of BC patients was analyzed. As presented in Table 3, both EZH2 and NSD2 protein expression correlated with pathological grade of tumor and lymph node metastasis; the high grade and positive lymph  node involvement were associated with increased rate of positive EZH2 and NSD2 (P < 0.001, P = 0.003). On the contrary, EZH2 and NSD2 protein expression were not associated with other clinicopathological parameters, including age, size of tumor and histological type of tumor. Kaplan-Meier Plotter database was used to investigate the relationship between EZH2 or NSD2 and survival rate of BC patients. 3,951 patients were classified into high-or low-expression group according to the median expression of EZH2 or NSD2. As shown in Figure 3, high EZH2 expression level was significantly correlated with lower RFS (P = 4.17e−13, Figure 3A), OS (P = 0.014, Figure 3B), DMFS (P = 0.0021, Figure 3C) and PPS (P = 0.00097, Figure  3D) rates. High NSD2 expression was associated with RFS (P = 4.6e−9, Figure 4A) and DMFS (P = 0.003, Figure 4C), but had no effect on OS (P = 0.65, Figure 4B) and PPS (P = 0.68, Figure  4D) in BC patients. By evaluating risk factors associated with patients using Cox proportional hazard models, it was determined that EZH2 and NSD2 were significant predictors of outcomes in BC patients (Figures 3 and 4).

EZH2 and NSD2 Expression Are Higher in TNBC Than That in Other Molecular Subtypes of BC
Interestingly, we observed that the positive rates of EZH2 and NSD2 protein expression in TNBC were higher than that in other molecular subtypes of BC, including luminal A, luminal B and HER2-enriched BC (P = 0.033, P = 0.026, Table 4). UALCAN database was used to explore the transcriptional levels of EZH2 and NSD2 in different molecular subtypes of BC. EZH2 expression in TNBC was significantly higher than that in luminal BC (P = 8.88e−16), while there was no statistical difference in EZH2 expression between TNBC and HER2enriched BC (P = 0.103, Figure 5A). Compared with luminal and HER2-enriched BC, NSD2 expression was higher in TNBC (P = 8.83e−3, P = 8.04e−5, Figure 5B).

EZH2 Upregulates NSD2 Expression and Histone Methylation and Promotes the Proliferation, Migration, and Invasion of TNBC Cells
As EZH2 and NSD2 expression were significantly upregulated in TNBC, a TNBC cell line MDA-MB-231 was used to examine the functional relationship between EZH2 and NSD2. We constructed a MDA-MB-231 cell line in which EZH2 expression was stably knocked down by lentivirus shRNA vector. The results showed that knockdown of EZH2 decreased NSD2 mRNA and protein expression remarkably in MDA-MB-231 cells ( Figures 6A, B). In addition, knockdown of EZH2 also reduced the levels of H3K27me3 and H3K36me2, which are histone methylation markers catalyzed by EZH2 and NSD2, respectively ( Figure 6B). By contrast, adenovirus-mediated EZH2 overexpression significantly increased NSD2 expression as well as the methylation levels of H3K27 and H3K36 in MDA-MB-231 cells ( Figures 6A, B). These results suggested that NSD2 functions as a downstream gene of EZH2. To assess the biological role of EZH2 in TNBC, CCK8 assay revealed that cell proliferation was significantly impaired in shEZH2 stably transfected cells ( Figure 6C). Next, we evaluated cancer cell migration and invasion abilities through wound healing assay and transwell assay. Decreased expression of EZH2 also impeded the migration and invasion of cells ( Figures 6D, E). On the contrary, overexpression of EZH2 in MDA-MB-231 cells displayed an inverse phenotype ( Figures 6C-E). Similar effects were observed in MDA-MB-468 cells with knockdown or overexpression of EZH2 ( Supplementary Figures 1A-E).

EZH2-Mediated Oncogenic Effects Require NSD2 Expression
Knockdown of NSD2 could inhibit the proliferation, migration, and invasion in MDA-MB-231 and MDA-MB-468 cells ( Supplementary Figures 2A-D). To determine the role of NSD2 in EZH2-mediated oncogenic effects, overexpression of EZH2 was performed in MDA-MB-231 cells with stable knockdown of NSD2 and control cells. As expected, overexpression of EZH2 upregulated H3K27me3 level in cells, but EZH2-induced NSD2 expression and H3K36me2 level were abolished by NSD2 shRNA treatment ( Figures 7A, B). Cells overexpressing EZH2 displayed high proliferation, migration,

Function and Pathway Enrichment Analysis of Co-Expressed Genes of EZH2/ NSD2 Axis
To analyze the functions and pathways that EZH2/NSD2 axis may be involved in, 94 co-expressed genes of EZH2 and NSD2 with a Pearson correlation (≥0.5) in BC were obtained using UALCAN database (Table S1). A PPI network comprising 96 nodes and 2,087 edges with PPI-enrichment P-value < 1.0e−16 was constructed using STRING online system to explore the interactions of EZH2, NSD2, and their co-expressed genes ( Figure 8A). GO enrichment analysis revealed that these genes were significantly enriched in the biological processes (BP) about cell proliferation, such as cell division, mitotic nuclear division, chromosome segregation and condensation, microtubule-based movement, cytokinesis, transition of mitotic cell cycle, and DNA repair ( Figure 8B and Table S2). These genes were widely distributed in various kinds of cellular component (CC) and mostly located in nucleus, nucleoplasm and cytosol ( Figure 8C and Table S3). In terms of molecular function (MF), these genes were responsible for ATP, microtubule, DNA, and protein binding, thus maintaining the activity of ATPase, microtubule motor and protein kinase ( Figure 8D and Table S4). KEGG pathway analysis showed that these genes were mainly enriched in the pathways involved in cell cycle, meiosis, DNA replication, and repair ( Figure 8E and Table S5).

Function and Pathway Enrichment Analysis of Interactive Genes of EZH2/NSD2 Axis
Based on the PPI analysis ( Figure 8A), 45 genes interacted with EZH2 and seven genes interacted with NSD2. Six genes   interacting with both EZH2 and NSD2 were further identified, that were cyclin A2 (CCNA2), cyclin dependent kinase 2 (CDK2), lysine demethylase 2B (KDM2B), kinesin family member 11 (KIF11), kinesin family member 23 (KIF23), and proliferating cell nuclear antigen (PCNA) ( Figure 9A). The eight genes in Figure 9A were subjected to GO and KEGG pathway analysis. As shown in Table 5, numerous biological processes were associated with the development of cancer, including mitotic nuclear division, cell division, cellular response to multiple stimuli, Ras protein signal transduction, and transition of mitotic cell cycle. The significant enriched GO terms in CC were nucleoplasm, cyclin A2-CDK2 complex, centrosome, nucleus, kinesin complex, and spindle ( Table 5). MF analysis indicated that these genes may function in RNA polymerase II core promoter sequence-specific DNA binding and microtubule motor activity ( Table 5). The enriched KEGG pathway of these genes is cell cycle pathway ( Table 5). The expression of CCNA2, CDK2, KDM2B, KIF11, KIF23, and PCNA in 116 cases of TNBC and 114 cases of normal breast tissue were evaluated using UALCAN database. Similar with the expression patterns of EZH2 and NSD2, the gene expression heatmap showed that these six genes were also significantly upregulated in TNBC ( Figure 9B). Compared with normal breast tissues, the median expression levels of CCNA2, CDK2, KDM2B, KIF11, KIF23, and PCNA in TNBC were increased by 14.5, 1.44, 1.55, 9.29, 12.47, and 2.87 folds, respectively (P<0.001, Figure 9C). Additionally, the expression levels of these genes were further verified by qRT-PCR in MDA-MB-231 and MDA-MB-468 cells with stable knockdown of EZH2 or NSD2. Knockdown of EZH2 decreased the expression of CCNA2,  CDK2, KIF11, KIF23, and PCNA, but had no effect on the expression of KDM2B ( Figure 9D). The expression levels of CCNA2, CDK2, and PCNA were also downregulated in cells with stable knockdown of NSD2, while no significant changes of KDM2B, KIF11and KIF23 expression were observed ( Figure 9E).

DISCUSSION
Generally, breast carcinogenesis arises from altered gene expression or gene mutations, which eventually leads to the dysregulation of numerous oncogenes, tumor suppressor gene, and non-coding RNAs. Epigenetic alterations are heritable changes which affect gene expression profiles but do not change the primary DNA sequence (23). Investigations about aberrant epigenetic factors including histone modifications and DNA methylation, mainly focus on the molecular mechanisms in cancer initiation and development, the new biomarkers for cancer progression and the potential of epigenetic therapy for BC (23). Histone methyltransferases are potential therapeutic targets for malignancies, but it is poorly understood that how they are linked together. Asangani (18) (18). Studies have highlighted the multiple roles of EZH2 in the pathogenesis of BC. The increased expression of EZH2 in histologically normal breast epithelium is associated with higher risk of developing cancer, suggesting that EZH2 may be a potential marker for identifying preneoplastic lesions of the breast as well as a possible target for preventative intervention (24). Previous studies showed that EZH2 is frequently expressed in 75.7% of inflammatory BCs and associated with unfavorable prognostic factors including larger size, high histological grade, increased risk for distant metastasis, and first degree family history of BC (25,26). Therefore, EZH2 may function as a useful biomarker of long-term metastatic risk and worse prognosis in women with BC, and warrant further validation studies. As a novel oncogene identified in recent years, NSD2 overexpression is associated with the clinical stage and TNM classification in BC (27). Gene expression profiling analysis revealed that NSD2 expression is highly elevated in endocrinotherapy-resistant BC tissues and cell lines (28). Consistent with previous reports, this study demonstrated that EZH2 and NSD2 proteins were expressed in 82.2% (120/146) and 79.5% (116/146) cases of BC, respectively. Gene mRNA expression analysis also showed that the transcriptional levels of EZH2 and NSD2 were upregulated in BC. Moreover, we found a high correlation between the expression of these two histone methyltransferases at both protein and mRNA levels. Both EZH2 and NSD2 expression strongly correlated with histological grade of tumor and lymph node metastasis. Although the effect of EZH2 or NSD2 on survival in BC patients has been previously described, this study is the first to offer a comprehensive analysis about the prognostic values of EZH2 and NSD2 for RFS, OS, DMFS, and PPS in a large group of BC patients. Increased NSD2 expression level was correlated with shorter RFS and DMFS, and EZH2 could be used for prognosis assessment in RFS, OS, DMFS and PPS. Notably, this study highlighted that EZH2 and NSD2 were highly expressed in TNBC when compared with other molecular subtypes of BC. These results implied that EZH2 and NSD2 are involved in the development and progression of BC, especially for TNBC.
Ectopic overexpression of EZH2 induces malignant transformation of the mammary gland cells by promoting cell invasion and anchorage-independent growth in vitro (29). EZH2 overexpression disrupts ductal morphogenesis and causes epithelial hyperplasia in vivo (30). In this study, knockdown of EZH2 significantly inhibited TNBC cell proliferation and impaired cell migration and invasion, whereas overexpression of EZH2 produced an inverse phenotype. Silencing endogenous expression of NSD2 significantly suppresses the proliferation and A B FIGURE 5 | The relative expression of EZH2 messenger RNA (mRNA) (A) and NSD2 mRNA (B) among different molecular subtypes of breast cancer (BC) were analyzed using UALCAN database. *P < 0.05, **P < 0.01, ***P < 0.001. metastasis of BC cells through inhibiting the Wnt/b-catenin signaling pathway (27). NSD2 induces resistance to tamoxifen in BC by methylating H3K36me2 at the promoters of key glucose metabolic enzyme genes and heightening the pentose phosphate pathway (28). NSD2 also mediates TNBC cell resistance to the EGFR inhibitor gefitinib via stimulating ADAM9-EGFR-AKT signaling (31). This is the first report that knockdown of NSD2 abolishes EZH2-mediated TNBC cell proliferation, migration and invasion, indicating that the oncogenic function of EZH2 depends on NSD2 expression.
A previous report suggested that the regulation of NSD2 expression by EZH2 may occur at the posttranscriptional level through microRNAs network (18). The authors found that a significant induction of NSD2 protein in normal epithelial cells when EZH2 is overexpressed, but NSD2 transcription level remains unchanged (18). In prostate cancer cells, Asangani (18) found that EZH2 catalyzes H3K27 trimethylation to inhibit miR-26a, miR-31, and miR-203, which target the 3'-UTR of NSD2. However, we testified that EZH2 functioned upstream of NSD2 and upregulated both the mRNA and protein expression of NSD2, as well as H3K27me3 and H3K36me2 levels in TNBC cells. It was reported that EZH2 has a dual nature, because it can act as a gene repressor or activator (3). On one hand, EZH2 can restrain gene expression by methylating both histone and non-histone proteins on promoter regions of several genes, such as CDKN1C, FOXC1, E-cadherin, RAD51, and RUNX3 (3,32). On the other hand, it can activate gene expression and pathways by interacting with transcription factors and co-factors, such as androgen receptor, b-catenin, STAT3, estrogen receptor a, RelA/RelB, and PCNA-associated factor (3). Based on the results of the present study, we speculated that EZH2 may enhance the transcription of NSD2 though interacting with some transcription factors or co-factors in TNBC. Thus, further researches will be needed to elucidate the direct or indirect regulation of NSD2 expression by EZH2. In 94 co-expressed genes with EZH2 and NSD2, six interactive genes with EZH2/NSD2 axis were further identified, that are CCNA2, CDK2, KDM2B, KIF11, KIF23, and PCNA. GO analysis of these genes revealed that the function of EZH2/NSD2 axis was associated with cell division, mitotic nuclear division and transition of mitotic cell cycle. These results suggested that EZH2/NSD2 axis may participate in the tumorigenesis and development of TNBC by impacting the process of cell mitosis. KEGG pathway analysis indicated that EZH2/NSD2 axis was connected to the cell cycle pathway. CCNA2-CDK2 complex acts as a core component of cell cycle and often aberrantly expressed in cancer (33). PCNA plays a crucial role in DNA replication during cell cycle progression (34). KDM2B promotes TNBC cell proliferation by repressing the transcription of the cell cycle inhibitors p15INK4B, p16INK4A, and p57KIP2 (35). KIF11 is an essential molecular motor protein for chromosome segregation and mitotic spindle formation. Knockdown of KIF11 causes G2/ M phase arrest and inhibits self-renewal in chemotherapy resistant TNBC stem cells (36). KIF23 is responsible for cytoplasm separation and axon elongation, and its high level is significantly correlated with unfavorable overall survival in BC (37). These genes, as essential regulators of the cell division cycle, are usually taken as markers of cell proliferation. Furthermore, recent data suggest that they also function in cell migration and invasion. CCNA2 participates in the modulation of RhoA activity and its dysregulation by autophagy in the late phase of mitosis appears to be associated with EMT, affecting the aggressiveness of the tumor (38). CCNA2 also facilitates EMT via the integrin avb3 signaling, representing a crucial regulator of tumor cells metastasis (39). Inhibition of CDK2 reduces invasion of prostate cancer cells and reintroduction of CDK2 rescues the invasion ability, indicating that CDK2 is a crucial factor in metastasis of cancer (40). Phosphorylation of SIRT2 at S331 site by CDK2 impairs the catalytic activity of SIRT2, influences microtubule dynamics and antagonizes SIRT2-mediated inhibition of cell motility, including adhesion and migration (41). Positive expression rate of PCNA mRNA in colorectal cancer with liver metastasis was higher than those without liver metastasis, so it may be useful for evaluating liver metastasis of cancer cells (42). KDM2B regulates actin cytoskeleton and epithelial-to-mesenchymal transition by upregulating Rho-GTPases and activating FAK/PI3K signaling, thus promoting tumor cell motility (43)(44)(45). KIF11 is a key downstream molecule that responds to directional cues in chemotaxis to govern the direction of migration (46). KIF11 knockdown inhibits the migration and invasion of BC cells via decreasing the phosphorylation of ERK, AMPK, AKT, and CREB (47). KIF23 promotes gastric cancer cell proliferation, migration, and invasion by activating Wnt/b-catenin signaling pathway through direct interaction with Amer1 (48). The expression levels of CCNA2, CDK2, and PCNA were reduced in cells with  knockdown of EZH2 or NSD2, while the expression levels of KIF11and KIF23 were only reduced in cells with knockdown of EZH2, suggesting that CCNA2, CDK2, and PCNA may be regulated by both EZH2 and NSD2. The proposed mechanism of cell proliferation, migration and invasion regulated by EZH2/ NSD2 axis was summarized in Figure 10. Our study provided a foundation for investigating the interaction between these six genes and EZH2/NSD2 axis and the mechanism of TNBC growth and metastasis in follow-up work.
Given the important role of EZH2 in tumor pathogenesis, developing small molecule inhibitors of EZH2 becomes popular and clinical trials are underway to test these molecules. Nevertheless, there are still some challenges in direct targeting of EZH2 in cancers due to its various roles. For example, EZH2 plays a major part in normal hematopoietic stem cell differentiation and self-renewal, B cell differentiation, thymopoiesis, and lymphopoiesis (49, 50). Additionally, inactivating mutations or deletion of EZH2 have been identified in several types of tumor, such as Kras-driven lung adenocarcinoma and T cell precursor acute lymphoblastic leukemia, suggesting the tumor-suppressive effect of EZH2 in certain contexts (51,52). Thus, identifying alternative downstream targets of EZH2 such as NSD2, may provide a new therapeutic strategy.
In summary, the expression of EZH2 and NSD2 are coordinately increased in BC, particularly in TNBC. These two histone methyltransferases act as prognostic indicators for BC patients. EZH2 promotes the proliferation, migration and invasion abilities of TNBC cells via upregulating NSD2 expression. EZH2/ NSD2 axis may contribute to the progression of TNBC by affecting the cell cycle pathway. Our findings unveiled the role of EZH2/ NSD2 histone methyltransferases axis in pathogenesis and progression of TNBC, which may facilitate the development of diagnosis and therapy against TNBC.

DATA AVAILABILITY STATEMENT
The original contributions presented in the study are included in the article/Supplementary Material. Further inquiries can be directed to the corresponding author.

ETHICS STATEMENT
The studies involving human participants were reviewed and approved by the Research Ethics Committee of First Affiliated Hospital of Dali University. The patients/participants provided their written informed consent to participate in this study.

AUTHOR CONTRIBUTIONS
BG and YP designed the study. BG, XL, ZL, and LZ acquired interpreted the data. BG and YP drafted and revised the manuscript. All authors contributed to the article and approved the submitted version.