A prognosis marker MUC1 correlates with metabolism and drug resistance in bladder cancer: a bioinformatics research

Background MUC1 is a type I transmembrane protein that plays an important role in tumor cell signal transduction. Although current studies have shown that MUC1 is upregulated in bladder cancer (BC), the specific mechanism is still unclear. Methods We performed expression analysis, gene set enrichment analysis, survival analysis, immune infiltration analysis, drug sensitivity analysis, and metabolism-related gene expression analysis on TCGA-BLCA, GES31684 and GSE13507. Results The expression of MUC1 in the tumor and lymphatic metastasis positive samples was significantly increased. Genes related to MUC1 expression were significantly enriched in immune response, ribosomes, exosomes, and energy metabolism. The results of the immune infiltration analysis showed that M1 macrophages in BC with high MUC1 expression were significantly decreased. Expression of MUC1 increases drug resistance in BC patients. In addition, MUC1 increases glycolysis, glucose uptake, and lactate production by inducing metabolic reprogramming. Conclusion MUC1 has a significant effect on the metabolism and immune cell infiltration of BC, which may be the cause of increased drug resistance, and can be used as a molecular target for the diagnosis and treatment of BC.


Introduction
Bladder cancer (BC) is a malignant tumor that seriously endangers human health. There are approximately 573,000 new cases of BC each year worldwide, and the number of deaths caused by BC is approximately 213,000 ranking tenth among all tumors [1]. BC is divided into muscle invasive bladder cancer (MIBC) and non-muscle invasive bladder cancer (NMIBC) according to the depth of invasion [2]. BC mainly includes transitional cell carcinoma, squamous cell carcinoma, and adenocarcinoma on the basis of the pathological type, among which transitional cell carcinoma accounts for more than 90% [3]. There are plenty of risk factors in the development of BC, including gender, smoking prevalence, occupational carcinogen exposure, Schistosoma haematobium infection, etc. [4]. Recently, researchers have considered the interaction of hypoxia and tumor, immune cell infiltration and immune microenvironment as the key step in the occurrence and development of BC [5,6]. Currently, early diagnosis of BC can rely on serum markers, urine markers, and cystoscopy. However, the lack of clear markers and poor patients' compliance make the early diagnosis of BC very difficult [7]. Despite the 5-year survival rate for BC is as high as 77%, the recurrence rate of BC is still high, and the 5-year survival rate for metastatic BC is less than 5% [8]. Therefore, it is particularly urgent to find highly sensitive biomarkers for BC.
MUC1 encodes a type I transmembrane protein with a polar distribution, called Mucin 1. N-terminal alpha subunit of Mucin1 functions in cell adhesion and the C-terminal beta subunit is involved in cell signaling [9]. Recent studies have found that MUC1 plays an important role in tumor hypoxia microenvironment and tumor metabolism [10]. MUC1 can not only mediate hypoxia-driven angiogenesis by regulating a variety of angiogenic factors; it can also stabilize and activate hypoxia inducible factor-1α (HIF-1α) to promote the metabolic reprogramming of cancer cells [11,12]. Recently, studies have developed antibody-drug conjugated (ADC) based on MUC1 employed in the treatment trastuzumab-resistant breast cancer patients [13]. Besides, a multitargeted recombinant Ad5 PSA/ MUC-1/brachyury-based immunotherapy vaccine developed by Bilusic et al. was recently introduced into phase I trials [14]. The types and abundance of immune cell infiltration in the tumor microenvironment have proved to be very different from normal tissues [14]. Studies have found that the expression of MUC1 induces immune suppression in colon cancer, and this suppression can be reversed by blocking the PD1/PDL1 pathway [15]. MUC1 has shown a strong potential in the diagnosis and treatment of tumors.
Although some studies claim that MUC1 has significance in BC, the specific mechanism remains unclear. In this study, we used gene expression data obtained from the the Cancer Genome Atlas (TCGA) and the Gene Expression Omnibus (GEO) database to explore the expression of MUC1 in human BC samples. We used R (version 4.0.3) to analyze the correlation between MUC1 expression and clinical characteristics. To better understand the mechanism of MUC1 in the occurrence and development of BC, we conducted Gene Set Enrichment Analysis (GSEA) and Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) analysis [16]. Next, we used Tumor Immune Estimation Resource (TIMER) and CIBERSORT to analyze the relationship between MUC1 expression and tumor immune cell infiltration. Cancer Cell Line Encyclopedia (CCLE) and Genomics of Drug Sensitivity in Cancer (GDSC) was used to analyze the relationship between MUC1 expression and drug sensitivity. Finally, we use Gene Expression Profile Interaction Analysis (GEPIA), Kaplan-Meier survival analysis (KM), and Human Protein Atlas (HPA) to analyze the relationship between MUC1 and patient prognosis.

Data retrieval and download
To explore the expression of MUC1 in BC and its relationship with clinical characteristics, gene expression data and clinical data of BC in TCGA (TCGA-BLCA cohort) were downloaded from the Genomic Data Commons Data Portal (https:// portal. gdc. cancer. gov). Besides, we downloaded the gene expression data and survival data of GSE31684 and GSE13507 as a validation dataset from the GEO (https:// www. ncbi. nlm. nih. gov/ geo) [17,18].

Gene expression analysis
We compared the expression of MUC1 in cancer and adjacent normal tissues, as well as in positive and negative lymphatic metastasis tissues. With the GSE13507 dataset, we compared MUC1 expression in normal, paraneoplastic, primary and recurrent tumour tissues. In addition, we used GEPIA to show the expression of MUC1 in bladder tumor tissues and adjacent tissues (http:// gepia. cancer-pku. cn/) [19]. Next, we removed samples with incomplete clinical data (age, gender, American Joint Committee on Cancer (AJCC) pathologic stage, AJCC pathologic N, AJCC pathologic M, AJCC pathologic T). The protein expression of MUC1 between normal and cancer tissues was retrieved and compared from the HPA database (www. prote inatl as. org).

Prognostic analysis
Univariate and multivariate analyses were conducted using the Cox Proportional Hazards Regression model on the remaining samples to obtain risk scores. To further understand the relationship between MUC1 expression and clinical characteristics, we conducted receiver operating characteristic (ROC) analysis and logistic regression analysis. After dividing the samples into high and low expression group according to the median expression of MUC1, KM survival analysis was conducted to analyze the overall survival (OS) rate of BC patients related to MUC1 expression. Furthermore, we conducted a similar analysis in GSE31684 and GES13507 to verify the results obtained previously.

Gene set enrichment analysis
To explore the mechanism of MUC1 in BC, we divided the samples into MUC1-H and MUC1-L groups according to the level of MUC1 expression. Next, the gene differential expression analysis was performed on two sets of samples, and the genes that were differentially expressed were obtained (logFC > 2, p value < 0.05). Incorporate differentially expressed genes in GSEA analysis. The number of permutations was set to 1000. To investigate the possible biological functions of MUC1, we used GSEA to analyze the GO pathway and the KEGG pathway. The absolute value of normalized enrichment score (|NES|) is greater than 2, nominal p value (NOM p value) is less than 0.05, and false discovery rate (FDR) is less than 0.05 of all enrichment pathways.

Immune infiltrate analysis
Li et al. developed a deconvolution algorithm (TIMER) to calculate the type and abundance of immune cell infiltration based on the gene expression profiles (http:// timer. cistr ome. org/) [20]. We evaluated MUC1 expression in BC and its correlation with the type and abundance of immune cell infiltration, including B cells, CD4 + T cells, CD8 + T cells, Dendritic cells, Neutrophils and Macrophages. Another deconvolution algorithm based on gene expression called CIBERSORT (http:// ciber sort. stanf ord. edu/) was used to further assess the relationship between MUC1 and immune cell [21,22]. And we performed a nonparametric Wilcox test to assess whether the infiltration of immune cells differed between MUC1-H and MUC-L groups. In addition, we also assessed the correlation between different immune cells.

Metabolism analysis
The samples in the TCGA, BLCA, and GSE13507 datasets were divided into two groups according to the expression of MUC1, and the expression levels of genes related to glucose metabolism (HK2, GLUT1, PKM, RWDD3, SLC16A3, SLC5A12, ENO1) in different groups were compared. Data were analyzed using t-test and p value < 0.05 was considered statistically significant.

Expression of MUC1 at the mRNA and protein levels and its relationship to clinical characteristics
We compared MUC1 mRNA expression between BC and normal tissues in TCGA database, and also between lymphphatic metastasis negative and positive samples. A total of 411 tumor tissues and 19 normal tissues were enrolled in the expression analysis, including 233 lymphatic metastasis negative and 81 lymphatic metastasis positive samples. The expression of MUC1 in different groups was represented by Scatterplot (Fig. 1A). The results showed that the expression of MUC1 in cancer tissues was significantly higher than that in normal tissues (p value = 2.2e−2), and it was significantly higher in lymph node metastasis positive samples than in lymph node metastasis negative (p value = 4.2e−5). The GSE13507 dataset contains normal, paraneoplastic, primary and recurrent tumour samples. We analysed the expression of MUC1 in different subgroups of GSE13507. The results showed significant differences in MUC1 between primary tumours and paraneoplastic tissues (p value < 0.002). Similarly, MUC1 was significantly different between primary tumours and normal tissues (p value = 0.0298) (Fig. 1B). Furthermore, we used GEPIA to compare the expression of MUC1 in normal and tumor tissues and found that the mRNA expression of MUC1 significantly increased in tumor tissues (p value < 0.01) (Fig. 1C).
We evaluated the connection between clinical characteristics and MUC1 expression levels of BC patients. There are 382, 382, 192, 350 samples with complete data of AJCC pathologic Stage, AJCC pathologic T stage, AJCC pathologic M stage, and AJCC pathologic N stage, respectively. The results showed that MUC1 has little correlation with AJCC pathological T stage (p value > 0.05, Fig. 1D) and AJCC pathological Stage (p value > 0.05, Fig. 1F), but is significantly correlated with AJCC pathological N stage (p value < 0.05, Fig. 1E) and AJCC pathological M stage (p value < 0.05, Fig. 1G). In addition, we conducted logistic regression to analyze the correlation between MUC1 expression and clinical characteristics ( Table 1). The increase of MUC1 expression level in BC was significantly correlated with lymphatic metastasis (N1 vs. N0, p value = 0.007). These indicate that compared with BC patients with low levels of MUC1 expression, BC patients with high levels of MUC1 expression are more likely to develop lymphatic metastasis. Next, to understand the difference in protein levels of MUC1, we searched HPA for the immunohistochemical data of MUC1 in BC. The results showed that MUC1 has higher levels of expression in tumor compared to normal tissues (Fig. 1H).

The role of MUC1 in the prognosis of bladder cancer
We used univariate Cox proportional hazards regression analysis to explore the relationship between MUC1 expression, clinical characteristics, and OS. The results showed that Age (HR  (Table 2). Next, multivariate Cox proportional hazards regression analysis showed that MUC1 expression is an independent prognostic factor ( Fig. 2A, Table 3).
The distribution of MUC1 expression, survival status of patients with BC, and expression profiles of MUC1 are shown in Fig. 2B. ROC curve analysis showed that the AUC values of 1, 3, and 5 years were 0.6, 0.634, and 0.634, respectively, indicating that the expression of MUC1 has a potential prognostic ability (Fig. 2C).
Next, we performed a survival analysis in TCGA based on the expression of MUC1. The results of TCGA show that MUC1 expression is significantly related to OS (p value = 1.47e−03) (Fig. 2D). To further verify the influence of MUC1 on the OS of BC patients, we further conducted similar analysis in GSE13507 and GSE31684. As expected, similar results appeared in GSE31684 (p value = 1.21e−03) and GSE13507 (p value = 2.66e−02) (Fig. 2E, F).

MUC1 related enrichment pathway
To explore the potential biological functions of MUC1, we performed GO term and KEGG pathway enrichment analyses of differentially expressed genes. Signal pathways that meet the following conditions (|NES| > 2, NOM p value < 0.05, FDR < 0.05) are considered enriched. As shown in Table 4, there are five positive GO Biological Process terms related to MUC1 (nuclear-transcribed mRNA catabolic process, nonsense-mediated decay; cotranslational protein targeting to the membrane; SRP-dependent cotranslational protein targeting to the membrane; protein targeting to ER; establishment of protein localization to endoplasmic reticulum) and five negative GO Biological Process terms related to MUC1 (cornification; keratinization; immune responseregulating signaling pathway; immune response-regulating cell surface receptor signaling pathway; defense response to other organism) (Fig. 3A). There are five positive GO Cellular Component terms related to MUC1 (cytosolic large ribosomal subunit; cytosolic ribosome;   (Fig. 3C). As shown in Table 5, there are two positive KEGG pathways related to MUC1 (Coronavirus disease, COVID-19 and Ribosome) and two negative KEGG pathways related to MUC1 (Staphylococcus aureus infection and Intestinal immune network for IgA production) (Fig. 3D). These results indicate that the pathways that regulate ribosomal and lipid metabolism are critical to BC patients and are closely related to the expression of MUC1.

Relationship between MUC1 expression and immune cell infiltration
We used TIMER to clarify the correlation between the expression of MUC1 and the level of immune infiltration of BC. The results have shown that MUC1 expression is significantly correlated with B cells (p value = 1.50e−03), macrophages (p value = 3.593e−02) and Neutrophils (p value = 3.94e−04), indicating that MUC1 plays a key role in the immune infiltration of BC (Fig. 4A). To explore the impact of MUC1 expression on the immune microenvironment of BC patients, we divided the samples into MUC1-H and MUC1-L according to MUC1 expression. CIBERSORT was used to calculate the infiltration abundance of 22 immune cells in each sample to assess the difference between different groups. The results indicate that M1 macrophages are affected by the expression of MUC1 (p value = 0.033) (Fig. 4B). We also assessed the potential correlation between immune cell infiltration and clinical features. The results showed a strong correlation between immune cell infiltration and clinical features (Fig. 4C).

Discussion
Recently, the role of MUC1 in tumors has gradually been emphasized. In pancreatic cancer, MUC1 reduces radiation-induced pancreatic cancer cell toxicity and DNA damage by enhancing glycolysis, the pentose phosphate pathway, and nucleotide biosynthesis [23]. In triple-negative breast cancer, MUC1-C increases PD-L1 transcription and promotes tumor immune escape by recruiting MYC and NF-κB p65 to the PD-L1 promoter region [24]. In prostate cancer, MUC1-C inhibits androgen receptor (AR) axis signaling and induces the expression of the neural BRN2 transcription factor, increasing the self-renewal capacity and tumorigenicity of prostate cancer cells [25].
In addition, MUC1 also plays an important role in ovarian cancer [26], non-small cell lung cancer [27] and pancreatic cancer [28]. One study found that the positive rates of MUC1 and MUC2 in urothelial carcinoma were 89.7% and 44.3%, respectively, and MUC1 expression was significantly correlated with tumor grade, while MUC2 was the opposite. Different staining patterns of MUC1 and MUC2 have important prognostic implications in patients with difficult-to-recognize low-grade papillary urothelial lesions and precancerous lesions [29]. Urothelial carcinoma patients expressing high levels of MUC1-C have poor survival and are strongly associated with cisplatin resistance [30]. However, the specific mechanism of MUC1 in BC remains unclear. We analyzed the potential mechanism of MUC1 in BC (Fig. 5D). In our current study, to explore the possibility of MUC1 as a prognostic biomarker for BC, we used TCGA BLCA cohort to evaluate the prognostic value of MUC1. First, we analyzed the expression of MUC1 (protein level and mRNA level), and the results showed that the expression of MUC1 in tumor tissues was significantly higher than that in adjacent normal tissues, and the expression of MUC1 in lymphatic metastasis positive samples was significantly higher than that in lymphatic metastasis negative samples. BC patients with high expression of MUC1 showed advanced N stage, M stage, and tumor status. We analyzed the role of MUC1 and clinical features in prognosis, and the results showed that high MUC1 expression was strongly associated with worse prognosis. These results suggest that MUC1 may serve as an independent prognostic factor for overall survival in BC.
MUC1 is closely related to metabolism and resulting drug resistance. Recent studies have found that there is a crosstalk between the MUC1 and HIF-1 signal pathways. The HIF-1 signaling pathway has been shown to be abnormally activated in numerous cancers, especially solid tumors [31]. HIF-1 consists of two subunits, HIF-1α and HIF-1β, of which the α subunit is the active unit. In oxygenated cells, HIF-1α will be hydroxylated by proline hydroxylase (PHD1), which will then be recognized and degraded by the ubiquitin-protease complex; while in hypoxic cells, the activity of PHD1 decreases and the lack of a substrate for the hydroxylation reaction leads to the accumulation of HIF-1α [32]. The transcriptional targets of HIF-1α include genes that play key roles in angiogenesis, erythropoiesis, metabolic reprogramming, vasomotor function, and apoptosis proliferation response [33]. In the hypoxic condition, MUC1 promotes the recruitment of HIF-1α and p300 to the promoter region of glycolysis genes, upregulates the expression of glycolysis-related genes, increases glucose uptake and lactate production, and assists tumor cells to survive in hypoxic environments [12]. In addition, MUC1 can increase the carbon flux of the non-oxidative pentose phosphate pathway (PPP) and pyrimidine synthesis pathway of pancreatic cancer cells by promoting the expression of HIF-1α, leading to gemcitabine resistance [34].
In this study, GLUT1, PKM, and SLC16A3 were significantly upregulated in the MUC1-H group, while RWDD3 was significantly downregulated. GLUT1 encodes the major glucose transporter in cell [35]. Pyruvate kinase is encoded by the PKM and is the main ratelimiting enzyme in the glycolysis process [36]. SLC16A3 removes lactate produced by glycolysis from the plasma membrane and plays a key role in the growth and proliferation of hypoxic cancer cells [37]. RWDD3 negatively regulates the HIF-1α signaling pathway by increasing the sumoylation of HIF-1α, promoting its stabilization, transcriptional activity, and the expression of its target gene VEGFA during hypoxia [38]. Therefore, we can draw the hypothesis: MUC1 increases the uptake of glucose and the production of lactate through the HIF-1α signaling pathway, so that tumor cells can survive under hypoxic conditions. In addition, GO and KEGG pathway enrichment analysis also showed that up-regulated MUC1 was mainly closely related to metabolism. Through the CCLE and GDSC databases, we found that the expression of MUC1 in BC was associated with drug resistance, which may be caused by increased pyrimidine metabolism, and the specific mechanism remains to be further studied. We used the TIMER database to explore the relationship between the expression of MUC1 and the level   [40]. Almost simultaneously, Maureen et al. successfully constructed a BCG-human interleukin-2 (hIL2)-MUC1 breast cancer vaccine, which significantly reduced the incidence and rate of tumorigenesis in mice, suggesting that co-expression of MUC1 and IL-2 plays a key role in the induction of tumor-specific immune responses [41]. Yuan et al. constructed a BCG-MUC1-based breast cancer vaccine (rBCG-MVNTR4-CD80), which significantly suppressed tumors, induced interferon production, and stimulated both CD4+ and CD8+ positive lymphocyte production [42]. In addition, the authors constructed two other rBCG-MUC1 vaccines (rBCG-MVNTR4-CSF and rBCG-MVNTR8-CSF, with similar effects [43]. In 2016, an anti-tumor vaccine based on MUC1 (MUC1maltose-binding protein (MBP)/BCG) entered preclinical trials and showed no significant toxicity in mice, rats or crab-eating monkeys.MUC1-MBP/BCG acts as an antitumor agent by inducing Th1 cell activation and MUC1specific IgG antibody secretion [44]. At present, most of the research on MUC1 and BCG is targeting MUC1 and building tumor vaccines based on BCG to better play the role of BCG. Therefore, high MUC1 expression can theoretically improve the responsiveness of tumors to such rBCG drugs, but there is no research on the relationship between traditional BCG and MUC1 expression. This needs further study.
MUC1 has shown great potential as a widely expressed membrane protein in the diagnosis and treatment of tumors. However, the mechanism of action of MUC1 in BC has not been reported in the literature. The role of MUC1 in BC was investigated by expression analysis, prognostic analysis, immunological analysis, metabolic analysis and drug sensitivity analysis. The potential mechanisms of MUC1-mediated drug resistance in bladder cancer were uncovered. However, there are some shortcomings in this study: no experimental study was conducted to verify the function of MUC1, and the role of MUC1 needs to be verified by cellular and animal experiments subsequently.

Conclusion
In short, this is the first report to explore the possibility of MUC1 as a marker of BC at the mechanistic level. This work helps to elucidate the role of cellular metabolism, immune cell infiltration, and genes regulating these processes in BC development and drug resistance. In addition, MUC1 can be used as a target for monoclonal antibodies, antibody-drug conjugate, cancer vaccine, Aptamer, and radiotherapy. Understanding the mechanism of MUC1 will help to develop corresponding treatment and diagnosis methods.