Six immune-related promising biomarkers may promote hepatocellular carcinoma prognosis: a bioinformatics analysis and experimental validation

Identification of differentially expressed miRNAs in HCC

The miRNA microarray data (GSE69580) containing five HCC tumors and five normal tissues were first quantile-normalized before data analysis (Additional file 1: Fig. S1). The filtering criteria were as described previously (|log2FC| > 0.5, adjusted p < 0.05). We identified seven upregulated and two downregulated miRNAs (Table 1). Volcano and heatmaps were drawn to show the differential expression of the nine miRNAs between the tumor and normal tissues (Fig. 2A, B). We then imported nine miRNAs into FunRich (3.1.3) software to perform miRNA GO enrichment analysis. The GO biological process terms (BP) showed that most miRNAs were involved in the regulation of nucleobase, nucleoside, nucleotide, and nucleic acid metabolism (19.7%, p < 0.001) (Fig. 2C). Regarding the GO cellular component terms, most miRNAs may be localized in the nucleus (48.2%, p < 0.001) and cytoplasm (45%, P < 0.001) (Fig. 2D). The GO molecular function terms showed that most of the miRNAs were associated with transcription factor activity (8%, p < 0.001), protein serine/threonine kinase (2.9%, P = 0.001), transcription regulator activity (6.4%, P = 0.007), ubiquitin-specific protease activity (3.3%, P = 0.01), and guanyl-nucleotide exchange factor (1.3%, P = 0.014) (Fig. 2E). For transcription factor analysis, we found that most miRNAs were related to transcription factors (TFs). We chose the top 10 TFs that were closely related to miRNAs. As shown in Fig. 2F, the TFs SP1, EGR1, POU2F1, SP4, MEF2A, FOXA1, SOX1, FOXO1, NKX6-1, and HOXD8 are associated with differentially expressed miRNAs.

Table 1 Differentially expressed miRNAs in miRNA microarray data (GSE69580)Fig. 2figure 2

Identification of differentially expressed miRNAs in HCC. A The volcano plot shows differentially expressed genes (DEGs) (|log2FC| > 0.5, adj P-value < 0.05). B The heatmap shows the different expression level of DEGs. C The biological process terms of miRNAs. D The cellular component terms of miRNAs. E The molecular function terms of miRNAs. F The transcription factors enrichment analysis of miRNAs

Identification of differentially expressed immune-related miRNA targeted-genes

We first uploaded nine miRNAs to three databases (miRecords, miRTarbase, and Tarbase), and then aggregated all validated interacting target genes in these three databases. After that, the final 393 immune-related target mRNA genes were obtained by screening overlapping genes between miRNA targeted genes and the list of genes downloaded from the ImmPort Portal database (Additional file 1: Fig. S2). We then extracted the expression matrix of these 393 genes from the TCGA database and normalized it before differential gene analysis (Additional file 1: Fig. S3). The filter criteria are |log2FC| > 1 and adjusted p-value < 0.05. We confirmed that 97 genes were upregulated and 37 genes were downregulated in HCC tumor tissues. Volcano and heatmaps were drawn to show the differential expression of these genes (Fig. 3A, B). Ultimately, we obtained the final genes by screening the overlapping genes between the upregulated/downregulated mRNAs in HCC tumor tissues and the targets of downregulated/upregulated miRNAs (Fig. 3C).

Fig. 3figure 3

Identification of differentially expressed immune-related miRNA targeted-genes. A The volcano plot shows the differentially expressed mRNAs (DGMis). B The heatmap shows the different expression level of DEMis. C The Venn diagrams show the overlapping genes between the up-/down-regulated mRNAs in HCC tumor tissues and the targets of down-/up-regulated miRNA

Establishment of the network of miRNA-mRNA interaction

The top 20 upregulated and downregulated mRNAs were chosen from the overlapping gene list for further analysis (Tables 2 and 3). As shown in Fig. 4A, B, we constructed the miRNA-mRNA network and utilized the MCODE plugin of Cytoscape software to obtain the sub-network from the whole network. Next, KEGG and GO analyses were performed using the “ClusterProfiler” R package. The results showed that the majority of the genes were enriched in the MAPK signaling pathway, cytokine-cytokine receptor interaction, and PI3K-Akt signaling pathway, etc. (Fig. 4C). Furthermore, most of the genes were related to GO:0050673 (epithelial cell proliferation), GO:0031649 (heat generation), and GO:0043434 (response to peptide hormone) (Fig. 4D). Subsequently, survival analyses of 40 genes were performed. Genes with a p-value < 0.05, were included in further analyses. Therefore, 15 miRNA targets were identified. We then retrieved studies and papers about the 15 miRNA targets in PubMed, and ultimately obtained four targets (NTF3, PSMD14, SORT1, and CD320) that had not been previously studied in-depth in HCC, especially the relationship between immune infiltration and HCC progression (Fig. 4E, Additional file 1: Fig. S4). Finally, we established a ceRNA network of hsa-miR-125b-5p-PSMD14/CD320/SORT1 and hsa-miR-21-5p-NTF3.

Table 2 The top 20 up-regulated mRNAs in 97 immune-related target upregulated genesTable 3 The top 20 down-regulated mRNAs in 37 immune-related target downregulated genesFig. 4figure 4

Establishment of the network of miRNA-mRNA interaction. A Network diagram composed of miRNA and top 20 up-/down-regulated. B The sub-network of miRNA-mRNA network. C, D The KEGG and GO analysis of top 20 up-/down-regulated mRNAs. E The overall survival of DEMs. (*P < 0.05, **P < 0.01, ***P < 0.001, ****P < 0.0001, n.s. not statistically significant)

Correlation between miRNA targeted-genes expression and immune infiltration

We further investigated whether certain links exist between genes and immune cells. The correlations were explored from the Tumor Immune Estimation Resource (TIMER) database, which showed that CD320 was positively associated with B cells, CD4+ T cells, macrophages, and dendritic cells immune infiltration level (Fig. 5A); PSMD14 was positively associated with B cells, CD8+ T cells, CD4+ T cells, macrophages, neutrophils, and dendritic cell infiltration (Fig. 5B); NTF3 was negatively associated with tumor purity and was related to CD4+ T cells, macrophage, and neutrophil infiltration (Fig. 5C), and SORT1 was positively associated with B cells, CD4+ T cells, macrophages, neutrophils, and dendritic cells immune infiltration level (Fig. 5D). Meanwhile, the average expression heatmaps of CD320, SORT1, NTF3, and PSMD14 in various immune cells are shown in Additional file 1: Fig. S5 using the TISCH database.

Fig. 5figure 5

Correlation between miRNA targeted-genes expression and immune infiltration. A Correlation of CD320 expression with immune infiltration level in HCC. B Correlation of PSMD14 expression with immune infiltration level in HCC. C Correlation of NTF3 expression with immune infiltration level in HCC. D Correlation of SORT1 expression with immune infiltration level in HCC. E The correlation between the four hub genes (CD320, PSMD14, SORT1, and NTF3) and the cell markers of corresponding immune cells by Cell Marker database

Based on the above results, we further investigated the correlation between these four hub genes and the cell markers of the corresponding immune cells. Immune marker data were obtained from the Cell Marker database (http://bio-bigdata.hrbmu.edu.cn/CellMarker/). Correlation heatmaps are presented in Fig. 4E.

Relationship between miRNA targeted-genes expression and methylation

Previous studies reported that tumor purity as a confounding factor affects gene expression and DNA methylation levels, and copy number affects gene expression levels, which in turn is related to tumor purity and immune cell infiltration levels [20]. Therefore, we obtained the methylation expression levels of four genes from the UALCAN database. We also in-depth investigated the relationship between gene expression and methylation levels in HCC using the MEXPRESS database. We found that the promoter methylation level of CD320 in HCC tissues was significantly lower than that in normal tissues, and CD320 expression was negatively related to its promoter methylation level (Fig. 6A, E). Moreover, similar to CD320, the promoter methylation levels of SORT1 and PSMD14 were considerably lower in tumor tissues than in normal tissues, and there was also a negative correlation between gene expression and promoter methylation levels (Fig. 6B, C, E). However, the analysis results of NTF3 differed from those of CD320, SORT1, and PSMD14. NTF3 expression was positively associated with promoter methylation levels (Fig. 6D, E). This partly indicates that NTF3 as a tumor suppressor gene and SORT1, CD320, and PSMD14 as oncogenes regulate HCC progression. Next, we performed correlation analyses between the four genes and related methyltransferase genes (DNMT1, DNMT3A, and DNMT3B). The results showed that HCC tumors with high levels of CD320, SORT1, and PSMD14 had high levels of methyltransferase genes (DNMT1, DNMT3A, and DNMT3B), whereas there was no significant correlation between NTF3 and methyltransferase genes (Fig. 6F). Moreover, we analyzed the relationship between miRNA (hsa-miR-21-5p and hsa-miR125b-5p) and methyltransferase genes (DNMT1, DNMT3A, and DNMT3B) in HCC. As shown in Additional file 1: Fig. S6A, B, the results showed that there was a positive relationship between hsa-miR-21-5p and DNMT1 (R = 0.190, p < 0.001) and DNMT3A (R = 0.170, p < 0.001), but there was no significant relationship between hsa-miR-21-5p and DNMT3B (R = 0.013, p = 0.811). Furthermore, a significant negative relationship was observed between hsa-miR125b-5p and DNMT1 (R = − 0.270, p < 0.001), DNMT3A (R = − 0.400, p < 0.001), and DNMT3B (R = − 0.190, p < 0.001). Together, these results provide important insights that altered methylation levels may contribute to the function of these hub genes.

Fig. 6figure 6

Relationship between miRNA targeted-genes expression and methylation. AE Correlation between gene expression (CD320, PSMD14, SORT1, and NTF3) and promoter methylation levels. F The correlation analyses between the four hub genes (CD320, PSMD14, SORT1, and NTF3) and related methyltransferase genes (DNMT1, DNMT3A, and DNMT3B)

Gene set enrichment analysis of miRNA targeted-genes in HCC tissues

To study the downstream pathways of these miRNA-targeted genes, we grouped the matrix of TCGA database according to the gene expression level for GSEA analysis. We chose five of all statistically significant analysis results. CD320 is related to RNA Polymerase (NES = 1.9784905, NOM p < 0.001, FDR = 0.06683415), Pyrimedine Metabolism (NES = 1.9404032, NOM p < 0.001, FDR = 0.063131504), Purine Metabolism (NES = 1.8211129, NOM p < 0 0.001, FDR = 0.15423618), Base Excision Repair (NES = 1.7193841, NOM p < 0 0.001, FDR = 0.13022694), and Proteasome (NES = 2.0096643, NOM p = 0.001953125, FDR = 0.10067051) (Fig. 7A). PSMD14 was associated with oocyte meiosis (NES = 1.9725554, NOM p < 0.001, FDR = 0.114430845), cell cycle (NES = 1.9679518, NOM p < 0.001, FDE = 0.057919133), vasopressin-regulated water reabsorption (NES = 1.9385664, NOM p < 0.001, FDR = 0.05596375), ubiquitin-mediated proteolysis (NES = 1.9356284, NOM p < 0.001, FDR = 0.04440799), and regulation of autophagy (NES = 1.9276263, NOM p < 0.001, FDR = 0.045232568) (Fig. 7B). NTF3 was associated with the calcium signaling pathway (NES = 2.3056374, NOM p < 0.001, FDR = 0.00468132), cytokine-cytokine receptor interaction (NES = 2.2316432, NOM p < 0.001, FDR = 0.003475299), chemokine signaling pathway (NES = 2.116517, NOM p < 0.001, FDR = 0.003489959), TGF-β signaling pathway (NES = 2.1123161, NOM p < 0.001, FDR = 0.003384397), and MAPK signaling pathway (NES = 2.0579696, NOM p < 0.001, FDR = 0.004097538) (Fig. 7C). SORT1 was associated with the mTOR signaling pathway (NES = 1.9278517, NOM p < 0.001, FDR = 0.021522397), pathways in cancer (NES = 1.9270834, NOM p < 0.001, FDR = 0.018652743), VEGF signaling pathway (NES = 1.916239, NOM p < 0.001, FDR = 0.017952878), lysosome (NES = 1.9507663, NOM p < 0.001, FDR = 0.026280008), and neurotrophin signaling pathway (NES = 1.9135792, NOM p < 0.001, FDR = 0.017722571) (Fig. 7D).

Fig. 7figure 7

Gene set enrichment analysis of miRNA targeted-genes in HCC tissues. AD The GSEA analysis of CD320, PSMD14, NTF3, and SORT1

Survival analysis and prognostic model

To investigate the impact of four key genes on the prognosis of HCC, an effective model was established for predicting prognostic status by univariate and multivariate Cox proportional hazards regression analysis. The area under the curve (AUC) of the ROC curve showed that CD320 (AUC = 0.922460), PSMD14 (AUC = 0.937861), SORT1 (AUC = 0.871925), and NTF3 (AUC = 0.966288) were significant predictors (Fig. 8A). In the univariate Cox proportional hazards regression analysis, four genes and tumor stage were identified as prognostic biomarkers (Fig. 8B). Multivariate Cox proportional hazards regression analysis showed that CD320 (HR = 2.484, 95% CI = 1.571–3.928, p < 0.001), PSMD14 (HR = 1.787, 95% CI = 1.153–2.767, p = 0.009), and SORT1 (HR = 1.743, 95% CI = 1.121–2.710, p = 0.014), with significant effects on prognosis, were identified (Fig. 8C-a–c). Meanwhile, NTF3 (HR = 0.649, 95% CI = 0.422–0.998, p = 0.049) was also statistically significant in the multivariate Cox proportional hazards regression analysis (Fig. 8C-d). A nomogram was used for prognostic judgment. “Points” is a scoring scale for each factor, and “Total points” is a scale for total score (Fig. 8D-a–d). Moreover, the results of the calibration analysis suggest that the four prognostic models were in good concordance with the outcomes of HCC patients (Additional file 1: Fig. S7A–D).

Fig. 8figure 8

Survival Analysis and Prognostic Model. A The ROC curve of CD320, PSMD14, NTF3, and SORT1 in HCC. B The univariate cox regression analysis of CD320, PSMD14, NTF3, SORT1 and clinicopathological characteristics. C The multivariate cox regression analysis of CD320, PSMD14, NTF3, SORT1. D The Nomogram shows the prognostic model of CD320, PSMD14, NTF3, SORT1

Validation of miRNA and their targeted genes expressions

We validated the expression of two miRNAs and their target genes in HCC cell lines and HCC samples using qRT-PCR. hsa-miR-21-5p expression was significantly upregulated in six HCC cell lines (HCCLM3, MHCC97H, Hep3B, Huh7, PLC/PRF/5, and HepG2) compared to that in L02 cells (Fig. 9A-a). In contrast, hsa-miR-125b-5p expression was significantly upregulated in L02 liver cells compared to that in the five HCC cell lines (Fig. 9A-b). The mRNA expression of CD320 and PSMD14 was significantly upregulated in six HCC cell lines (HCCLM3, MHCC97H, Hep3B, Huh7, PLC/PRF/5, and HepG2) compared to liver cells L02 (Fig. 9B-a, b). The mRNA level of SORT1 was higher in the five HCC cell lines than in the L02 cells (Fig. 9B-c). In contrast, NTF3 mRNA levels in five HCC cell lines (HCCLM3, MHCC97H, Hep3B, PLC/PRF/5, and HepG2), except for huh7 cell lines, were remarkably lower than those in liver cells L02 (Fig. 9B-d). Interestingly, we found that the expression of miRNAs and mRNAs did not significantly change in Huh7 cells. Then, we selected 15 pairs of human liver cancer and para-cancerous tissues to further detect the expression levels of the four hub genes. The results showed that NTF3 was downregulated in 9 of 15 HCC tissues, PSMD14 was upregulated in 11 of 15 HCC tissues, and SORT1 was upregulated in 9 of 15 HCC tissues (Fig. 9C-a–c).

Fig. 9figure 9

Validation of miRNA and their targeted-genes expressions. A-a, b and B-ad The qRT-PCR results of hsa-miR-125b-5p, hsa-miR-21-5p, CD320, PSMD14, SORT1, and NTF3 in different cell lines. C-ad The qRT-PCR results of CD320, NTF3, PSMD14 and SORT1 in 15 pairs of human liver cancer and para-cancerous tissues. D Western blot of CD320, NTF3, PSMD14 and SORT1 in 12 pairs HCC tumor tissues and matched adjacent normal tissues. E-ad IHC of CD320, NTF3, PSMD14 and SORT1 in 15 pairs HCC tumor tissues and matched adjacent normal tissues. (*P < 0.05, **P < 0.01, ***P < 0.001, ****P < 0.0001, n.s. not statistically significant)

CD320 expression was not significantly different between HCC tumor tissues and matched adjacent normal tissues (Fig. 9C-d). Next, western blot analysis was used to validate the protein levels of PSMD14, SORT1, and NTF3 in 12 pairs of HCC tissues and normal liver tissues (Fig. 9D). Moreover, we further explored the changes of 4 genes protein expression from the Human Protein Atlas database. The results indicated that CD320, PSMD14, and SORT1 protein levels were higher in tumor tissues than in normal tissues (Additional file 1: Fig. S8A, B, D), while the level of NTF3 in the normal tissues was higher than that in the tumor tissues (Additional file 1: Fig. S8C). Similarly, immunohistochemistry assay indicated that the levels of CD320, PSMD14, and SORT1 were higher than those in matched adjacent normal tissues in 15 HCC patients (Fig. 9E-a–c), whereas NTF3 showed the opposite result (Fig. 9E-d). We also utilized the GSE14520, GSE76427, and TCGA expression matrices to investigate the expression levels of the four hub genes (Additional file 1: Fig. S9A–C).

Validation of the potential mechanisms in HCC

According to the previous PCR, MHCC97H and HCCLM3 cell were selected for the further experiments. The PSMD14 and SORT1 expression in MHCC97H and HCCLM3 cells transfected with siRNA were confirmed using qRT-PCR (Fig. 10A). And the CCK-8 and transwell assays suggest that the downregulation of PSMD14 and SORT1 can slow down the growth and migration and invasion of MHCC97H and HCCLM3 cells (Fig. 10B–D). Furthermore, the results of cellular function assay are consistent with recent studies, suggesting that the oncogenic role of PSMD14 and SORT1 in various cancer [21,22,23,24]. However, the potential mechanisms of them in HCC are still unclear. Above GSEA reveals that PSMD14 is related to the autophagy process and SORT1 is association with the mTOR signaling pathway. Therefore, we detected the protein expression levels of these pathways’ biomarker by performing the western blot assay. And we found that the interference of SORT1 downregulates the p-mTOR (Ser2448) expression, leading to the inactivation of mTOR signaling pathway (Fig. 10E). Then, the autophagy process in HCC cells was inhibited after the downregulation of PSMD14. The LC3B and ATG5 expression were significantly decreased and p62 protein was remarkably increased when interfering the PSMD14 (Fig. 10F). Consistent with GSEA, SORT1 can activate the mTOR pathway to enhance the HCC progression. Meanwhile, we conjectured that elevated PSMD14 may maintain tumor cell survival by stimulating autophagy enabling then to ensure their own energy metabolism as well as reduce damage under specific circumstances. All in all, these results reveal the oncogenic role of SORT1 and PSMD14 in HCC cells.

Fig. 10figure 10

Validation of the potential mechanisms in HCC. A The interference of PSMD14 and SORT1 were confirmed using qRT-PCR. B The proliferation of HCC cells was assessed using CCK-8 assay. C, D The migration and invasion of HCC cells assessed by transwell assay. E, F The protein levels of mTOR pathway and autophagy-related biomarkers in HCC cells with si-SORT1 and si-PSMD14 were detected using western blot

Prediction of DEIRGs-related drugs

The miRNA targeted genes screened before were uploaded to DSigDB for drug prediction enrichment analysis (Additional file 1: Table S1). We then chose the top 10 drugs that were related to miRNA-targeted genes. The screen criteria were adjusted with p-value < 0.05. The results of the drug prediction enrichment analysis are shown in Additional file 1: Fig. S10 and Table S1. We found that pinaflavol TTD 00010236, harmaline CTD 00006074, GNF-Pf-3464 TTD 00008265, sulfuretin TTD 00011132, chloroxine TTD 00007143, Tyrphostin B48 CTD 00003485, chloroxine, and Redoxal TTD 00010526 may be associated with PSMD14; AlphaRedisol BOSS may be associated with CD320, and selenium methyl cysteine CTD 00000103 may be associated with SORT1. Some of these drugs have been reported to have anti-cancer effects. For example, harmaline can suppress the growth of liver cancer cells by inducing the p53/p21 and Fas/FasL signaling pathways [25], and may have therapeutic potential for controlling breast cancer invasiveness [26], and chloroxine can facilitate platinum-induced DNA damage to induce cancer cell death in high-grade serous cancer [27]. Therefore, the prediction of DERIGs-related drugs as a good reference advances our future scientific research.

留言 (0)

沒有登入
gif