Identification of DELRGs in BLCA
The RNA-seq levels of 857 LRGs in BLCA (n = 414) and normal bladder tissues (n = 19) in TCGA were acquired. We identified 113 DELRGs with FDR< 0.05 and | log2 FC| >1 as the screening criteria, including 49 downregulated and 64 upregulated genes (Figures 1A and B).
Functional enrichment analysis of DELRGs
Functional enrichment analysis was performed with the 113 DELRGs. In the biological processes (Figure 2A), the DELRGs were mainly enriched in steroid metabolic process, fatty acid metabolic process, lipid catabolic process, steroid biosynthetic process, and fatty acid derivative metabolic process. In the cellular components (Figure 2A), the DELRGs were mainly enriched in lipid droplet, peroxisomal matrix, myelin sheath, microbody lumen, peroxisome, and microbody. In the molecular functions (Figure 2A), the DELRGs were mainly enriched in cofactor binding, monooxygenase activity, iron ion binding, phospholipase A2 activity, and NADP binding. In the KEGG pathways, the results showed that the DELRGs were mainly enriched in glycerophospholipid metabolism, PPAR signaling pathway, fatty acid metabolism, and AMPK signaling pathway (Figure 2B).
Analysis of small molecular drugs
To identify candidate small molecular drugs for treating BLCA, all DELRGs were divided into up-regulated and down-regulated groups, which were uploaded to the CMAP database. We identified five small molecular drugs with anti-cancer functions in the BLCA progression (enrichment score < 0) with p< 0.01 and n >2 as the screening criteria. Table 1 shows the 5 small molecule drugs.
Construction of LRGs signature for predicting OS
29 DELRGs were associated with the prognosis of BLCA patients through univariate Cox regression analysis (Figure 1C). To ensure the stability and feasibility of clinical prognosis based on 29 genes, we obtained 19 DELRGs associated with the prognosis of BLCA patients by lasso survival analysis (Figure 3A and B). After multivariate Cox
regression analysis, 11 genes were identified and used to construct a prognostic signature for OS (Figure 4E). 11 genes, including FASN, MBOAT7, SERPINA6, PPARGC1B, FADS1, CPT1B, HSD17B1, OSBPL10, AKR1B1, CCDC58, and PLA2G2F, were selected to construct LRGs signature for predicting OS. We developed a 11 genes signature-based risk score based on their COX coefficient as follows:
Risk score = (0.2728 × FASN expression) +(-0.2101 × MBOAT7 expression) + (0.8956 × SERPINA6 expression) +(-0.3559 × PPARGC1B expression) + (0.1154 × FADS1 expression) +(-0.2844 × CPT1B expression) + (0.1872 × HSD17B1 expression) +(0.193 × OSBPL10 expression) + (0.1156 × AKR1B1 expression) + (-0.2434 × CCDC58 expression) +(0.064 × PLA2G2F expression)
Patients were then divided into high- and low-risk groups on the basis of the median value. Patients in the high-risk group showed poorer prognosis than those in the low-risk group (p < 0.05) (Figure 3C, E and G). Time‐dependent ROC analysis indicated that the prognostic accuracy of the 11 LRGs signatures in the discovery set was 0.716 at 5- year (Figure 3I).
Validation of the LRGs signature in GEO dataset
To ensure the prediction value of the LRGs signature, we used GSE13507 as validation sets to validate our results. According to the LRGs-based classifier identified above, BLCA patients in the validation sets were divided into a high‐ and a low‐risk group by the median risk score. In accord with the results from above, significantly higher survival rates were observed in low‐risk groups compared with the high‐risk ones in the validation set (Figure 3D, F and H). Time‐dependent ROC analysis indicated that the prognostic accuracy of the LRGs signatures was 0.721 at 5-year (Figure 3J).
Independent Prognostic Role of LRGs Signature
To explore the independence of the LRGs signature, we compared the clinical features including gender, age, grade, stage, T, N, and LRGs signature risk score by performing univariate and multivariate Cox regression analysis. Age (HR = 1.037, 95% CI = 1.020–1.055; p < 0.001), stage (HR = 1.783, 95% CI = 1.44–2.207, p < 0.001), T (HR = 1.569, 95% CI = 1.233–1.819, p < 0.001), N (HR = 1.548, 95% CI = 1.317–1.819, p < 0.001), and risk score (HR = 1.615, 95% CI = 1.424–1.832, p < 0.001) were significantly associated with OS in the univariate analysis (Figure 4A). Multivariate analysis suggested that age (HR = 1.032, 95% CI = 1.014–1.05, p < 0.001) and risk score (HR = 1.479, 95% CI = 1.298–1.685, p < 0.001) were also significantly associated with OS (Figure 4B). Therefore, this indicated that the risk score was an independent prognostic predictor. To better forecast the prognosis of BLCA patients, we constructed a nomogram from the variables associated with OS (Figure 4C). The calibration curve suggested that the nomogram showed good performance in consistent with between the nomogram’s 3‐ or 5‐year OS estimates and the Kaplan–Meier estimates (Figure 4D and E).
Analysis of immune infiltration
In order to define the relationships between the risk scores and immune infiltration, we investigated the enrichment scores of 16 immune cells and their participated pathways by using ssGSEA R package. The results revealed that immune cells (such as DCs, aDCs, pDCs, mast cells, B cells, CD8+ T cells, macrophages, Th2 cells, Th1 and Tfh cells, TIL cell, Neutronphils, Treg cells, T helper cells,) were remarkably elevated in the high-risk groups compared with low-risk group (Figure 5A). In addition, the scores of the immune-related pathways (such as the type I IFN response, parainflamation, type II IFN response, MHC class I, CCR, T cell co-inhibition, cytolytic activity, APC inhibition, inflammation-promoting and check-point) were also generally higher in high-risk group, indicating that their immunological functions related with lipid metabolism were more active in the high-risk group (Figure 5B).