Characterization of the Lipid Metabolism in Bladder Cancer to Guide Clinical Therapy

Background Bladder cancer is one of the most common malignancies of the urinary system with an unfavorable prognosis. More and more studies have suggested that lipid metabolism could influence the progression and treatment of tumors. However, there are few studies exploring the relationship between lipid metabolism and bladder cancer. This study aimed to explore the roles that lipid metabolism-related genes play in patients with bladder cancer. Methods TCGA_BLCA cohort and GSE13507 cohort were included in this study, and transcriptional and somatic mutation profiles of 309 lipid metabolism-related genes were analyzed to discover the critical lipid metabolism-related genes in the incurrence and progression of bladder cancer. Furthermore, the TCGA_BLCA cohort was randomly divided into training set and validation set, and the GSE13507 cohort was served as an external independent validation set. We performed the LASSO regression and multivariate Cox regression in training set to develop a prognostic signature and further verified this signature in TCGA_BLCA validation set and GSE13507 external validation set. Finally, we systematically investigated the association between this signature and tumor microenvironment, drug response, and potential functions and then verified the differential expression status of signature genes in the protein level by immunohistochemistry. Results A novel 6-lipidmetabolism-related gene signature was identified and validated, and this risk score model could predict the prognosis of patients with bladder cancer. In addition, the prognostic model was tightly related to immune cell infiltration and tumor mutation burden. Gene set variation analysis (GSVA) and gene set enrichment analysis (GSEA) showed that mTOR signaling pathway, G2M checkpoint, fatty acid metabolism, and hypoxia were enriched in patients in the high-risk score groups. Furthermore, 3 therapies specific for bladder cancer patients in different risk scores were identified. Conclusion s. In conclusion, we investigated the lipid metabolism-related genes in bladder cancer through comprehensive bioinformatic analysis. A novel 6-gene signature associated with lipid metabolism for predicting the outcomes of patients with bladder cancer was conducted and validated. Furthermore, the risk score model could be utilized to indicate the choice of therapy in bladder cancer.


Introduction
Bladder cancer is one of the most common malignancies of the urinary system; it has the 13 th highest mortality among all cancers, and the mortality is still rising despite tremendous efforts that have been made for the treatment [1]. e progression of bladder cancer is a multistage process including environmental and genetic factors [2]. Previous studies indicate that tobacco smoking and occupational exposure to carcinogens are major factors [3]. e primary treatment for bladder cancer is transurethral resection of the bladder (TURB) accompanied by intravesical chemotherapy or immunotherapy [3]. However, the prognosis of patients still remains unfavorable. erefore, it is essential to identify prognostic biomarkers to guide the selection of treatment for improving the curative effects.
Accumulating evidence demonstrates that clinical outcome, epigenetic status, and treatment resistance are associated with tumor metabolism [4,5]. Metabolic reprogramming has been considered to be a new hallmark of malignant tumors [6]. Cancer cells usually require more energy to meet their biological needs than normal cells [7], while fatty acid oxidation is an important energy source for cancer cells, so lipid metabolism in cancer cells has been recognized as playing an important role in cancer progression [8].
Recent studies have shown that fatty acid metabolism has a close connection with bladder cancer [9]. Epidemiologic studies have shown that free fatty acid (FFA) level was increased in bladder cancer compared with paracancerous tissues. Cheng et al. demonstrated that inhibition of fatty metabolism is important in sustaining tumor growth through PPAR-c-mediated pathway [10]. Besides, the molecular mechanism of drug resistance of cancer therapy might include lipid metabolism reprogramming [11]. Rysman E et al. indicated that altered lipid composition of cellular membranes could disrupt the uptake of chemotherapeutic agents and lead to chemotherapeutic resistance [12]. What is more, the tumor microenvironment (TME) plays an important role in the progression of bladder cancer, while metabolic disorders including FA metabolism changes have a crucial impact on cancer [13]. Fatty acid secreted in the microenvironment could affect the function and phenotype of infiltrating immune cells [14].
In this study, we explored the lipid metabolism-related genes in bladder cancer to conduct a model to predict patient prognosis. e prognostic risk score model independently predicted the survival outcome of bladder cancer patients. What is more, the relationship between risk score model and TME cell-infiltrating characteristics was investigated and suitable therapy treatment could be selected through the risk score model. ese findings can provide a new sight into exploring the metabolic mechanism and treatment for bladder cancer.

Data Acquisition.
We included two datasets in this study, and they are TCGA_BLCA cohort and GSE13507 cohort. e transcriptome profiles, somatic mutation profiles, and clinical information about TCGA_BLCA were downloaded from the GDC_Portal (https://portal.gdc.cancer.gov/), and the transcriptome profiles and corresponding clinical information of GSE13507 cohort were downloaded from the GEO database (https://www.ncbi.nlm.nih.gov/gds/?term�GSE13507) and served as independent external validation sets. Notably, the batch effects between TCGA_BLCA and GSE13507 cohorts were normalized by the R package "sva."

Landscape of Lipid Metabolism-Related Genes in Bladder
Cancer. We firstly investigated the role of lipid metabolismrelated genes in TCGA_BLCA cohort, performed differential expression analysis between tumor and normal samples, and screened differentially expressed lipid metabolism-related genes with a threshold of FC (fold change, FC) > 1.5 and adjusted p value < 0.05. Following this, we performed univariate cox regression of these differentially expressed lipid metabolism-related genes in TCGA_BLCA cohort to further identify the critical differentially expressed lipid metabolism-related genes with significant prognostic value. Subsequently, we explored the somatic mutation of these critical genes and draw the mutation atlas.

Identification of the Prognostic Signature.
Having systematically summarized the role of lipid metabolism-related genes in TCGA_BLCA cohort, we would like to establish a prognostic signature by these critical genes. us, we first randomly split all the patients in TCGA_BLCA cohort with a ratio of 1 : 1 that one for the construction and the other for the verification. Moreover, the GSE13507 cohort was served as an external independent validation cohort. en, we separately conducted LASSO regression to screen appropriate variables and multivariate cox regression to establish the prognostic signature in the training set, and a formula of risk score based on these lipid metabolism-related genes was established: risk score � coef(i). exp(i). (1) i represents each gene in the prognostic signature, coef (i) represents the coefficient of this gene, and exp (i) is the expression value of it. us, each sample in training set, validation set, and GSE13507 acquired a risk score according to this formula. Moreover, we set the medium value of the risk score in training set as the threshold, and each patient received a risk level that the higher is high risk and the lower is low risk.

Further Verification of the Prognostic Signature.
We first conducted survival analysis in all three sets to confirm the prognostic value of this signature, then performed univariate Cox regression to calculate the hazard ratio of the risk score, and carried out a meta-analysis to summarize the HR of risk score in three different sets. Besides, the differences in clinical information between high-risk patients and low-risk patients were explored by the chi-square test. Also, the ROC curves of the risk score in three sets were plotted and the area under the curves was calculated. Furthermore, we combined the indicator of tumor mutation burden and our risk score to predict the survival of patients together and explored the difference in TMB between high-/low-risk patients and the correlation between TMB and risk score. Finally, we further investigated the mutation differences between high-risk patients and low-risk patients in all TCGA_BLCA patients.

Gene Set Enrichment Analysis (GSEA) and Gene Set
Variation Analysis (GSVA). Having constructed and verified the prognostic signature, we wonder about the further potential mechanisms behind the risk score.
us, we separately conducted the GSEA and GSVA in all TCGA_BLCA patients by the R package "clusterProfiler" and "GSVA." e gene ontology (GO) gene sets, KEGG gene sets, Hallmarks gene sets, metabolism-related gene sets, and cell death-related gene sets were used to do the corresponding analysis.
2.6. Tumor Microenvironment, Drug Response, and Immunohistochemistry. Notably, we carried out eight different algorithms to quantify the immune cells or immunerelated function of each patient according to its transcriptome profile, and they were XCELL, TIMER, QUANTISEQ, MCPCOUNTER, EPIC, CIBERSORT-ABS, CIBERSORT, and ssGSEA. en, both differential infiltration analysis and correlation tests were conducted. Following this, we separately estimate the drug response of each patient to the commonly used drugs including cisplatin, gemcitabine, and others small molecule drugs by R package "proPhetics." Also, the TIDE score was compared between high-risk patients and low-risk patients to predict the drug response to immunotherapy. Finally, we searched the immunohistochemistry profiles of these signature genes in the Human Protein Atlas database (HPA, https://www.proteinatlas.org/) and compared and verified the differential expression status in the protein expressed level.

Statistical Analysis.
All the data were processed or analyzed by the R program version 4.1.1 and Microsoft Office Excel. P < 0.05 was regarded with significant statistical differences.

Identification of 18 Vital Differentially Expressed Lipid
Metabolism-Related Genes in Bladder Cancer. Expression data of 309 lipid metabolism-related genes (LMRGs) were collected from the GEO and TCGA cohorts. e flow chart of this research is presented in Figure 1. 89 differently expressed LMRGs were found between normal and bladder cancer tissues, with 57 genes upregulated and 32 genes downregulated in cancer samples when the cutoff was set to FC (fold change, FC) > 1.5 and FDR<0.05 (Figures 2(a) and 2(b)). Univariate Cox regression analysis was conducted on 89 differently expressed LMRGs in TCGA_BLCA cohort. A total of 18 genes with prognostic value were identified with a p value < 0.05 (Figure 2(c)). e somatic mutation profile of 18 LMRGs associated with prognosis was first summarized. A total of 54 of 412 bladder cancer samples experienced mutations of LMRGs, with a frequency of 13.11%. ACOX2, SLC27A2, and ACLY had the highest mutation frequency (Figure 2(d)). Further analyses demonstrated a mutation co-occurrence relationship between ACSF2 and PTGIS, ACOX2 and GPX1, SLC27A2 and CYP1B1, and ACLY and DHCR24 (Figure 2(e)).

Construction and Verification of the Prognostic Index.
To construct a prognostic index of lipid metabolism-related genes, we obtained 404 samples of bladder cancer from TCGA database and divided them into two groups as training cohort (N � 204) and testing cohort (N � 200). e basic characteristics of the patients included are shown in Table 1.
en, the least absolute shrinkage and selection operator (LASSO) Cox regression analysis was conducted to narrow the number of genes. Finally, six genes (CYP1B1, ACOT13, NUDT19, SCD, IL4I1, and DECR1) were used for the construction of prognostic risk score model (Figures 3(a) and 3(b)). en, six genes were shifted from them through multivariate Cox regression (Figure 3(d)). e coefficient of these genes is displayed in Table 2.
e chord diagram showed the correlation between the six genes filtered by LASSO Cox regression analysis (Figure 3(c)). To estimate the effectiveness of this model, we also collected 165 samples from GEO database (GSE13507) as validation. We then calculated the risk score of samples and divided the patients into high-and low-risk score groups according to median risk score acquired from training cohort. We also found that patients in the high-risk score group had the worse survival status than those in the low-risk group (Figures 3(f )-3(k)). We found that clinical characteristics such as stage and status were positively related to high-risk group (Figure 3 ). Furthermore, we conducted metaanalysis based on the three cohorts to demonstrate that the risk score was in good validity (HR � 1.49, 95% CI � 1.15-1.94, p � 0.01) (Figure 3(r)).

e Relationship between Tumor Mutation Burden (TMB)
and Risk Score. We compared the tumor mutation burden (TMB) in high-and low-risk score patients and found the TMB of patients in the low-risk score group was higher (Figures 4(a) and 4(b)). e mutation spectrum of patients with high-risk score (Figure 4(c)) and low-risk score ( Figure 4(d)) was mapped, and 14 significantly mutant genes including (ATAD5, PIK3CA, FGFR3, HUWE1, SPTAN1, GON4L, ALMS1, RELN, STAG2, AHNAK, MED13, UTRN, C2orf16, and ATR) were obtained (Table 3). Survival curves demonstrated that the prognosis of patients with high TMB was better (Figure 4(e)); furthermore, patients with high TMB and low-risk score had the best prognosis ( Figure 4(f )).

Gene Set Variation Analysis (GSVA) and Gene Set Enrichment Analysis (GSEA).
We performed GSEA and GSVA to verify the correlation between risk score and pathways involved in the formation of bladder cancer. We conducted GSVA using Hallmark gene sets, and the results showed that risk score was positively associated with 31 hallmark pathways including mTOR signaling pathway, G2M checkpoint, fatty acid metabolism, and hypoxia ( Figure 6(a)). We performed GSEA in different risk score groups using the Kyoto Encyclopedia of Genes and Genomes (KEGG) signaling pathways (Figures 6(b) and 6(c)). Besides, we confirmed that risk score was associated with pathways in cell death including autophagy and ferroptosis ( Figure 6(d)). Also, we performed metabolism analysis and discovered risk score was positively related to fat-soluble vitamins ( Figure 6(e)).

Response to Chemotherapy and Immunotherapy for
Bladder Cancer of High-and Low-Risk Score Patients. We estimated the response to chemotherapy drugs in different risk score groups and found that patients in the low-risk score were more sensitive to methotrexate sensitivity; furthermore, metformin could be used in bladder cancer patients of low-risk score (Figure 7(a) and 7(b)). We also conducted immunotherapy analysis between the high-and low-risk score groups and found that anti-CTLA4 and anti-PD-1 therapy could make a difference in therapeutic effect between two groups no matter whether used alone or in combination, and patients in the high-risk score groups showed better sensitivity compared with those in the lowrisk score groups (Figure 7(c)).

Immunohistochemistry (IHC) Verification of ACOT13, CYP1, DECR1, IL4I1, and SCD.
We obtained the results of immunohistochemical staining of ACOT13, CYP1, DECR1, IL4I1, and SCD in both normal tissues and bladder cancer tissues and found that the expression of CYP1 was higher in normal tissues, while other genes were higher in tumor tissues than in normal tissues (Figure 8).

Discussion
Bladder cancer is a common malignant tumor with high rates of recurrence [15]. e treatments of bladder cancer have advanced a lot during the past decade; however, the high morbidity and mortality remain unchanged [16]. Effective therapies and prognostic markers still need to be identified. Accumulating studies have indicated that lipid metabolism is a crucial step in metabolic reprogramming, while metabolic reprogramming is a new hallmark of malignant tumors [17]. Despite the importance of lipid metabolism in bladder cancer, few studies were conducted to explore the association between lipid metabolism and bladder cancer.
In this study, we identified the potential mechanism and prognostic value of lipid metabolism-related genes in bladder cancer via bioinformatic analysis. A 6-gene prognostic risk model was constructed by LASSO Cox analysis. In addition, these lipid metabolism-relatedgene-based signatures were tightly associated with TNM stage, T stage, N stage, and status. e patients with low-risk scores were found better outcomes than those with high-risk scores. Our findings showed that this risk model was an independent prognostic prediction for survival and was tightly associated with tumor mutation burden (TMB) and tumor microenvironment (TME).
Recently, more and more studies have suggested that the progression and prognosis of bladder cancer were tightly related to immune cell infiltration [18]. Furthermore, the functions of immune cell could be influenced by metabolic reprogramming [19]. erefore, there must exist a relationship between risk scores and immune cell infiltration. Our results showed that patients in the high-risk score group had higher expression of M0 macrophages, while Tregs and CD8+T cells were upregulated in the low-risk group, demonstrating a differential infiltration pattern between the subgroups. In addition, patients in high-risk groups had higher expression of NRP1, CD44, CD276, and TNFSF9 than those in low-risk groups. ese results indicated that the unfavorable prognosis of patients in the high-risk score groups might be due to the immunosuppressive environment and elevated expression of immune checkpoint genes [20]. We also found that tumor mutation burden was negatively associated with the risk score, and the prognosis of patient with high tumor mutation burden and low score was the best. It might be due to that tumor mutations could generate immunogenic neoantigens, thus leading to immune checkpoint blockade [21].
GSEA and GSVA of hallmarks suggested that G2M checkpoint, EMT pathway, and hypoxia were highly expressed in the high-risk score groups. MT Dillon et al. suggested that inhibition of G2M checkpoint could slow down the progression of tumors [22]. Epithelialmesenchymal transition (EMT) is a process of epithelial cells acquiring mesenchymal features. It is associated with tumor initiation, invasion, metastasis, and resistance to therapy [23]. Hypoxia has emerged as a crucial factor in tumor pathophysiology, and microenvironment promotes altered cellular metabolism including lipid metabolism. Furthermore, reports suggested that hypoxia could trigger EMT in bladder cancer [24]. e pathways of autophagy and ferroptosis were also involved in patients of high-risk score groups, and Enyong Dai et al. suggested that autophagy-  [25]. e poorer prognosis of patients in the high-risk score groups might be tightly associated with mechanisms above. e risk signature was constructed with ACOT13, CYP1, DECR1, IL4I1, NUDT19, and SCD ACOT13, and a member of acy1-CoA thioesterase (ACOT) enzymes can catalyze hydrolysis of fatty acyl-CoA into free fatty acids. It is usually enriched in oxidative tissues and tightly related to mitochondria [26]. Previous studies found that it was associated with many diseases, including lung cancer, pheochromocytomas, and paragangliomas [27]. CYP1 enzymes could catalyze the metabolic activation of procarcinogens and deactivation of certain anticancer drugs. Inhibition of CYP1   Survival probability (f ) Figure 4: Relationship between tumor mutation burden (TMB) and risk score. (a, b) Comparison of TMB between patients in the high-and low-risk score groups, and high tumor mutation burden was associated with low-risk score. (c, d) Mutation spectrum of high-risk and lowrisk patients. (e) Survival analyses for low-and high-TMB patient groups using the Kaplan-Meier curves (P < 0.001, log-rank test). (f ) Survival analyses for four groups of TMB and risk score using the Kaplan-Meier curves. e high-TMB and low-risk score groups indicated better overall survival than the other three groups (P < 0.001, log-rank test). Journal of Oncology * * * *** ** ** *** *** * ** * * *** ** * * * ** ** *** * *   2   NRP1  ICOSLG  CD48  TNFRSF14  KIR3DL1  TNFRSF4  CD44  TMIGD2  IDO2  PDCD1  CD276  TNFSD15  CD40LG  CD40  CD200R1  ICOS  TNFSF9  TIGIT  CTLA4 LGALS9 CD27  BTLA  CD160  CD200  CD200R1  CD244  CD27  CD276  CD28  CD40  CD40LG  CD44  CD48  CD70  CD80  CD86  CTLA4  HHLA2  ICOS  ICOSLG  IDO2  KIR3DL1  LAG3  LAIR1  LGALS9  NRP1  PDCD1  TIGIT  TMIGD2  TNFRSF14  TNFRSF18  TNFRSF4  TNFRSF8  TNFRSF9  TNFSF14  TNFSF15  TNFSF4  TNFSF9   is an effective approach for chemoprevention, and many studies have suggested that inhibitors and prodrug target CYP1 are promising anticancer strategies [28]. DECR1 is an auxiliary enzyme of beta-oxidation, and it participates in redox homeostasis by controlling the balance between saturated and unsaturated phospholipids. Deletion of DECR1 can impair lipid metabolism and reduce tumor growth; therefore, DECR1 is important in the progression of tumor growth and treatment resistance [29]. NUDT19 has been identified to promote the proliferation, migration, and EMT process of tumor via mTORC1/P70S6K signaling pathway [30]. IL4I1 frequently associates with AHR (aryl hydrocarbon receptor) activity and activates the AHR through the generation of indole metabolites and kynurenic acid. In summary, it associates with reduced survival in patients with tumor and enhances the progression of tumor [31]. e previous study indicated that upregulation of SCD could proliferate cancer cells in a lipid-depleted environment for it could synthesize monounsaturated fatty acids. Decreased tumor SCD activity could slow tumor growth [32]. Furthermore, we identified the molecule drugs highly related to lipid metabolism genes for the treatment of bladder cancer. Methotrexate, a well-established antimetabolite, has been used separately or in combination for antitumoral activity for a while [33], and we found it suitable to treat patients with low-risk score. Excepting chemotherapy, immunotherapy is another important treatment for bladder cancer; in our study, we found that anti-CTLA4 and anti-PD1 were sensitive to patients with high-risk score of bladder cancer no matter whether used separately or in combination. ese two drugs are immune checkpoint inhibition and have been licensed for the treatment of bladder cancer [34]. Metformin, often used for diabetes, is known to induce apoptosis in many types of cancers and has the feasibility as a drug repositioning used for the treatment of bladder cancer [35]. Our study indicated that patients in the low-risk score groups were suitable for the treatment of metformin compared with patients in the high-risk score groups.
Our study constructed and validated a prognostic signature model based on lipid metabolism genes, which could predict the prognosis of patients with bladder cancer well and guide the treatment for patients with bladder cancer. Our study also has limitations, and the main limitation of the study is that we do not have experimental studies in vivo and in vitro. Further studies would be conducted to validate what roles the lipid metabolism-related genes play in bladder cancer.

Conclusions
In conclusion, we investigated the lipid metabolism-related genes in bladder cancer through comprehensive bioinformatic analysis. A novel 6-gene signature associated with lipid metabolism for predicting the outcomes of patients with bladder cancer was conducted and validated. Furthermore, the risk score model could be utilized to indicate the choice of therapy in bladder cancer.

Data Availability
e source data of this study were derived from the public repositories, as indicated in the section of "Methods" of the manuscript. All data that support the findings of this study are available from the corresponding author upon reasonable request.

Conflicts of Interest
e authors declare that they have no conflicts of interest.