Abstract
Cancer stem cells (CSCs) represent a subset of tumor cells that are responsible for recurrence and metastasis of tumors. These cells are resistant to radiotherapy and chemotherapy. Immunotherapeutic strategies that target CSCs specifically have provided initial results; however, the mechanism of action of these strategies is unclear. The data were requested from The Cancer Genome Atlas and Genotype-Tissue Expression, followed with the survival analysis and weighted gene co-expression network analysis to detect survival and stemness related genes. Patients were divided into three groups based on their immune status by applying single sample GSEA (ssGSEA) with proven dependability by ESTIMATE analysis. The filtered key genes were analyzed using oncomine, GEPIA, HPA, qRT-PCR, and functional analysis. Patients in a group with a higher stemness and a lower immune infiltration showed a worse overall survival probability, stemness and immune infiltration characteristics of breast cancer progressed in a non-linear fashion. Thirteen key genes related to stemness and immunity were identified and the functional analysis indicated their crucial roles in cell proliferation and immune escape strategies. The qRT-PCR results showed that the expression of PIMREG and MTFR2 differed in different stages of patients. Our study revealed a promising potential for CSC-target immunotherapy in the early stage of cancer and a probable value for PIMREG and MTFR2 as biomarkers and targets for immunotherapy.
Introduction
Breast cancer (BC) has the highest incidence rate and mortality rate among female malignant tumors, which impacts women’s health significantly (1). It is considered to be heterogeneous depending on molecular subtype and on different stages of cancer progression (2). This heterogeneity poses are changed during treatment, even through various treatment strategies have been developed based on different pathological types (3), especially for the triple-negative BC (4). Thus, further explorations are needed to identify new markers for guiding individualized treatment.
Growing evidence has established the presence of a subpopulation of cancer cells with stem-like properties in most human malignancies, frequently referred to as ‘cancer stem cells’ (CSCs), which possess the long-term ability to initiate and repopulate tumors (5, 6). Diverse mechanisms by which CSCs manage to survive through various strategies including tumor initiation, metastatic reactivation, oncogene- and immune-targeted therapy resistance have been unvealed (4, 7).
Immuno-resistance is one of the main features of tumors that helps them escape immunosurveillance and evade eradication by resisting immunosuppression (8, 9). A myriad of strategies have been discovered in the tumor cells that allow them to circumvent the immune attack, including genetic and epigenetic alterations in the genome of tumor cells that reduce immune recognition and promote protective microenvironment (5, 10).
Evidence has emerged that CSCs have a potential role in regulating their immune characteristics (9, 10), while the molecular mechanism is unclear. In this study, we focused on CSCs in BC. To this end, differentially expressed genes (DEGs) were screened using The Cancer Genome Atlas (TCGA) and Genotype-Tissue Expression (GTEx) databases, mRNAsi index and WGCNA analysis were used in turns to profile the association between stemness of tumor and clinical characters and to identify genes that are related to stemness. Next, we applied immune infiltration analysis to the filtered immune-related genes from stemness and survival-related genes. Finally, the qRT-PCR analysis was used to verify that the expression of PIMREG and MTFR2 was aberrant and diverse in different clinical stages.
Materials and Methods
Design and Data Processing of This Study
The schematic flowchart of the present study is shown in Figure 1. Gene expression and clinical information of patients and healthy individuals was collected from the common database by applying strict filters. Then the computational biology tools were applied, to evaluate stemness and immune characteristics, and biomarkers for patients with BC were revealed. The expression of key genes was evaluated by qRT-PCR.
Figure 1
Data Acquisition
The comprehensive data of BC including 1,104 tissues of patients and 113 non-tumor tissues was downloaded from the TCGA database (https://portal.gdc.cancer.gov September 26, 2019). Due to the lack of corresponding samples from non-tumor tissues, we utilized GTEx project from where 80 samples of healthy individuals were obtained. RNA-seq results of healthy tissue samples and cancer samples were combined into a matrix file using a merge script in the Perl language (http://www.perl.org/). The mRNAsi were acquired from Tathiane M. Malta et al. (11). The 29 marker gene sets for immune cell types and signaling pathway activation were obtained from Bao X et al. and Liu Z et al. (12, 13).
Differentially Expressed Genes
The “edgeR” R package was used to screen DEGs of normal breast and cancer samples with the following selection criteria: FDR < 0.05, and |log2 FC| > 1. The values of genes calculated by limma R package, and genes with expression of < 0.5 were deleted. The heatmap and Kaplan-Meier (K-M) curves were drawn.
WGCNA and Module Merge
Co-expression network was established using the WGCNA R package based on the DEGs. First, RNA-seq data were matched with corresponding mRNAsi and filtered by hierarchical cluster analysis to detect outliers, the Pearson correlation matrix was constructed by correlating coefficient of genes and then transforming them into weighted adjacency matrix using the power function: apq=|cpq|β (cpq=Pearson’s correlation between gene p and gene q, apq=adjacency between gene p and gene q, and β = soft threshold). The best soft threshold (soft threshold = 4) which was selected according to scale independence and mean connectivity for achieving a scale-free co-expression network. The weighted adjacency matrix was then transformed into the topological overlap matrix, DEGs were allocated into different modules for average linkage hierarchical clustering and similar genes were grouped into one module, module dendrogram was drawn with the minimum size (genome), modules with similar heights (cutoff < 0.3) were merged. As a correlation value between genes and sample traits, gene significance (GS) was calculated based on statistical significance, which was determined using the relevant p-values in the linear regression between gene expression and clinical phenotypes (mRNAsi and EREG-mRNAsi). The significant modules related to mRNAsi were selected according to Module significance which was defined as the average GS within the module and revealed the correlation between the module and sample traits.
Identification of mRNAsi- and Survival-Associated Genes
Module membership (MM) was defined as the correlation between the module’s own genes and gene expression profiles. The key genes associated with mRNAsi were screened by GS and MM defined as cor. gene MM > 0.8 and cor. gene GS > 0.5. The interactive genes between survival-related DEGs and mRNAsi-related genes were finally identified and depicted using venn charts.
Immune Infiltration Grouping of BC Patients
BC gene set was prepared as a gmt file for further quantitative measurements of the immune activation status. The RNA-seq data of individual cancer samples was transformed into enrichment scores of each immune-related term by Single-Sample Gene Set Enrichment Analysis (ssGSEA) in the R package gsva. Tumors samples with qualitatively different enrichment scores were divided into low, median, high infiltration clusters by using hierarchical clustering analysis in sparcl R package. Results were presented as a color dendrogram and heatmap.
Tumor Microenvironment (TME) Analysis of BC
The TME scores were calculated based on BC gene set and tumor purity was predicted by using estimate R package. The assessment of TME were divided into four clusters (stromal score, immune score, estimate score and tumor purity) (14, 15). The heatmap and the violin plot were conducted to further explore the relationship between immune groups and TME.
Identification of Key Genes
All interactive genes associated with survival time and mRNAsi of Venn chart were enrolled for correlation analysis with three infiltration clusters. Genes with no significant difference between immune clusters and mRNAsi were excluded (P ≥ 0.05). Correlation analysis and survival analysis of these genes were conducted. The results were shown in box plots, heatmaps, and K-M curves. Oncomine (http://www.oncomine.org), GEPIA (http://gepia.cancer-pku.cn/), and The Human Protein Atlas database (HPA; https://www.proteinatlas.org/) was used to inspect differential expression of key genes between BC and healthy tissues and between different tumors.
Functional Analysis of Key Genes
Protein-protein interaction (PPI) network analysis was applied to show the relationship between different proteins of these key genes by using String (https://string-db.org/). The functional enrichment analysis was also conducted using clusterProfiler R package to investigate the biological functions of key genes. Gene ontology (GO) functional annotations and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment were used in the study with the threshold values: P < 0.05, and FDR < 0.05. The results were shown in PPI network graph and dot plots. The relation of co-expressed key genes and clinical characters are analyzed by houyunhuang/ggcor R package.
Cell Culture
The human breast cancer cell lines, MCF-7 and MDA-MB-231 were obtained from the Chinese Academy of Sciences (Shanghai, China), cell lines were maintained in DMEM (Invitrogen, Carlsbad, CA, USA) with 10% fetal bovine serum (FBS; Sigma-Aldrich, St. Louis, MO, USA), penicillin (100 U/ml), and streptomycin (100 g/ml; Life Technologies, Grand Island, NY, USA) at 37°C with 5% CO2.
Patients
BC samples and corresponding para-cancerous tissues were collected from patients who underwent modified radical mastectomy from January 2020 to March 2020 from the department of endocrine and breast surgery in the First Affiliated Hospital of Chongqing Medical University. Tissues were excised and immediately transferred into liquid nitrogen. All patients were informed and written informed consent were provided. The study was conducted according to the clinical practice guidelines of the International Conference on Harmonization and the Declaration of Helsinki. All patients were diagnosed with triple-negative BC with no evidence of distant metastasis, only patients with Tumor Node Metastasis (TNM) stages II-III were included in the study, certified by two pathologists. BC was divided into two subtypes based on their TNM stages. Requisite clinical data were acquired from hospital records and pathology reports. This study protocol was approved by the ethical committee of Affiliated Hospital of Chongqing Medical University. This study protocol registered with Chinese Clinical Trial Registry (http://www.chictr.org.cn/showproj.aspx?proj=19710; Date of registration: 25/09/2017; Registration number: ChiCTR-PDN-17012784) and was approved by the ethical committee of Affiliated Hospital of Chongqing Medical University (approval number: 2020-119).
Total RNA Extraction
Total RNA from BC tissues was extracted using an UNIQ-10 column Total RNA Extractio Kit (Sangon Biotech). The RNA concentration and purity were assessed using a SMA4000 microspectrophotometer (Merinton Instrument, Inc) and by RNA electrophoresis with DYY-6C electrophoresis apparatus (Liuyi. Beijing).
Reverse Transcription and qRT-PCR Quantification
RNA with concentration ranging from 91.84 ng/µl to 1325.94 ng/µl from human BC tumor and para-carcinoma tissues and were reversed-transcribed using a RR047A cDNA synthesis kit (TaKaRa, China). Quantitative PCR was performed for CEP55, MTFR2, and PIMREG using a 2X SG Fast qPCR Master Mix (High Rox, B639273, BBI) with the Step One Plus fluorescence quantitative PCR instrument (ABI, Foster, CA, USA), GAPDH is used for internal control gene. The primers designed by Primer Premier 5.0 are listed as below: MTFR2-F: CTCCTCCACCACTTCCTCCTCAG; MTFR2-R: CGCTCAATTGCACGAAGCTTAACC; PIMREG-F: GAGTGCTTTGGGTGCCGTGTC; PIMREG-R: CCGCCTTGATCGCCGTAATGG; CEP55-F: GTGGGGATCGAAGCCTAGTA; CEP55-R: TCATACACGAGCCACTGCTG.
Immunohistochemical
The stained images of IHC and the corresponding information in breast tissues and tumors were downloaded from The Human Protein Atlas (https://www.proteinatlas.org). Samples from two patients were included in each group, the IHC procedure were conducted as descripted in the official IHC protocol by Human protein Atlas (https://www.proteinatlas.org/download/IHC_protocol.pdf), in brief, the 4-μm-thick tissue sections from the tissues were pretreated and then stained with the antibodies (CEP55: Sigma-Aldrich, HPA023430, MTFR2: Sigma-Aldrich, HPA029792, and PIMREG: Sigma-Aldrich, HPA043783). Heat Induced Epitope Retrieval (HIER) was used for antigen retrieval. In the Immunohistochemical staining process, slides were first rinsed in wash buffer. then they were Incubated with Ultra V Block, primary antibody, and labeled HRP polymer in order. Next, they were developed in DAB solution, counterstained in hematoxylin, rinsed, dehydrated, and coversliped. The information of stain, intensity and quantity of the photo description were displayed as the results that are shown in The Human Protein Atlas, eg: https://www.proteinatlas.org/ENSG00000129195-PIMREG/tissue/breast#img, to be noted, the phenomenon of the information that is not fully matched with the photo might be caused by the process of manually adjustment: https://www.proteinatlas.org/about/help#2).
Statistical Analysis
The K-M analysis was applied to show the survival difference of mRNAsi groups and immune clusters. Wilcox analysis was utilized to analyze the difference of mRNAsi between normal and cancer groups. Krudkal-Wallis H test was conducted to reveal the relationship between mRNAsi, infiltration clusters, TME scores, tumor purity and clinical information like TNM stage (P<0.05). The beeswarm plots and violin plots were made to visualize the results.
The differences in gene expression between tumor and normal tissues, and para comparison between subgroups were expressed as mean with standard deviations and compared through Student’s paired t-test and Student’s t-test accordingly by applying GraphPad Prism 8, P < 0.05 indicated statistically significant difference.
Results
Screening for Survival-Related DEGs
Samples for which with RNA-seq and comprehensive information was available were enrolled. A total of 2261 different expression genes were screened with 989 up-regulated genes and 1272 down-regulated genes and compared with normal samples. The results were visualized by heat map (Figure 2), in which 212 survival-related DEGs were filtered (including 89 upregulated genes and 123 down-regulated genes).
Figure 2
Correlation Analysis of mRNAsi and Clinical Characteristics
mRNAsi was revealed by Tathiane M. Malta et al. as an index for evaluating the stemness of tumor cells, thus it is also considered as a quantitative description of CSCs. Wilcox analysis showed that the mRNAsi level in cancer group was significantly higher than that of normal group (Figure 3A). The K-M curve showed that patients in low mRNAsi group had a longer overall-survival time within 5 year-follow-up (Figure 4A). As shown in Figures 3C, D, the score for mRNAsi had a rising tendency in a non-linear manner in the progression of BC, to be noted, the dynamic change of stemness of BC tissue from stage II to stage III were inverse comparing to the overall tendency (stage I to stage IV).
Figure 3
Figure 4
Identification of Significant Modules and mRNAsi-Associated Genes
In order to identify biologically significant genes that were associated with BC stemness, a gene co-expression network was conducted by performing WGCNA. DEGs with similar expression patterns were allocated into one module by WGCNA (Figure 5A), soft threshold in this study was set at 4 in order to ensure a scale-free network (Figures 5C, D), and nine significant gene modules were obtained for subsequent analysis as shown in Figure 5B. Next, we explored the correlation between these nine modules and the clinical traits (including mRNAsi and EREG-mRNAsi). The result showed that three gene modules were correlated with mRNAsi (the value of correlation > 0.6), among these, the turquoise module with a positive correlation at 0.83 was most significantly related to mRNAsi, the blue and brown module had significantly negative correlation with mRNAsi (Figure 5B). Next, 100 key genes were screened from turquoise module with a threshold defined as cor.MM > 0.8 and cor.GS > 0.5 (Figure 5E). Finally, 18 interactive genes were selected from the intersection of 212 survival-related DEGs as well as 100 mRNAsi-related genes.
Figure 5
Immune Phenotype Landscape in the TME of BC
The immune cell infiltration status of each sample was assessed by ssGSEA and transformed into enrichment scores based on its transcriptomes according to 29 marker gene sets. All 1104 BC samples were allocated into three hierarchical clusters (low infiltration: 11 patients; intermediate infiltration: 618 patients; and high infiltration: 475 patients) (Figure 6A). ESTIMATE method was then applied to testify the reliability of this model, which evaluated the conditions of cellular composition in each sample by multiple indicators, including stromal score, immune score, tumor purity, and an ESTIMATE score that is regarded as sum of the first three scores. The results showed that patients in high infiltration cluster had higher stromal score, higher immune score, lower tumor purity, and higher overall ESTIMATE score (Figures 6B–E). Survival analysis was also applied to assess the of prognosis value of ESTIMATE scores. As shown in Figures 4B–E, high immune score, low tumor purity, and high ESTAMATE score had a longer overall-survival time. As shown in Figure 7, stromal score and estimate score showed a downward tendency in a non-linear manner in the tumor progression, while estimate score showed a rising tendency in a non-linear manner in the progression of Lymph Node. To be noted, the change of stomal score from stage II to stage III and ESTAMATE score from T2 to T3 showed inverse correlation compared with the global trend. The Stromal score of N2 reaches the top.
Figure 6
Figure 7
Key Genes Identification and Correspondence Functional Analysis
Of 18 mRNAsi-related genes that were analyzed for correlation with three immune cell infiltration clusters, 13 key genes were finally enrolled (Figure 8). There was statistical significance of the expression of ANLN, CCNB2, CEP55, FAM83D, MELK, MTFR2, PIMREG, RAD54L, RAD51, SHCBP1, SKA1, UBE2T, and SKA3 within three immune cell infiltration clusters. As shown in Figure 9, the expression levels for all 13 key genes were higher in tumor samples than in the normal samples. The K-M curves showed that patients with higher level of expression of these genes had a poor overall survival time (Figure 10). The results indicated that these 13 key genes were associated with mRNAsi and immunization as good prognostic factors.
Figure 8
Figure 9
Figure 10
In order to further explore the molecular function of the 13 key genes identified, PPI and enrichment analysis were conducted. As shown in Figures 11B, C, CEP55 proteins had one of the highest node number (node number is 12). The clinical significance of these co-expressed key genes is shown in Figure 11. GO functional enrichment analysis indicated these genes were related to cell proliferation, such as organelle fission, chromosome segregation, spindle, midbody, and chromosomal region (Figure 12A). The results of KEGG showed the top enriched terms were homologous recombination and Fanconi anemia pathway and indicated a close correlation with DNA damage and repair (Figure 12B). The result of TIMER showed that most of the key genes had strong correlation with B lymphocyte and dendritic cells’ function (Figure 3E). Through differentially expression analysis of pan-cancer and normal samples by Oncomine and protein expression scan by HPA, it is noted that the evaluation of protein expression score was based on immunohistochemical, which were manually scored with regard to staining intensity including negative, weak, moderate or strong, and fraction of stained cells including <25%, 25–75%, and >75%, next the combination of intensity and fractions was converted into an protein expression level score and listed as follows: weak combined with either 25 - 75% or 75% and moderate <25% were shown as low, which included PIMREG and MTFR2 in breast tissue; moderate combined with either 25 - 75% or 75% and strong <25% were shown as medium, which contained the rest of the images. Protein expression values were also manually adjusted by expert annotators when necessary. It was found that PIMREG, MTFR2, and CEP55 were overexpressed in BC and also in many other cancers (Figures 13A, H–M).
Figure 11
Figure 12
Figure 13
Validation of Key Genes Expression Using qRT-PCR
The expression levels of PIMREG, MTFR2, and CEP55 were assessed by qRT-PCR, which showed a significantly higher level of expression for these genes in BC samples compared with that in corresponding adjacent breast samples. We also investigated the relationship between gene expressions and TNM stage. Significantly different expression levels of PIMREG and MTFR2 were identified between TNM stage II group and stage III group as shown in Figures 13B–G. Next, we investigated whether expression levels of PIMREG, CEP55, and MTFR2 is related to the expressions of BCSC markers in breast cancer cell lines. Two breast-cancer cell lines were compared, including MCF7 and MDA-MB-231, which are respectively known as negative and positive cell lines for BCSC markers [47]. Two well-known breast cancer stem-cell markers, CD44 and ALDH1 family member A1, were applied to verify this result. RT-PCR indicated that the expression level of PIMREG, CEP55, and MTFR2 are relatively high in BCSC-marker-positive cells (MDA-MB-231) in comparing with BCSC-negative cells (MCF-7).
Discussion
Typical treatments for BC besides surgery include endocrine therapy, chemical therapy, and radiotherapy; these have led to increased survival rates in majority of patients since most of the tumor cell express and respond to receptors for estrogen and progesterone receptor. However, as the most aggressive subtype, triple-negative BC have limited treatment options and worse prognosis. The presence of CSCs was discovered as one of the main causes for treatment resistance and there is lack of targets for CSCs in advanced cancers. In fact, one of the main manifestations in the cancer progression is the gradual loss of the differentiated phenotype and gaining stem cell-like characteristics (5). Evidence showed that the frequency of CSC is increasing with tumor progression in multiple solid cancer types (16); however, this process is non-linear in BC due to the dynamic negative feedback control of tumor cell populations and high heterogeneity within cells in BC tissue (17). Therefore, it is crucial to observe the dynamic change of the stemness characteristics in the progression of BC. In this study, we discovered that the mRNAsi score was the lowest in stage I and reaches the highest in stage IV, which matched with the data published by earlier researchers. To be noted, tumor reached the first peak in stage II (and T2) and decreased in stage III (and T3) of the stem cell characteristics, this provided an insight into a specific timeline for detecting the vitality of breast CSC (17).
Diverse mechanisms have been unraveled that could help CSCs to survive under hostile environment, including enhancing DNA-repair capacity (18), among which are overexpression of multifunctional efflux transporters and aberrant activation of developmental pathways (19, 20). We applied mRNAsi-based WGCNA and identified the module that was associated with BCSC characteristics. Functional annotations of the module were primarily associated with cell cycle control, cell proliferation characteristics, and DNA-damage repair.
Immunotherapy has been developed rapidly in recent years due to the deepened acknowledgement of the immune system evasion by tumors. Developments, such as immune checkpoint inhibitors and receptors, have been made for the treatment of aggressive cancer. In order to observe the potential immune target in the CSC, we depicted the immune landscape of BC by clustering all patients into three clusters of different immune infiltration degrees. Our immune landscape revealed that T cell mediated immune response and innate immune response may lead the anti-tumor effect, further results from TIMER have confirmed these findings and provided a more specific profile of immune include dendritic cell and B lymphocyte. This provided a prospect and reminder in the direction of CSC-target immunotherapy since there have been successful attempts of generating CSC-primed T cells in vitro that showed targeting of CSCs after adoptive transfer in vivo (21). It has also been demonstrated that significant anti-CSC immunity was induced by dendritic cell vaccine basing on CSC both in vitro and in immune-competent hosts (8, 22). B lymphocytes were recently found to play a pivotal role in the cancer immune by being an assistant for T lymphocyte (23).
We also testified the reliability of the immune landscape by applying ESTIMATE analysis, as shown by someone, the relevance of breast CSC markers may vary according to the heterogeneity of the TME (24), where contains multiple types of cells, such as fibroblasts, endothelial cells, and mesenchymal stem cells, they are thought to have an affection to breast CSCs survival by secreting signaling molecules. Our results showed that patients in the high-infiltration cluster had the highest immune cell, highest stromal cell proportion, and lowest tumor purity compared with those in patients in the low- and intermediate-infiltration clusters, which is in agreement with previous studies. The change the stromal score from stage II to stage III and the estimate score from T2 to T3 also showed inversion, which synchronized with the stemness character. We found samples in stage IIIA tend to have highest stromal scores (N2), higher estimate scores (T3) and lower mRNAsi while samples in stage II tended to have lower stromal scores (N0-1), lower estimate scores (T2) and higher mRNAsi. These results indicated that the immune-like activity enhanced and likely validly regulated the stemness of BC cells in stage IIIA compared to stage II. The relationship between immune characteristics and stemness needs to be explored and the genes related to both stemness and immune properties are promising research targets.
As Pan Q et al. stated, identification of specific antigens or genetic alterations in CSCs may provide more specific targets for immunotherapy (25). In order to uncover the epigenetic regulations of BCSC have utilized to survive under immune surveillance, the stemness-related genes were assessed under different immune infiltration clusters, 13 stemness-related genes were found to expressed differently, expression level of these genes were significantly upregulated in tumor samples and associated with poorer overall survival and progression in BC. Among these, RAD51 is a protein encoding gene that is recruited to the perturbed replication DSBs and forks sites and respectively blocks the exonuclease activity of MRE11 on DSB repair and on the replicated genome and eventually limit self-DNA accumulation in the cytosol, this process also prevents the initiation of innate immune signaling mediated by STING (26), which brought an insight to the potential relationship between stemness and immune escape mechanisms in the tumor cell.
In order to observe the relationship between the expression of key genes and clinical characteristics of patients, we applied qRT-PCR to evaluate the expression of PIMREG, MTFR2, and CEP55 (Figure 14). Results have shown that these genes were highly expressed in BC tissue, these results indicated that the selected key genes may become therapeutic targets of BC, to be noted, PIMREG, MTFR2 were shown to have higher expression in TNM stage II compared with stage III, this indicate that the selection of preferable timing is rather crucial, for both the stemness of the tumor cells and immune infiltration characters of tumor tissue are in the phases of competition. PIMREG is commonly known as a cell cycle promoter of hypoxic fetal cardiomyocytes (27). It is known as a key gene involved in exogenous infections include suppurative periapical periodontitis (28) and paragonimiasis (29). It has been discovered that PIMREG promoted the aggressiveness of BC by disrupting the NF-κB/IκBα feedback loop (30). To be noted, NF-κB activation was proven to be involved in the tumor associated macrophages-mediated tumor growth in human pancreatic ductal adenocarcinoma (31). PIMREG is known to positively regulate STAT3 activity to promote cell differentiation and shown to be associated with poor survival in the BC and pancreatic cancer (31–34). It has been proven that both anti- and proinflammatory cytokines signal through associated receptor/JAK complexes and result in the phosphorylation of STAT3 (35). CEP55 protein is required for membrane fusion of cytokinesis through the inhibitory of cyclin-dependent kinase 1 phosphorylation (36, 37), CEP55 was found to inhibit apoptosis of human glioma cells via the PI3K/Akt signaling pathway, the PI3K/Akt pathway is a crucial pathway in the immune escape, it promotes the progression of BC growth by suppressing NK cell cytotoxicity through eIF2b (38), also, the inhibition of PI3K in lung cancer downregulated PD-L1 expression (39). MTFR2 is known as a protein that belongs to the MTFR family, it has been found that the MTFR2 was upregulated in BC and associated with poor survival of BC patients (40). MTFR2 alters glucose metabolism through activating HIF1α and HIF2α (EPAS1) in the BC cell lines (40, 41), the stability of which is also controlled by PI3K/AKT signaling axis (42). MTFR2 was also identified as an activator of TTK promoter in glioma stem-like cells (43). It is also reported to play a role in mitochondrial aerobic respiration essentially (43) and the promotion of mitochondrial fission (41). One of the strengths of the immune-therapies is the ability to target multiple antigens, which makes these approaches perfectly suited for targeting heterogenous CSC populations (44). PIMREG and MTFR2 are considered as potential target genes in CSC population that are required for effective immune targeting, the pan-cancer analysis implies that these genes may also play important roles in other tumors.
Figure 14
Conclusion
We identified 13 key genes related to stemness and immune escape and the stemness and immune escape degree of BC were increased non-linearly in the early stage of BC progression. PIMREG and MTFR2 are considered as effective diagnostic markers and potential targets for therapy. The development of stemness is involved in the regulation of both innate- and acquired-immune microenvironment.
Statements
Data availability statement
The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation.
Ethics statement
The studies involving human participants were reviewed and approved by The First Affiliated Hospital of Chongqing Medical University. The patients/participants provided their written informed consent to participate in this study. 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
Conception and design: JXT and LW. Acquisition of data: JXT and YW. Analysis and interpretation of: LW. Writing, review and/or revision of manuscript: LW, JL, and JJT. Administrative, technical or material support: XY. Study supervision: LW and JL. All authors contributed to the article and approved the submitted version.
Acknowledgments
This manuscript has been released as a pre-print at [https://www.researchsquare.com/article/rs-21080/v1], ([Linbang Wang] et al.) (43).
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.
Abbreviations
BC, Breast cancer; CSCs, Cancer stem cells; DEGs, Differentially expressed genes; TCGA, The Cancer Genome Atlas; GTEx, Genotype-Tissue Expression; K-M, Kaplan-Meier; GS, Gene significance; MM, Module membership; ssGSEA, Single-Sample Gene Set Enrichment Analysis; TME, Tumor microenvironment; PPI, Protein-protein interaction; GO, Gene ontology; KEGG, Kyoto Encyclopedia of Genes and Genomes; TNM, Tumor Node Metastasis.
References
1
DeSantisCEMaJGoding SauerANewmanLAJemalA. Breast cancer statistics, 2017, racial disparity in mortality by state. CA: Cancer J Clin (2017) 67(6):439–48. doi: 10.3322/caac.21412
2
LiuZZhangXSZhangS. Breast tumor subgroups reveal diverse clinical prognostic power. Sci Rep (2014) 4:4002. doi: 10.1038/srep04002
3
StinglJCaldasC. Molecular heterogeneity of breast carcinomas and the cancer stem cell hypothesis. Nat Rev Cancer (2007) 7(10):791–9. doi: 10.1038/nrc2212
4
BianchiniGBalkoJMMayerIASandersMEGianniL. Triple-negative breast cancer: challenges and opportunities of a heterogeneous disease. Nat Rev Clin Oncol (2016) 13(11):674–90. doi: 10.1038/nrclinonc.2016.66
5
BaiXNiJBeretovJGrahamPLiY. Cancer stem cell in breast cancer therapeutic resistance. Cancer Treat Rev (2018) 69:152–63. doi: 10.1016/j.ctrv.2018.07.004
6
SampieriKFoddeR. Cancer stem cells and metastasis. Semin Cancer Biol (2012) 22(3):187–93. doi: 10.1016/j.semcancer.2012.03.002
7
LiuSWichaMS. Targeting breast cancer stem cells. J Clin Oncol (2010) 28(25):4006–12. doi: 10.1200/jco.2009.27.5388
8
PellegattaSPolianiPLCornoDMenghiFGhielmettiFSuarez-MerinoBet al. Neurospheres enriched in cancer stem-like cells are highly effective in eliciting a dendritic cell-mediated immune response against malignant gliomas. Cancer Res (2006) 66(21):10247–52. doi: 10.1158/0008-5472.can-06-2048
9
Dianat-MoghadamHRokniMMarofiFPanahiYYousefiM. Natural killer cell-based immunotherapy: From transplantation toward targeting cancer stem cells. J Cell Physiol (2018) 234(1):259–73. doi: 10.1002/jcp.26878
10
HsuJMXiaWHsuYHChanLCYuWHChaJHet al. STT3-dependent PD-L1 accumulation on cancer stem cells promotes immune evasion. Nat Commun (2018) 9(1):1908. doi: 10.1038/s41467-018-04313-6
11
BindeaGMlecnikBTosoliniMKirilovskyAWaldnerMObenaufACet al. Spatiotemporal dynamics of intratumoral immune cells reveal the immune landscape in human cancer. Immunity (2013) 39(4):782–95. doi: 10.1016/j.immuni.2013.10.003
12
BaoXShiRZhangKXinSLiXZhaoYet al. Immune Landscape of Invasive Ductal Carcinoma Tumor Microenvironment Identifies a Prognostic and Immunotherapeutically Relevant Gene Signature. Front Oncol (2019) 9:903:903. doi: 10.3389/fonc.2019.00903
13
LiuZLiMJiangZWangX. A Comprehensive Immunologic Portrait of Triple-Negative Breast Cancer. Trans Oncol (2018) 11(2):311–29. doi: 10.1016/j.tranon.2018.01.011
14
JiaQWuWWangYAlexanderPBSunCGongZet al. Local mutational diversity drives intratumoral immune heterogeneity in non-small cell lung cancer. Nat Commun (2018) 9(1):5361. doi: 10.1038/s41467-018-07767-w
15
CharoentongPFinotelloFAngelovaMMayerCEfremovaMRiederDet al. Pan-cancer Immunogenomic Analyses Reveal Genotype-Immunophenotype Relationships and Predictors of Response to Checkpoint Blockade. Cell Rep (2017) 18(1):248–62. doi: 10.1016/j.celrep.2016.12.019
16
FedeleMCerchiaLChiappettaG. The Epithelial-to-Mesenchymal Transition in Breast Cancer: Focus on Basal-Like Carcinomas. Cancers (2017) 9(10):134–53. doi: 10.3390/cancers9100134
17
LiuXJohnsonSLiuSKanojiaDYueWSinghUPet al. Nonlinear growth kinetics of breast cancer stem cells: implications for cancer stem cell targeted therapy. Sci Rep (2013) 3:2473. doi: 10.1038/srep02473
18
VitaleIManicGDe MariaRKroemerGGalluzziL. DNA Damage in Stem Cells. Mol Cell (2017) 66(3):306–19. doi: 10.1016/j.molcel.2017.04.006
19
SteinbichlerTBDudasJSkvortsovSGanswindtURiechelmannHSkvortsovaII. Therapy resistance mediated by cancer stem cells. Semin Cancer Biol (2018) 53:156–67. doi: 10.1016/j.semcancer.2018.11.006
20
ChoiYJParkJHHanJWKimEJae-WookOLeeSYet al. Differential Cytotoxic Potential of Silver Nanoparticles in Human Ovarian Cancer Cells and Ovarian Cancer Stem Cells. Int J Mol Sci (2016) 17(12):2077–94. doi: 10.3390/ijms17122077
21
RappCWartaRStamovaSNowrouziAGeisenbergerCGalZet al. Identification of T cell target antigens in glioblastoma stem-like cells using an integrated proteomics-based approach in patient specimens. Acta Neuropathologica (2017) 134(2):297–316. doi: 10.1007/s00401-017-1702-1
22
Arias-PulidoHCimino-MathewsAChaherNQuallsCJosteNColpaertCet al. The combined presence of CD20 + B cells and PD-L1 + tumor-infiltrating lymphocytes in inflammatory breast cancer is prognostic of improved patient outcome. Breast Cancer Res Treat (2018) 171(2):273–82. doi: 10.1007/s10549-018-4834-7
23
WangHYangMLinLRenHLinCLinSet al. HepG2 cells acquire stem cell-like characteristics after immune cell stimulation. Cell Oncol (Dordrecht) (2016) 39(1):35–45. doi: 10.1007/s13402-015-0249-1
24
DirkseAGolebiewskaABuderTNazarovPVMullerAPoovathingalSet al. Stem cell-associated heterogeneity in Glioblastoma results from intrinsic tumor plasticity shaped by the microenvironment. Nat Commun (2019) 10(1):1787. doi: 10.1038/s41467-019-09853-z
25
PanQLiQLiuSNingNZhangXXuYet al. Concise Review: Targeting Cancer Stem Cells Using Immunologic Approaches. Stem Cells (Dayton Ohio) (2015) 33(7):2085–92. doi: 10.1002/stem.2039
26
MukherjeeSAbdisalaamSBhattacharyaSSrinivasanKSinhaDAsaithambyA. Mechanistic link between DNA damage sensing, repairing and signaling factors and immune signaling. Adv Protein Chem Struct Biol (2019) 115:297–324. doi: 10.1016/bs.apcsb.2018.11.004
27
HashimotoKKodamaAHondaTHanashimaAUjiharaYMurayamaTet al. Fam64a is a novel cell cycle promoter of hypoxic fetal cardiomyocytes in mice. Sci Rep (2017) 7(1):4486. doi: 10.1038/s41598-017-04823-1
28
SousaELMartinhoFCNascimentoGGLeiteFRGomesBP. Quantification of endotoxins in infected root canals and acute apical abscess exudates: monitoring the effectiveness of root canal procedures in the reduction of endotoxins. J Endodontics (2014) 40(2):177–81. doi: 10.1016/j.joen.2013.10.008
29
ChengNXuXNZhouYDongYTBaoYFXuBet al. Cs1, a Clonorchis sinensis-derived serodiagnostic antigen containing tandem repeats and a signal peptide. PloS Neglected Trop Dis (2018) 12(8):e0006683. doi: 10.1371/journal.pntd.0006683
30
JiangLRenLZhangXChenHChenXLinCet al. Overexpression of PIMREG promotes breast cancer aggressiveness via constitutive activation of NF-kappaB signaling. EBioMedicine (2019) 43:188–200. doi: 10.1016/j.ebiom.2019.04.001
31
YeHZhouQZhengSLiGLinQWeiLet al. Tumor-associated macrophages promote progression and the Warburg effect via CCL18/NF-kB/VCAM-1 pathway in pancreatic ductal adenocarcinoma. Cell Death Dis (2018) 9(5):453. doi: 10.1038/s41419-018-0486-0
32
YaoZZhengXLuSHeZMiaoYHuangHet al. Knockdown of FAM64A suppresses proliferation and migration of breast cancer cells. Breast Cancer (Tokyo Japan) (2019) 26(6):835–45. doi: 10.1007/s12282-019-00991-2
33
JiaoYFuZLiYZhangWLiuY. Aberrant FAM64A mRNA expression is an independent predictor of poor survival in pancreatic cancer. PloS One (2019) 14(1):e0211291. doi: 10.1371/journal.pone.0211291
34
LiZZhangYZhangZZhaoZLvQ. A four-gene signature predicts the efficacy of paclitaxel-based neoadjuvant therapy in human epidermal growth factor receptor 2-negative breast cancer. J Cell Biochem (2019) 120(4):6046–56. doi: 10.1002/jcb.27891
35
OwenKLBrockwellNKParkerBS. JAK-STAT Signaling: A Double-Edged Sword of Immune Regulation and Cancer Progression. Cancers (2019) 11(12):2002–28. doi: 10.3390/cancers11122002
36
ZhouCHancockJLKhannaKKHomerHA. First meiotic anaphase requires Cep55-dependent inhibitory cyclin-dependent kinase 1 phosphorylation. J Cell Sci (2019) 132(18):132–43. doi: 10.1242/jcs.233379
37
MondalGRowleyMGuidugliLWuJPankratzVSCouchFJ. BRCA2 localization to the midbody by filamin A regulates cep55 signaling and completion of cytokinesis. Dev Cell (2012) 23(1):137–52. doi: 10.1016/j.devcel.2012.05.008
38
JinFWuZHuXZhangJGaoZHanXet al. The PI3K/Akt/GSK-3beta/ROS/eIF2B pathway promotes breast cancer growth and metastasis via suppression of NK cell cytotoxicity and tumor cell susceptibility. Cancer Biol Med (2019) 16(1):38–54. doi: 10.20892/j.issn.2095-3941.2018.0253
39
GaoYYangJCaiYFuSZhangNFuXet al. IFN-gamma-mediated inhibition of lung cancer correlates with PD-L1 expression and is regulated by PI3K-AKT signaling. Int J Cancer (2018) 143(4):931–43. doi: 10.1002/ijc.31357
40
LuGLaiYWangTLinWLuJMaYet al. Mitochondrial fission regulator 2 (MTFR2) promotes growth, migration, invasion and tumour progression in breast cancer cells. Aging (2019) 11(22):10203–19. doi: 10.18632/aging.102442
41
MonticoneMPanfoliIRaveraSPuglisiRJiangMMMorelloRet al. The nuclear genes Mtfr1 and Dufd1 regulate mitochondrial dynamic and cellular respiration. J Cell Physiol (2010) 225(3):767–76. doi: 10.1002/jcp.22279
42
DeyNSunYCarlsonJHWuHLinXLeyland-JonesBet al. Anti-tumor efficacy of BEZ235 is complemented by its anti-angiogenic effects via downregulation of PI3K-mTOR-HIF1alpha signaling in HER2-defined breast cancers. Am J Cancer Res (2016) 6(4):714–46.
43
WangJXieYBaiXWangNYuHDengZet al. Targeting dual specificity protein kinase TTK attenuates tumorigenesis of glioblastoma. Oncotarget (2018) 9(3):3081–8. doi: 10.18632/oncotarget.23152
44
EltoukhyHSSinhaGMooreCASandifordOARameshwarP. Immune modulation by a cellular network of mesenchymal stem cells and breast cancer cell subsets: Implication for cancer therapy. Cell Immunol (2018) 326:33–41. doi: 10.1016/j.cellimm.2017.07.011
Summary
Keywords
breast cancer, cancer stem cell, tumor immune infiltration, PIMREG, MTFR2
Citation
Wang L, Liu W, Liu J, Wang Y, Tai J, Yin X and Tan J (2021) Identification of Immune-Related Therapeutically Relevant Biomarkers in Breast Cancer and Breast Cancer Stem Cells by Transcriptome-Wide Analysis: A Clinical Prospective Study. Front. Oncol. 10:554138. doi: 10.3389/fonc.2020.554138
Received
21 April 2020
Accepted
31 December 2020
Published
24 February 2021
Volume
10 - 2020
Edited by
Sofia Diana Merajver, University of Michigan, United States
Reviewed by
Zhifang Zhang, Beckman Research Institute, City of Hope, United States; Giovanni Porta, University of Insubria, Italy
Updates
Copyright
© 2021 Wang, Liu, Liu, Wang, Tai, Yin and Tan.
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: Jinxiang Tan, tjx1202@163.com
This article was submitted to Cancer Molecular Targets and Therapeutics, a section of the journal Frontiers in Oncology
Disclaimer
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.