Qingshui Wang1,2,3†
Wenting Zhong2†Lin Deng4†Qili Lin5†Youyu Lin2Hongxia Liu2Luyun Xu2Lingfang Lu2Yajuan Chen2Jianping Huang5Meichen Jiang5Han Xiao5Jie Zhang6,7,8He Li4Yuxiang Lin6,7,8*
Chuangui Song6,7,8*
Yao Lin1,2*- 1Central Laboratory at the Second Affiliated Hospital of Fujian Traditional Chinese Medical University, Collaborative Innovation Center for Rehabilitation Technology, Fujian University of Traditional Chinese Medicine, Fuzhou, China
- 2Key Laboratory of Optoelectronic Science and Technology for Medicine of Ministry of Education, College of Life Sciences, Fujian Normal University, Fuzhou, China
- 3Fujian Provincial Key Laboratory of Hepatic Drug Research, Fuzhou, China
- 4Department of General Surgery, The 900th Hospital of the Joint Logistics Support Force, Fuzhou, China
- 5Department of Pathology, Fujian Medical University Union Hospital, Fuzhou, China
- 6Department of Breast Surgery, Fujian Medical University Union Hospital, Fuzhou, China
- 7Department of General Surgery, Fujian Medical University Union Hospital, Fuzhou, China
- 8Breast Cancer Institute, Fujian Medical University, Fuzhou, China
Background: Triple-negative breast cancer (TNBC) is the most invasive and metastatic subtype of breast cancer. SUMO1-activating enzyme subunit 1 (SAE1), an E1-activating enzyme, is indispensable for protein SUMOylation. SAE1 has been found to be a relevant biomarker for progression and prognosis in several tumor types. However, the role of SAE1 in TNBC remains to be elucidated.
Methods: In the research, the mRNA expression of SAE1 was analyzed via the cancer genome atlas (TCGA) and gene expression omnibus (GEO) database. Cistrome DB Toolkit was used to predict which transcription factors (TFs) are most likely to increase SAE1 expression in TNBC. The correlation between the expression of SAE1 and the methylation of SAE1 or quantity of tumor-infiltrating immune cells was further invested. Single-cell analysis, using CancerSEA, was performed to query which functional states are associated with SAE1 in different cancers in breast cancer at the single-cell level. Next, weighted gene coexpression network (WGCNA) was applied to reveal the highly correlated genes and coexpression networks of SAE1 in TNBC patients, and a prognostic model containing SAE1 and correlated genes was constructed. Finally, we also examined SAE1 protein expression of 207 TNBC tissues using immunohistochemical (IHC) staining.
Results: The mRNA and protein expression of SAE1 were increased in TNBC tissues compared with adjacent normal tissues, and the protein expression of SAE1 was significantly associated with overall survival (OS) and disease-free survival (DFS). Correlation analyses revealed that SAE1 expression was positively correlated with forkhead box M1 (FOXM1) TFs and negatively correlated with SAE1 methylation site (cg14042711) level. WGCNA indicated that the genes coexpressed with SAE1 belonged to the green module containing 1,176 genes. Through pathway enrichment analysis of the module, 1,176 genes were found enriched in cell cycle and DNA repair. Single-cell analysis indicated that SAE1 and its coexpression genes were associated with cell cycle, DNA damage, DNA repair, and cell proliferation. Using the LASSO COX regression, a prognostic model including SAE1 and polo-like kinase 1 (PLK1) was built to accurately predict the likelihood of DFS in TNBC patients.
Conclusion: In conclusion, we comprehensively analyzed the mRNA and protein expression, prognosis, and interaction genes of SAE1 in TNBC and constructed a prognostic model including SAE1 and PLK1. These results might be important for better understanding of the role of SAE1 in TNBC. In addition, DNA methyltransferase and TFs inhibitor treatments targeting SAE1 might improve the survival of TNBC patients.
Introduction
Triple-negative breast cancer (TNBC) is a specific type of breast cancer which lacks the expression of estrogen receptor (ER), progesterone receptor (PR), and human epidermal growth factor receptor (HER2) (Jhan and Andrechek, 2017; Vagia et al., 2020). Although TNBC only comprises approximately 15–20% of all breast cancers, it has a higher recurrence and mortality rate than other types of breast cancer (Sporikova et al., 2018; da Silva et al., 2020). Patients with TNBC generally develop distant metastasis within the first 3 years after initial treatment, with the mortality rate could reach about 40% in the first 5 years (Wilson et al., 2019; Chang-Qing et al., 2020). Generally, TNBC do not respond to hormone or HER2-targeted therapies due to the absence of ER, PR, and HER2. Therefore, a combination of surgery, chemotherapy, and radiotherapy remains the most common treatment strategy for TNBC patients. While significant advances have been made in the treatment of TNBC, the prognosis of these patients remains extremely poor. Therefore, it is still necessary to identify new predictive biomarkers for treatment response and promise therapeutic targets for TNBC.
SUMOylation is a key post-translational modification with the form of conjugation of SUMO molecule to an acceptor lysine of target proteins, which is catalyzed by E1-activating enzyme, E2-conjugating enzyme, and several E3-ligating enzymes (Hickey et al., 2012; Pichler et al., 2017; Wang Q. et al., 2019; Chang and Yeh, 2020). SUMO1-activating enzyme subunit 1 (SAE1), an E1-activating enzyme, is indispensable for protein SUMOylation. SAE1 has a molecular weight of 110 kDa, predominantly localizes to the nucleus and can be found in a variety of tissues. SAE1 forms a thioester bond between SUMO protein and UBA2/SAE2 by regulating the activation of ATP-dependent SUMO protein, and then participates in the process of SUMO protein modification (Amente et al., 2012; Yang et al., 2019; Wang et al., 2020; Ong et al., 2021). SAE1 mainly mediates acetylation and phosphorylation and plays an important role in chromosome division and apoptosis. It has been found that SUMOylation is an important mechanism in cellular responses to stress and one that is usually abnormal in many cancers (Song et al., 2014; Yang et al., 2019). For instance, SUMO E3 protein inhibitor of activated STAT3 (PIAS3) is upregulated in various cancer types, like prostate, lung, colon, and breast cancer and brain tumor (Abbas et al., 2015; Li et al., 2015; He et al., 2018; Polimeno et al., 2020; Tseng et al., 2020). Recently, research and function of SAE1 have been increasingly analyzed in organism. SAE1 has been found to play a significant role in various tumorigenesis and development, such as glioma, colon cancer, and non-small cell lung cancer (Yang et al., 2019; Zhang et al., 2019). Nevertheless, its functions of SAE1 in TNBC are indistinct by far. Whether SAE1 can be used as a therapeutic target for TNBC and its biological roles remain to be further studied.
Materials and Methods
Extraction of SUMO1-Activating Enzyme Subunit 1 Expression and Clinical Data From Triple-Negative Breast Cancer Datasets
The mRNA expression of SAE1 and clinical data for TNBC patients were extracted from the cancer genome atlas (TCGA)1 and gene expression omnibus (GEO) database2 including GSE21653 (Sabatier et al., 2011b), GSE31448 (Sabatier et al., 2011a), GSE45827 (Gruosso et al., 2016), GSE53752, GSE65216 (Maire et al., 2013), GSE36295 (Merdad et al., 2014), GSE37751 (Terunuma et al., 2014), GSE38959 (Komatsu et al., 2013), and GSE61724 (Mathe et al., 2015). The method for extracting microarray gene expression values is based on our previous research (Mathe et al., 2015; Wang Q. S. et al., 2019). First, the probe ID was converted into a gene symbol. When a gene was mapped to different probes, the genic expression value was calculated by the average expression value. Next, the median normalization was performed using the robust multichip averaging method. Kaplan-Meier survival curve was conducted to estimate the overall survival (OS) and disease-free survival (DFS) by using TCGA and GSE21653, respectively.
Patients and Specimens
Two hundred seven patients from Fujian Medical University Union Hospital between August 2013 and October 2017 were retrospectively enrolled in this study. All 207 patients were histologically confirmed TNBC, with a median age of 55 years (range, 27–72 years). Relevant clinicopathological information, including tumor size, nodal status, tumor grade, lymphovascular invasion (LVI), type of surgery, chemotherapy, and radiotherapy status, were obtained from patients’ medical records. OS was defined as the time from the date of diagnosis until death from any cause. DFS was defined as the time from the date of diagnosis to the date of clinical relapse (with histopathology confirmation or radiological evidence of tumor recurrence). The last follow-up was December 31, 2020. All patients met the following inclusion criteria: (1) no history of other malignant tumors or bilateral breast cancer; (2) received post-operative adjuvant chemotherapy with at least three cycles of anthracycline-based and three cycles of taxane-based regimens; (3) primary tumor size was pT1c-pT2 (1 cm < T ≤ 5 cm); and (4) demographic, clinicopathological, and follow-up information were complete. This study was approved by the Research Ethics Committee of Fujian Medical University Union Hospital and informed consent was obtained from each participant.
Immunohistochemistry Staining Analysis
We conducted immunohistochemistry (IHC) staining analysis to measure the protein expression of SAE1 in TNBC tissues and adjacent normal breast tissues according to the standard immunoperoxidase staining procedure. Slides were incubated with anti-SAE1 (ab97523, Abcam, Cambridge, United Kingdom, diluted 1:300) according to the manufacturer’s instructions. A negative control was also prepared via substitution of a primary antibody with phosphate-buffered saline (PBS). The IHC staining scores of SAE1 were evaluated by two independent pathologists. The percentage of stained positive cells was calculated from 1 to 4: 1, 0–25%; 2, 26–50%; 3, 51–75%; and 4, 75–100%. The staining intensity score was ranged from 0 to 3: 0, no staining; 1, weak staining; 2, moderate staining; and 3, strong staining. The percentage of positive tumor cells and the staining intensity were multiplied to produce a weighted score for each case. A score of 8–12 was defined as high expression level and a score of 0–7 was defined as low expression.
Methylation Analysis
The methylation site and beta values of SAE1 promoter in TNBC was obtained from UCSC Xena Browser3. The beta values of CpG islands had been acquired using the Illumina HumanMethylation450K platform and had been preprocessed by TCGA by using standard protocols. One of the most widely used techniques to measure DNA methylation is the Illumina Infinium HumanMethylation450 BeadChip array, which covers approximately 450,000 CpG sites. At each CpG site, methylation is quantified by the beta value beta.
Weighted Gene Correlation Network Analysis
Weighted gene coexpression network (WGCNA) is a common algorithm for building gene coexpression networks and performed through WGCNA R package (Langfelder and Horvath, 2008). The WGCNA hypotheses that the coexpression gene network follows the scale less distribution firstly, defines the adjoining function of the gene coexpression correlation matrix and gene network formation, calculates the difference coefficient of different nodes, and constructs a hierarchical cluster tree accordingly finally.
Gene Ontology and Kyoto Encyclopedia of Genes and Genomes Pathway Enrichment Analysis
Gene Ontology (GO) analysis and kyoto encyclopedia of genes and genomes (KEGG) pathway enrichment were calculated by functional enrichment tool DAVID4 (Huang da et al., 2009). DAVID bioinformatics resources provide an integrated biological database and a repository of analytic tools for systematic exploration of biological meaning of gene set DAVID. The default parameters in the tool were used, and enriched pathways were ranked according to their enrichment scores. A P-value of <0.05 was identified as enriched functions.
Single-Cell Analysis
CancerSEA5 (Yuan et al., 2019) depicts single-cell functional status maps that contain 14 functional states (including angiogenesis, apoptosis, invasion, EMT, differentiation, proliferation, DNA damage, metastasis, hypoxia, inflammation, cell cycle, DNA repair, stemness, and quiescence) from 25 types of cancer including 41,900 individual cells. CancerSEA was used to explore the potential roles of SAE1 and its coexpression genes in breast cancer.
Transcription Factors Identification
Cistrome DB Toolkit database6 (Zheng et al., 2019) is a website tool that allows users to query transcription factors (TFs) that might regulate genes. We used the Cistrome DB Toolkit to predict which TFs are most likely to increase SAE1 expression in breast cancer. The correlation coefficient R > 0.3 is considered relevant.
Gene Set Variation Analysis
Gene set variation analysis (GSVA) provides increased power to detect subtle pathway activity changes over a sample population in comparison with corresponding methods. In the study, the pathway activity of protein SUMOylation in TNBC were analyzed using R package “GSVA.”
Statistical Analysis
In the study, t-test was applied to calculate the statistical correlation. Correlations between SAE1 expression and clinicopathological characteristics were performed by the chi-squared test. DFS and OS were calculated by the Kaplan-Meier method, and survival curves were compared via log-rank test. We conducted Cox regression analysis for the univariate and multivariate survival analyses. All P-values < 0.05 were considered statistically significant.
Results
The Level of SUMO1-Activating Enzyme Subunit 1 mRNA Was Highly Expressed in Triple-Negative Breast Cancer and Predicted Poor Prognosis
To examine the mRNA expression pattern of SAE1 in TNBC, publicly available datasets from TCGA and nine related GEO databases, including GSE21653, GSE31448, GSE45827, GSE53752, GSE65216, GSE36295, GSE37751, GSE38959, and GSE61724, were applied (Figures 1A–J). The mRNA level of SAE1 was found to be significantly increased in TNBC tissues (TCGA, p < 0.001; GSE21653, p < 0.001; GSE31448, p < 0.001; GSE45827, p < 0.001; GSE53752, p < 0.01; GSE65216, p < 0.001; GSE36295, p < 0.01; GSE37751, p < 0.001; GSE38959, p < 0.001; and GSE61724, p < 0.001). Moreover, meta-analysis containing 847 tissues from 10 TNBC databases mentioned above further demonstrated that the level of SAE1 mRNA was increased in TNBC (p < 0.001; Figure 1K). Next, we found that the average mRNA expression of SAE1 was higher in TNBC patients compared with luminal A or luminal B breast cancer subtypes based on GEO database GSE31448, GSE45827, and GSE65216 (Supplementary Figure 1). However, no significant differences were detected between TNBC- and HER2-enriched subtypes. These results could indicate that the increased expression of SAE1 might be more associated with ER-negative and PR-negative status.
Figure 1. The mRNA expression of SAE1 in TNBC. (A) The mRNA expression of SAE1 for TNBC patients in TCGA database. (B–I) The mRNA expression of SAE1 for TNBC patients in GSE21653 (B), GSE31448 (C), GSE45827 (D), GSE53752 (E), GSE65216 (F), GSE36295 (G), GSE37751 (H), GSE38959 (I), and GSE6172 (J). Meta-analysis of SAE1 expression in TNBC patients (K). **p < 0.01; ***p < 0.001.
In order to explore whether the mRNA expression of SAE1 could predict poor survival for TNBC patients, Kaplan-Meier survival analysis was performed for the OS and DFS. TNBC patients with high SAE1 expression had a significantly detrimental prognosis compared with those with low SAE1 expression (OS, p > 0.05; DFS, p < 0.001) (Figures 2A,B). Meanwhile, we performed GSVA to conduct GO analysis and assign protein SUMOlyation pathway activity estimates to individual samples and found that the activity of SUMOylation pathway was higher in TNBC patients compared with all other breast cancer subtypes (Supplementary Figure 2). The result of survival analysis indicated that high activity of SUMOylation pathway had a significantly detrimental DFS prognosis compared with those with low activity of SUMOylation pathway (Supplementary Figure 3).
Figure 2. The OS and DFS of SAE1 for TNBC patients. The prognostic value of SAE1 mRNA expression for OS (A) and DFS (B) in TNBC patients by Kaplan-Meier analysis. The log-rank test was performed to evaluate survival difference with the best cutoff value.
SUMO1-Activating Enzyme Subunit 1 mRNA Expression Is Positively Correlated With Forkhead Box M1 and Negatively Correlated With Methylation Site (cg14042711)
To identify members of a putative molecular network that might affect SAE1 expression, we examined TFs that might regulate SAE1 gene transcription using the Cistrome DB Toolkit. We checked the TNBC cell line MDA-MB-231 and found that jun proto-oncogene (JUN), tumor protein p53 (TP53), methyl-CpG-binding domain protein 3 (MBD3), E2F transcription factor 1 (E2F1), v-myc avian myelocytomatosis viral oncogene homolog (MYC), and forkhead box M1 (FOXM1) possessed the regulatory potential in MDA-MB-231 cell (Figure 3A). Subsequent analyses revealed that SAE1 mRNA and FOXM1 mRNA expression were positively correlated in TNBC tissues using TCGA database (Figure 3B, r = 0.34, p < 0.001). However, the expression of SAE1 mRNA was not significantly correlated with JUN, TP53, MBD3, E2F1, and MYC (Figures 3C–G). The transcription factor FOXM1 coordinates the expression of cell cycle-related genes and plays a pivotal role in tumorigenesis and cancer progression. To understand the role of FOXM1 in breast cancer, gene expression data (GSE55204) from breast cancer cells before and after FOXM1 knockdown to highlight specific FOXM1 transcriptional networks were performed by Bergamaschi et al. (2014). Based on the GSE55204 database, we found that when FOXM1 was knocked down (Supplementary Figure 4A), SAE1 expression decreased (Supplementary Figure 4B). Meanwhile, the expression of FOXM1 and SAE1 was positively correlated (Supplementary Figure 4C). These experimental results further support our conclusion.
Figure 3. Relationships of SAE1 with TFs and DNA methylation. (A) TFs with high regulatory potential in MDA-MB-231 cell lines (10k distance to TSS). (B) Correlation between SAE1 and FOXM1 mRNA expression. (C) Correlation between SAE1 and E2F1 mRNA expression. (D) Correlation between SAE1 and TP53 mRNA expression. (E) Correlation between SAE1 and JUN mRNA expression. (F) Correlation between SAE1 and MBD3 mRNA expression. (G) Correlation between SAE1 and MYC mRNA expression. (H) Relationships of SAE1 with DNA methylation in TNBC. *p < 0.05; **p < 0.01; ***p < 0.001.
The methylation of promoter is one of the critical mechanisms regulating gene expression which has been shown to control transcription in mammals. Previous studies have revealed that increased expressions of many genes were associated with hypomethylation of promoter. The association between the SAE1 expression and methylation of SAE1 for TNBC patient was detected via UCSC Xena. The associations between methylation of 11 SAE1 CpG sites and SAE1 expression are displayed in Figure 3H. The Pearson’s correlation coefficient between methylation site (cg14042711) and SAE1 is −0.32 (p < 0.05). However, we did not detect any significant correlation between 10 other SAE1 CpG sites and SAE1 expression. These results indicate that TFs and DNA methylation modifications might play a considerable role in TNBC processes by regulating SAE1 expression.
SUMO1-Activating Enzyme Subunit 1 Expression Was Not Correlated With Immune Factors
It is well-known that tumor-infiltrating immune cells could play an important role during cancer development and are of prognostic value in TNBC patients. Therefore, we performed a correlation analysis of SAE1 expression with the quantity of 10 tumor-infiltrating immune cells (T cells, CD8 T cells, cytotoxic lymphocytes, B lineage, NK cells, monocytic lineage, myeloid dendritic cells, neutrophils, endothelial cells, and fibroblasts). However, it is indicated that the expression of SAE1 was not correlated with different immune cell infiltrates (Figures 4A–J).
Figure 4. Relationships of SAE1 with tumor-infiltrating immune cells. (A) Correlation between SAE1 expression and fibroblasts. (B) Correlation between SAE1 expression and T cells. (C) Correlation between SAE1 expression and CD8 T cells. (D) Correlation between SAE1 expression and cytotoxic lymphocytes. (E) Correlation between SAE1 expression and B lineage. (F) Correlation between SAE1 expression and NK cells. (G) Correlation between SAE1 expression and monocytic lineage. (H) Correlation between SAE1 expression and myeloid dendritic cells. (I) Correlation between SAE1 expression and neutrophils. (J) Correlation between SAE1 expression and endothelial cells.
The Genes Coexpressed With SAE1 Enriched in the Cell Cycle, DNA Repair, DNA Damage, and Proliferation
To better understand the relevance and underlying mechanisms of SAE1 expression in breast cancer, we investigated the functional state of SAE1 in single-cell level via the CancerSEA database. Figure 5A displays the correlation between SAE1 and 14 functional states (including angiogenesis, apoptosis, invasion, EMT, differentiation, proliferation, DNA damage, metastasis, hypoxia, inflammation, cell cycle, DNA repair, stemness, and quiescence). Results indicated that SAE1 was positively correlated with cell cycle, DNA repair, DNA damage, and proliferation (Figures 5B–F).
Figure 5. The function of SAE1 in breast cancer cells. (A) Single-cell analysis displayed the correlation between SAE1 and 14 functional states (including angiogenesis, apoptosis, invasion, EMT, differentiation, proliferation, DNA damage, metastasis, hypoxia, inflammation, cell cycle, DNA repair, stemness, and quiescence). (B) Single-cell analysis indicated that SAE1 is primarily involved in regulation of the cell cycle, DNA damage and repair, and proliferation in breast cancer. (C) Correlation between SAE1 expression and cell cycle. (D) Correlation between SAE1 expression and DNA repair. (E) Correlation between SAE1 expression and cell proliferation. (F) Correlation between SAE1 expression and DNA damage.
To reveal the highly correlated genes and coexpression networks of SAE1 gene in TNBC patients, WGCNA analysis was conducted. Fifty-one TNBC samples obtained from GSE53752 database was used to build the WGCNA network (Figures 6A–E). We calculated the network topology for soft-thresholding powers from 1 to 30 to choose the best threshold. One of the most critical parameters in WGCNA network construction is the power value, which affects the average connectivity and independence of each coexpression module. Power value 8 was the lowest power for the scale-free topology. The coexpression similarity matrix was transformed into adjacency matrix, and the topological overlap matrix (TOM) was calculated. The dynamic tree cut analysis gave rise to 46 modules with different colors, among which the genes coexpressed with SAE1 belong to the green module containing 1,176 genes.
Figure 6. Identification of coexpression module genes associated with SAE1 using the WGCNA. (A) Clustering dendrogram of TNBC tissues. (B) Relationship between scale-free topology model fit and soft thresholds (powers). (C) Relationship between the mean connectivity and various soft thresholds. (D) Dendrogram of modules identified by WGCNA. (E) Topological overlap matrix among detected genes from RNAseq. (F,G) GO-BP (F) and KEGG pathway (G) network for the target genes in green model.
In order to gain further insights into the function of 1,176 genes in green module, we used DAVID database to analyze GO and KEGG pathway (Figures 6F,G). The terms of biological processes (BP) were cell-cell adhesion, intracellular protein transport, covalent chromatin modification, vesicle-mediated transport, cell division, protein autophosphorylation, protein transport, chromatin remodeling, DNA repair, mRNA processing, and protein phosphorylation. The terms of KEGG pathway terms were acute myeloid leukemia, protein processing in endoplasmic reticulum, insulin signaling pathway, endocytosis, non-small cell lung cancer, progesterone-mediated oocyte maturation, regulation of autophagy, chronic myeloid leukemia, thyroid hormone signaling pathway, and cell cycle.
Meanwhile, we analyzed the differentially expressed genes (DEGs) between TNBC tissues and adjacent normal tissues by using GSE45827, GSE65216, GSE31448, and GSE21653 databases. Figure 7A shows that the 1,176 genes in green module and DEGs in GSE45827, GSE65216, GSE31448, and GSE21653 have 10 overlap genes including minichromosome maintenance complex component 5 (MCM5), polo-like kinase 1 (PLK1), PHD finger protein 19 (PHF19), ubiquitin-conjugating enzyme E2S (UBE2S), secretion associated, Ras related GTPase 1A (SAR1A), extra spindle pole bodies homolog 1 (ESPL1), coactosin-like F-actin binding protein 1 (COTL1), antisilencing function 1B histone chaperone (ASF1B), apolipoprotein O (APOO), and transcription factor 19 (TCF19). Next, single-cell analysis was performed to investigate the functions of the 11 genes including the 10 overlap genes and SAE1 gene in breast cancer via CancerSEA website. Results indicated that these 11 genes were positively correlated with cell cycle, DNA repair, DNA damage, and proliferation (Figures 7B–F).
Figure 7. The function of SAE1 coexpression genes in breast cancer cells. (A) Venn diagrams showing relationships between DEGs and green model genes. (B) Single-cell analysis indicated that those 11 genes are primarily involved in regulation of the cell cycle, DNA damage and repair, and proliferation in breast cancer. (C) Correlation between 11 gene expression and cell cycle. (D) Correlation between 11 gene expression and DNA repair. (E) Correlation between SAE1 expression and DNA damage. (F) Correlation between 11 gene expression and cell proliferation.
Construction of a Prognostic Scoring Model Based on SUMO1-Activating Enzyme Subunit 1 and Polo Like Kinase 1 Expression
The prognostic analysis of the 10 overlap genes for TNBC patients were further evaluated (Figure 8). TNBC patients with low APOO expression exhibited significantly decreased DFS. TNBC patients with low MCM5, PHF19, PLK1, and UBE2S expression exhibited significantly increased DFS. However, there was not any significant relationship between the expression of ASF1B, COTL1, ESPL1, SAR1A, and TCF19 expression and the DFS of TNBC patients.
Figure 8. The DFS of SAE1 coexpression genes for TNBC patients. The prognostic value of (A) APOO, (B) ASF1B, (C) COTL1, (D) ESPL1, (E) MCM5, (F) PHF19, (G) PLK1, (H) SAR1A, (I) TCF19, and (J) UBE2S mRNA expression for DFS in TNBC patients by Kaplan-Meier analysis. The log-rank test was performed to evaluate survival difference with the best cutoff value.
In order to construct a risk score model for predicting DFS of TNBC, LASSO Cox regression model was used to build a prognostic classifier, which included SAE1, APOO, MCM5, PHF19, PLK1, and UBE2S (Figures 9A,B). Using the LASSO Cox regression models, we calculated a risk score for each patient based on those six gene expression: risk score = (PLK1 × 0.1122) + (SAE1 × 0.3116). Survival analysis revealed that the survival time of TNBC patients in the high-risk group was significantly shorter than those patients in the low-risk group (Figure 9C). Then, the model reliability was verified through the ROC curves, and the results exhibited that the AUC value was 0.76 (Figure 9D). The above results indicated that the SAE1 and PLK1 risk assessment model had predictive value for the prognosis of TNBC patients.
Figure 9. Construction of a classifier to predict disease-free survival in TNBC patients. (A) Partial likelihood deviance of DFS for the LASSO coefficient profiles. (B) LASSO coefficient profiles of the SAE1, APOO, MCM5, PHF19, PLK1, and UBE2S expression for DFS. (C) Kaplan-Meier survival analysis showing the impact of the classifier on DFS. (D) ROC curves compare the prognostic accuracy of the classifier in TNBC patients using AUCs at 5 years to assess prognostic accuracy.
The Level of SUMO1-Activating Enzyme Subunit 1 Protein Was Highly Expressed in Triple-Negative Breast Cancer and Predicted Poor Prognosis Through Immunohistochemical
Finally, the protein expression of SAE1 in TNBC was detected via IHC staining analysis. SAE1 protein expression was demonstrated significantly increased in 207 TNBC tissues compared with the paired adjacent normal breast tissues (Figures 10A–C). In our study, we found that SAE1 mainly localized to the cytoplasm but not the nucleus in TNBC patients (Supplementary Figure 5). The clinicopathological characteristics of 207 TNBC patients are presented in Table 1, and the correlations between relevant clinicopathological factors and SAE1 protein expression were analyzed. However, no associations between clinicopathological features and SAE1 expression were observed.
Figure 10. Immunohistochemical (IHC) analyses of SAE1 protein expression in TNBC. Representative IHC images of SAE1 protein expression in TNBC tissues (A) and normal breast tissues (B). (C) Histogram of IHC score for SAE1 in TNBC patients. (D,E) The prognostic value of SAE1 protein expression for OS (D) and DFS (E) in TNBC patients by Kaplan-Meier analysis. ***p < 0.001.
Table 1. Characteristics of triple-negative breast cancer patients and their SAE1 expression level.
Survival analysis was also applied to analyze the relationship between SAE1 protein expression and the prognosis of these 207 TNBC patients. The results of Kaplan-Meier curve implied that TNBC patients with lower protein expression of SAE1 have better prognosis in both OS and DFS than those with higher expression (Figures 10D,E). Meanwhile, univariate and multivariate Cox analyses were performed to identify the independent prognostic factor for poor survival in TNBC patients (Tables 2, 3). Lymph node metastasis, tumor grade, LVI, and SAE1 expression were determined to be independent prognostic factors for OS of TNBC patients by multivariate Cox analysis. For DFS, tumor size, lymph node metastasis, tumor grade, LVI, and SAE1 expression were identified to be independent prognostic factors in TNBC patients.
Table 2. Univariate Cox proportional hazard model for disease-free survival (DFS) and overall survival (OS) in TNBC patients with adjuvant chemotherapy.
Table 3. Multivariate Cox proportional hazard model for disease-free survival (DFS) and overall survival (OS) in TNBC patients with adjuvant chemotherapy.
Discussion
SUMO1-activating enzyme subunit 1 is a component of SUMO which are small protein modifiers capable of regulating cellular localization and function of target proteins that act as linchpins in SUMOylation, a post-translational modification involved in cellular processes such as regulation of gene transcription, DNA repair, cell-cycle regulation, and apoptosis. Potential relevancy between malignant diseases and SUMOylation pathway has been put forward and analogously the aberrant SAE1 was proved to be a partaker in cancer development. Reported overexpression of SAE1 in a stage-dependent manner in glioma by activating AKT SUMOylation-mediated signaling pathways suggested a probable role of SAE1 in cancer initiation and progression. While the reported opposite tendency of prognostic influence on gastric cancer and hepatocellular carcinoma by differential SAE1 expression indicated the uncertain association between SAE1 and malignant disease.
To study the potential mechanism of malignancy that related to SAE1, SAE1 expression analysis and prognostic analysis by Kaplan-Meier method were operated to evaluate the preliminary clinical significance based on the publicly available TCGA and GEO databases. Compared with adjacent normal tissues, the mRNA expression of SAE1 was increased in TNBC tissues, and the upregulation of SAE1 in TNBC tumor tissues was revealed with a poor prognosis (DFS) for TNBC patients.
Forkhead box M1 is a key transcription factor that plays an important role in the development of embryos, the homeostasis of mature tissues, and the occurrence and progression of malignant tumors (Bella et al., 2014). It regulates many important physiological functions, including promoting the repair of DNA damage; regulating the cell cycle G1/S phase, G2/S phase transition, and M phase process to ensure the progress of cell mitosis. FOXM1 is upregulated in many malignant tumors (Kim et al., 2006). Tumor treatment methods usually include radiotherapy or drug treatment that can damage DNA, such as platinum compounds, anthracyclines, topoisomerase inhibitors, alkylating agents, etc., which are the main drugs for clinical treatment of tumors (Johnson and Seidman, 2005; Palmieri et al., 2010). However, most DNA-damaging drugs often develop drug resistance, which is the main reason for the failure of tumor treatment (Colak and Medema, 2014). Studies have found that FOXM1 plays a key role in the resistance of DNA damaging drugs. If it is activated abnormally or its expression is elevated, it may promote the development of resistance to DNA-damaging drugs (Myatt and Lam, 2007, 2008). Here, we confirmed that FOXM1 expression was positively correlated with SAE1 expression, suggesting that SAE1 may be involved in the resistance to DNA damaging drugs. TF inhibitors might also help to improve the treatments of TNBC patients.
TCGA database was quoted synchronously to determine whether the SAE1 level would be manipulated by adjusting DNA methylation in TNBC. The result presented a significant negative correlation between SAE1 expression and SAE1 methylation site (cg14042711) level. Simultaneously, correlation analysis between SAE1 and the quantity of tumor-infiltrating immune cells in the TNBC tumor microenvironment displayed a negative outcome.
To further examine the effects of elevated SAE1 levels, CancerSEA was performed for analysis. The results indicate that SAE1 may affect the progression of breast cancer by regulating the cell cycle, DNA damage, DNA repair, and cell proliferation. Next, we implemented the WGCNA for the determination of SAE1 correlated genes as well as the enrichment analysis for gene functional study. The result showed the SAE1 coexpression genes were enriched into the green module which contains 1,176 genes that mainly clustered in cell cycle and DNA repair. The 1,176 genes in green module and DEGs in GSE45827, GSE65216, GSE31448, and GSE21653 have 10 overlap genes. Results indicated that these 11 genes including SAE1 and 10 overlap genes were also positively correlated with cell cycle, DNA repair, DNA damage, and proliferation. This is consistent with contemporary knowledge that SUMOs such as SAE1 are essential for the regulation of several cellular processes, including transcriptional regulation, transcript processing, genomic replication, and DNA damage repair (Sacher et al., 2006; Hang et al., 2011; Psakhye and Jentsch, 2012; Sarangi et al., 2014; Schimmel et al., 2014; Sarangi and Zhao, 2015). The expressions of these 10 genes were further examined by survival analysis, and the results indicated that SAE1, APOO, MCM5, PHF19, PLK1, and UBE2S were significantly correlated with the prognosis. Meanwhile, we identified two core genes including SAE1 and PLK1 closely related to the DFS of TNBC and constructed a survival prediction model based on SAE1 and PLK1 expression. The AUC value of the ROC curves in this model was 0.76, which proved that the model was reliable.
Finally, our IHC findings suggested that the protein expression of SAE1 was upregulated in TNBC tissues of 207 patients. According to survival analysis, SAE1 protein was associated with both OS and DFS separately. These observations collectively provide a new idea for the prognosis of TNBC patients.
Conclusion
To generalize, detection on expression patterns of SAE1 and its prognostic significance in TNBC was primarily implemented by bioinformatic analysis and IHC which showed a positive correlation between SAE1 downregulation and better prognosis. SAE1 expression correlated with DNA methylation and TF. It is therefore of great clinical importance to investigate DNA methyltransferase and TF inhibitors that might decrease the expression of SAE1.
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 authors.
Ethics Statement
Written informed consent was obtained from the individual(s) for the publication of any potentially identifiable images or data included in this article.
Author Contributions
QW and WZ: conceptualization. LD: methodology. QL, YYL, and HXL: software. LX, LL, and YC: validation. JH and MJ: formal analysis. HX and JZ: data curation. QW and YXL: writing—original draft preparation. CS and YaL: writing—review and editing. All authors contributed to the article and approved the submitted version.
Funding
This research was funded by the National Science Foundation for Young Scientists of China (82003095), the Natural Science Foundation of Fujian Province (2020J01311402), the Joint Funds for the Innovation of Science and Technology, Fujian Province (2019Y9054), the Natural Science Foundation of Fujian Province (2019J01278), the Fujian Young Teacher Education Research Project, Fujian Province (JAT170228), the National Natural Science Foundation of China (81672817), the Fujian Provincial Key Laboratory of Hepatic Drug Research (KFLX2020001), the Scientific Research Innovation Program “Xiyuanjiang River Scholarship” of College of Life Sciences, Fujian Normal University, and Special Project of Clinical Application Research in Fuzong Clinical Medical College of Fujian Medical University (2019L10) corresponded by LD.
Conflict of Interest
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
Publisher’s Note
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.
Supplementary Material
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fcell.2021.729211/full#supplementary-material
Supplementary Figure 1 | The mRNA expression of SAE1 in breast cancer. The mRNA expression of SAE1 in breast cancer patients based on GSE31448 (A), GSE45827 (B), and GSE65216 (C). ns, p > 0.05; ∗p < 0.05; ∗∗∗p < 0.001.
Supplementary Figure 2 | The pathway activity of SUMOylation in breast cancer. The pathway activity of SUMOylation in breast cancer based on GE21635. ns, p > 0.05; ∗p < 0.05.
Supplementary Figure 3 | The OS and DFS of protein SUMOlyation pathway for TNBC patients. The prognostic value of protein SUMOlyation pathway for OS (A) and DFS (B) in TNBC patients by Kaplan-Meier analysis.
Supplementary Figure 4 | Forkhead box M1 (FOXM1) regulates the expression of SAE1. The expression of FOXM1 (A) and SAE1 (B) in breast cancer cells transfected with control or FOXM1 siRNA base on GSE55204. Correlation between SAE1 expression and FOXM1 expression.
Supplementary Figure 5 | Forkhead box M1 (FOXM1) regulates the expression of SAE1. Another two representative IHC images of SAE1 protein expression in TNBC tissues.
Footnotes
- ^ http://www.cbioportal.org/
- ^ http://www.ncbi.nlm.nih.gov/geo/
- ^ http://xena.ucsc.edu/
- ^ https://david.ncifcrf.gov/
- ^ http://biocc.hrbmu.edu.cn/CancerSEA/
- ^ http://dbtoolkit.cistrome.org
References
Abbas, R., McColl, K. S., Kresak, A., Yang, M., Chen, Y., Fu, P., et al. (2015). PIAS3 expression in squamous cell lung cancer is low and predicts overall survival. Cancer Med. 4, 325–332. doi: 10.1002/cam4.372
Amente, S., Lavadera, M. L., Palo, G. D., and Majello, B. (2012). SUMO-activating SAE1 transcription is positively regulated by Myc. Am. J. Cancer Res. 2, 330–334.
Bergamaschi, A., Madak-Erdogan, Z., Kim, Y. J., Choi, Y. L., Lu, H., Katzenellenbogen, B. S. (2014). The forkhead transcription factor FOXM1 promotes endocrine resistance and invasiveness in estrogen receptor-positive breast cancer by expansion of stemlike cancer cells. Breast Cancer Res. 16:436. doi: 10.1186/s13058-014-0436-4
Bella, L., Zona, S., Nestal de Moraes, G., and Lam, E. W. (2014). FOXM1: a key oncofoetal transcription factor in health and disease. Semin. Cancer Biol. 29, 32–39. doi: 10.1016/j.semcancer.2014.07.008
Chang-Qing, Y., Jie, L., Shi-Qi, Z., Kun, Z., Zi-Qian, G., Ran, X., et al. (2020). Recent treatment progress of triple negative breast cancer. Prog. Biophys. Mol. Biol. 151, 40–53.
Colak, S., and Medema, J. P. (2014). Cancer stem cells–important players in tumor therapy resistance. FEBS J. 281, 4779–4791. doi: 10.1111/febs.13023
da Silva, J. L., Cardoso Nunes, N. C., Izetti, P., de Mesquita, G. G., and de Melo, A. C. (2020). Triple negative breast cancer: a thorough review of biomarkers. Crit. Rev. Oncol. 145:102855.
Gruosso, T., Mieulet, V., Cardon, M., Bourachot, B., Kieffer, Y., Devun, F., et al. (2016). Chronic oxidative stress promotes H2AX protein degradation and enhances chemosensitivity in breast cancer patients. EMBO Mol. Med. 8, 527–549. doi: 10.15252/emmm.201505891
Hang, L. E., Liu, X., Cheung, I., Yang, Y., and Zhao, X. (2011). SUMOylation regulates telomere length homeostasis by targeting Cdc13. Nat. Struct. Mol. Biol. 18, 920–926. doi: 10.1038/nsmb.2100
He, T., McColl, K., Sakre, N., Chen, Y., Wildey, G., and Dowlati, A. (2018). Post-transcriptional regulation of PIAS3 expression by miR-18a in malignant mesothelioma. Mol. Oncol. 12, 2124–2135. doi: 10.1002/1878-0261.12386
Hickey, C. M., Wilson, N. R., and Hochstrasser, M. (2012). Function and regulation of SUMO proteases. Nat. Rev. Mol. Cell Biol. 13, 755–766. doi: 10.1038/nrm3478
Huang da, W., Sherman, B. T., and Lempicki, R. A. (2009). Systematic and integrative analysis of large gene lists using DAVID bioinformatics resources. Nat. Protoc. 4, 44–57. doi: 10.1038/nprot.2008.211
Jhan, J. R., and Andrechek, E. R. (2017). Triple-negative breast cancer and the potential for targeted therapy. Pharmacogenomics 18, 1595–1609. doi: 10.2217/pgs-2017-0117
Johnson, M. L., and Seidman, A. D. (2005). Emerging targeted therapies for breast cancer. Oncology 19, 611–618.
Kim, I. M., Ackerson, T., Ramakrishna, S., Tretiakova, M., Wang, I. C., Kalin, T. V., et al. (2006). The Forkhead Box m1 transcription factor stimulates the proliferation of tumor cells during development of lung cancer. Cancer Res. 66, 2153–2161. doi: 10.1158/0008-5472.can-05-3003
Komatsu, M., Yoshimaru, T., Matsuo, T., Kiyotani, K., Miyoshi, Y., Tanahashi, T., et al. (2013). Molecular features of triple negative breast cancer cells by genome-wide gene expression profiling analysis. Int. J. Oncol. 42, 478–506. doi: 10.3892/ijo.2012.1744
Langfelder, P., and Horvath, S. (2008). WGCNA: an R package for weighted correlation network analysis. BMC Bioinform. 9:559. doi: 10.1186/1471-2105-9-559
Li, J., Cui, J., Zhang, J., Liu, Y., Han, L., Jia, C., et al. (2015). PIAS3, an inhibitor of STAT3, has intensively negative association with the survival of gastric cancer. Int. J. Clin. Exp. Med. 8, 682–689.
Maire, V., Baldeyron, C., Richardson, M., Tesson, B., Vincent-Salomon, A., Gravier, E., et al. (2013). TTK/hMPS1 is an attractive therapeutic target for triple-negative breast cancer. PLoS One 8:e63712. doi: 10.1371/journal.pone.0063712
Mathe, A., Wong-Brown, M., Morten, B., Forbes, J. F., Braye, S. G., Avery-Kiejda, K. A., et al. (2015). Novel genes associated with lymph node metastasis in triple negative breast cancer. Sci. Rep. 5:15832.
Merdad, A., Karim, S., Schulten, H. J., Dallol, A., Buhmeida, A., Al-Thubaity, F., et al. (2014). Expression of matrix metalloproteinases (MMPs) in primary human breast cancer: MMP-9 as a potential biomarker for cancer invasion and metastasis. Anticancer Res. 34, 1355–1366.
Myatt, S. S., and Lam, E. W. (2007). The emerging roles of forkhead box (Fox) proteins in cancer. Nat. Rev. Cancer 7, 847–859. doi: 10.1038/nrc2223
Ong, J. R., Bamodu, O. A., Khang, N. V., Lin, Y. K., Yeh, C. T., Lee, W. H., et al. (2021). SUMO-activating enzyme subunit 1 (SAE1) is a promising diagnostic cancer metabolism biomarker of hepatocellular carcinoma. Cells 10:178. doi: 10.3390/cells10010178
Palmieri, C., Krell, J., James, C. R., Harper-Wynne, C., Misra, V., Cleator, S., et al. (2010). Rechallenging with anthracyclines and taxanes in metastatic breast cancer. Nat. Rev. Clin. Oncol. 7, 561–574. doi: 10.1038/nrclinonc.2010.122
Pichler, A., Fatouros, C., Lee, H., and Eisenhardt, N. (2017). SUMO conjugation - a mechanistic view. Biomol. Concepts 8, 13–36. doi: 10.1515/bmc-2016-0030
Polimeno, L., Francavilla, A., Piscitelli, D., Fiore, M. G., Polimeno, R., Topi, S., et al. (2020). The role of PIAS3, p-STAT3 and ALR in colorectal cancer: new translational molecular features for an old disease. Eur. Rev. Med. Pharmacol. Sci. 24, 10496–10511.
Psakhye, I., and Jentsch, S. (2012). Protein group modification and synergy in the SUMO pathway as exemplified in DNA repair. Cell 151, 807–820. doi: 10.1016/j.cell.2012.10.021
Sabatier, R., Finetti, P., Adelaide, J., Guille, A., Borg, J. P., Chaffanet, M., et al. (2011a). Down-regulation of ECRG4, a candidate tumor suppressor gene, in human breast cancer. PLoS One 6:e27656. doi: 10.1371/journal.pone.0027656
Sabatier, R., Finetti, P., Cervera, N., Lambaudie, E., Esterni, B., Mamessier, E., et al. (2011b). A gene expression signature identifies two prognostic subgroups of basal breast cancer. Breast Cancer Res. Treat. 126, 407–420. doi: 10.1007/s10549-010-0897-9
Sacher, M., Pfander, B., Hoege, C., and Jentsch, S. (2006). Control of Rad52 recombination activity by double-strand break-induced SUMO modification. Nat. Cell Biol. 8, 1284–1290. doi: 10.1038/ncb1488
Sarangi, P., Bartosova, Z., Altmannova, V., Holland, C., Chavdarova, M., Lee, S. E., et al. (2014). Sumoylation of the Rad1 nuclease promotes DNA repair and regulates its DNA association. Nucleic Acids Res. 42, 6393–6404. doi: 10.1093/nar/gku300
Sarangi, P., and Zhao, X. (2015). SUMO-mediated regulation of DNA damage repair and responses. Trends Biochem. Sci. 40, 233–242. doi: 10.1016/j.tibs.2015.02.006
Schimmel, J., Eifler, K., Sigurðsson, J. O., Cuijpers, S. A., Hendriks, I. A., Verlaan-de Vries, M., et al. (2014). Uncovering SUMOylation dynamics during cell-cycle progression reveals FoxM1 as a key mitotic SUMO target protein. Mol. Cell 53, 1053–1066. doi: 10.1016/j.molcel.2014.02.001
Song, N., Gu, X., Wang, Y., Chen, Z., and Shi, L. (2014). Lentivirus-mediated siRNA targeting sae1 induces cell cycle arrest and apop- tosis in colon cancer cell RKO. Mol. Biol. 48, 107–116.
Sporikova, Z., Koudelakova, V., Trojanec, R., and Hajduch, M. (2018). Genetic markers in triple-negative breast cancer. Clin. Breast Cancer 18, e841–e850.
Terunuma, A., Putluri, N., Mishra, P., Mathé, E. A., Dorsey, T. H., Yi, M., et al. (2014). MYC-driven accumulation of 2-hydroxyglutarate is associated with breast cancer prognosis. J. Clin. Invest. 124, 398–412.
Tseng, J. C., Huang, S. H., Lin, C. Y., Wang, B. J., Huang, S. F., Shen, Y. Y., et al. (2020). ROR2 suppresses metastasis of prostate cancer via regulation of miR-199a-5p-PIAS3-AKT2 signaling axis. Cell Death Dis. 11:376.
Vagia, E., Mahalingam, D., and Cristofanilli, M. (2020). The landscape of targeted therapies in TNBC. Cancers 12:916. doi: 10.3390/cancers12040916
Wang, C., Xiao, Y., Lao, M., Wang, J., Xu, S., Li, R., et al. (2020). Increased SUMO-activating enzyme SAE1/UBA2 promotes glycolysis and pathogenic behavior of rheumatoid fibroblast-like synoviocytes. JCI Insight 5:e135935.
Wang, Q., Zhang, X., Chen, L., Weng, S., Xia, Y., Ye, Y., et al. (2019). Regulation of the expression of DAPK1 by SUMO pathway. Biomolecules 9:151. doi: 10.3390/biom9040151
Wang, Q. S., Li, F., Liao, Z. Q., Li, K., Yang, X. L., Lin, Y. Y., et al. (2019). Low level of Cyclin-D1 correlates with worse prognosis of clear cell renal cell carcinoma patients. Cancer Med. 8, 4100–4109. doi: 10.1002/cam4.2313
Wilson, T. R., Udyavar, A. R., Chang, C. W., Spoerke, J. M., Aimi, J., Savage, H. M., et al. (2019). Genomic alterations associated with recurrence and TNBC subtype in high-risk early breast cancers. Mol. Cancer Res. 17, 97–108. doi: 10.1158/1541-7786.mcr-18-0619
Yang, Y., Liang, Z., Xia, Z., Wang, X., Ma, Y., Sheng, Z., et al. (2019). SAE1 promotes human glioma progression through activating AKT SUMOylation-mediated signaling pathways. Cell Commun. Signal. 17:82.
Yuan, H., Yan, M., Zhang, G., Liu, W., Deng, C., Liao, G., et al. (2019). CancerSEA: a cancer single-cell state atlas. Nucleic Acids Res. 47, D900–D908.
Zhang, M., Jiang, D., Xie, X., He, Y., Lv, M., and Jiang, X. (2019). miR-129-3p inhibits NHEJ pathway by targeting SAE1 and represses gastric cancer progression. Int. J. Clin. Exp. Pathol. 12, 1539–1547.
Keywords: TNBC, SAE1, WGCNA, prognosis, cell cycle
Citation: Wang Q, Zhong W, Deng L, Lin Q, Lin Y, Liu H, Xu L, Lu L, Chen Y, Huang J, Jiang M, Xiao H, Zhang J, Li H, Lin Y, Song C and Lin Y (2021) The Expression and Prognostic Value of SUMO1-Activating Enzyme Subunit 1 and Its Potential Mechanism in Triple-Negative Breast Cancer. Front. Cell Dev. Biol. 9:729211. doi: 10.3389/fcell.2021.729211
Received: 22 June 2021; Accepted: 25 August 2021;
Published: 21 September 2021.
Edited by:
Erika Ruiz-Garcia, National Institute of Cancerology (INCAN), MexicoReviewed by:
Jennifer Diamond, University of Colorado Cancer Center, United StatesChang Li, University of Kentucky, United States
Copyright © 2021 Wang, Zhong, Deng, Lin, Lin, Liu, Xu, Lu, Chen, Huang, Jiang, Xiao, Zhang, Li, Lin, Song and Lin. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Yao Lin, yaolin@fjtcm.edu.cn; yaolin@fjnu.edu.cn; Chuangui Song, Songcg1971@hotmail.com; Yuxiang Lin, lyx0520@fjmu.edu.cn
†These authors have contributed equally to this work