Seven Key Genes Correlated with Immune Infiltration Reveal Prognostic Significance of Hepatocellular Carcinoma

Research Article

Austin J Gastroenterol. 2021; 8(1): 1115.

Seven Key Genes Correlated with Immune Infiltration Reveal Prognostic Significance of Hepatocellular Carcinoma

Tian M and Dong L*

Department of Gastroenterology and Hepatology, Shanghai Institute of Liver Diseases, Zhongshan Hospital of Fudan University, Shanghai, China

*Corresponding author: Ling Dong, Department of astroenterology and Hepatology, Shanghai Institute of Liver Diseases, Zhongshan Hospital of Fudan University, Shanghai, 200032, China

Received: June 22, 2021; Accepted: August 09, 2021; Published: August 16, 2021

Abstract

Background: The poor prognosis of Hepatocellular Carcinoma (HCC) is mainly due to late diagnosis, rapid progression and high recurrence rate. Reliable biomarkers for the diagnosis and prognosis of HCC are urgently needed.

Methods: We selected four public datasets from the GEO database to identify Differentially Expressed Genes (DEGs) and Differentially Expressed miRNAs (DE miRNAs). GO functional annotation, KEGG pathway enrichment analysis, and a Protein-Protein Interaction (PPI) network were constructed to explore the functions and importance of DEGs. To determine the target genes of 33 DE miRNAs screened from GSE10694, the miRNet WebServer was utilized. The key genes were screened out by mRNA-microRNA interaction analysis. Then, those highly expressed genes were verified in parallel databases (ONCOMINE, GEPIA and HPA databases). Further prognostic analysis by Kaplan Meier and diagnostic analysis based on TCGA were conducted. Furthermore, we investigated the association between key genes and immune infiltration in HCC tissues using the TIMER database.

Results: We identified seven key genes, including CCNA2, DLGAP5, MAD2L1, MELK, NCAPG, PRC1, and RRM2 based on public databases. The overexpression of these key genes has been demonstrated in HCC tissues and is associated with advanced disease and poor prognosis of patients with HCC. Furthermore, we found these key genes affect various of infiltrating immune cells and be positively correlated with the expression of gene markers of immune escape in HCC.

Conclusion: These seven key genes might be used as biomarkers for the diagnosis, prognosis, and prediction response to immunotherapy for patients with HCC, as well as the therapeutic targets of HCC.

Keywords: Hepatocellular carcinoma; Key genes; Immune infiltration; Immunotherapy; Biomarker

Abbreviations

AFP: a-Fetoprotein; ANTTs: Adjacent Non-Tumor Tissues; AUC: Area Under the Curve; BP: Biological Process; BRD9: Bromodomain-Containing Protein 9; CC: Cellular Component; DCs: Dendritic Cells; DEGs: Differentially Expressed Genes; DE miRNAs: Differentially Expressed miRNAs; DFS: Disease-Free Survival; FC: Fold Change; GEO: Gene Expression Omnibus; GO: Gene Ontology; GPC3: Glypican-3; HCC: Hepatocellular Carcinoma; HPA: Human Protein Atlas; LIHC: Liver Hepatocellular Carcinoma; MCODE: Molecular Complex Detection; MF: Molecular Function; OS: Overall Survival; PFS: Progression-Free Survival; PPI: Protein-Protein Interaction; RNAi: RNA interference; ROC: Receiver Operating Characteristic; TAMs: Tumor-Associated Macrophages; TCGA: The Cancer Genome Atlas; TILs: Tumor-Infiltrating Lymphocytes; TIMER: Tumor Immune Estimation Resource

Introduction

Liver cancer, of which Hepatocellular Carcinoma (HCC) accounts for almost 85-90 %, is the fourth leading cause of cancer-related deaths and the second most lethal tumor with 42,810 estimated new cases and 30,160 deaths in 2020 [1]. Despite clinical improvements in the therapeutic strategies for HCC, such as chemotherapy, targeted therapy, radiotherapy and immunotherapy, the prognosis remains unsatisfactory, with the 5-year-relative survival rate as low as 18%, mainly due to the concealed symptoms and late diagnosis [2-4]. Currently, the clinical screening and diagnosis of HCC mainly relies on imaging methods (computed tomography and ultrasound) and serum biomarkers. Nevertheless, imaging diagnosis of Small HCC (S-HCC) remains challenging. Besides, a-fetoprotein (AFP), the most used serological biomarker, has a sensitivity of only 65% for the diagnosis of clinical HCC, while its prediction for preclinical HCC is less than 40% [5]. Thus, there is an urgent need to identify more accurate molecular biomarkers for prognosis and diagnosis of HCC.

Progression of HCC is considered a multi-step process, in which gene aberrations, cellular context, and environmental factors play a role in the initiation, evolution, and metastasis of the tumor [6]. Studies have been performed to define the genetic mechanisms and pathways involved in HCC, revealing that many molecules such as HepatoScore-14 [7], Bromodomain-Containing Protein 9 (BRD9) [8], Glypican-3 (GPC3) [9], and METTL14 [10] have a significant association with the initiation or/and progression of HCC, and could serve as diagnostic or prognostic biomarkers. However, most of the existing studies are single-center and have been limited by the sample size and validation.

The liver can stimulate immune responses to prevent the invasion of harmful pathogens and the origination of tumors. On the other hand, immune escape plays an indispensable role in HCC [11]. The proliferation of Tumor-Infiltrating Lymphocytes (TILs) is considered representative of the host’s defense to cancer. However, a growing body of evidence reveals that high levels of TILs cannot control the development of the tumor. [12]. Notably, immune checkpoint pathways, including enrichment of many inhibitory co-stimulatory molecules, such as PD-1, LAG-3, CTLA-4, and TIM-3, have been shown to reflect Cytotoxic T Lymphocyte (CTLs) exhaustion [13- 15]. Clinical trials using immune checkpoint inhibitors in HCC, such as tremelimumab (anti-CTLA4), nivolumab (anti-PD1), and durvalumab (anti-PDL1), have shown positive responses. Unfortunately, the therapeutic effect is very limited in some patients with HCC [16]. Therefore, markers and drug targets are urgently needed to predict and improve the efficacy of immunotherapy for HCC.

In recent years, with the development of genome-wide spectrum analysis, our knowledge of the genetic background and driving pathways resulting in HCC has dramatically improved [6]. Microarray technology, the monitoring of genome-wide expression levels of genes in specific organisms, has become an indispensable tool for the classification, diagnosis, prognosis, and detection of therapeutic targets for different types of tumors. In addition, combined with bioinformatics, it has been widely used to describe overall gene expression patterns of multiple types of tumors to obtain a better understanding of the mechanisms behind carcinogenesis [17-19]. In our study, four public datasets consisting of 498 HCC tissues and 357 paired adjacent non-tumor tissues (ANTTs) were selected from the Gene Expression Omnibus (GEO) database [20] to identify reliable biomarkers for HCC diagnosis and prognosis prediction. We also added a series of in-depth analyses based on parallel databases (2,481 samples) to validate the expression of the key genes. The association of key genes with immune infiltration levels and the expression of gene markers in T cell exhaustion was analyzed in HCC.

Materials and Methods

Data source and identification of DEGs and DE miRNAs

We searched ‘Hepatocellular carcinoma’ in the GEO database, a global public library for archiving and free distribution of highthroughput gene expression and other functional genomics datasets [20], with filters activated for ‘DataSets’, ‘Homo sapiens’, and ‘Expression profiling by array,’ to obtain the expression data of mRNA and miRNA of HCC tissues. Finally, three gene expression profiles (GSE54236 based on GPL6480, GSE14520 based on GPL3921, GSE76427 based on GPL10558) and a miRNA expression profile (GSE10694 based on GPL6542) were retrieved after a systematic review. All data were freely accessible.

GEO2R was applied to recognize the DEGs and DE miRNAs between HCC tissues and paired ANTTs in each dataset. The adjusted P values (adj. P) and log Fold Change (FC) in expression were determined. To select DEGs and DEmiRNAs, we set the adj. P and |logFC| for selection as <0.05 and >1.0, respectively. We obtained the overlap of DEGs identified from the three datasets using the Funrich webserver [21] (http://funrich.org/index.html).

Functional pathway enrichment analysis of DEGs and target key genes selection

To explore the functions of DEGs, we conducted GO functional annotation and KEGG pathway enrichment analysis utilizing the Database for annotation, visualization, and integrated discovery database (DAVID database) [22]. Functional annotation consists of a Biological Process (BP), Cellular Component (CC), and Molecular Function (MF). An adj. P<0.05 was statistically significant.

To investigate the importance of these DEGs, a proteinprotein interaction (PPI) network was constructed by the online database STRING (https://string-db.org/) [23], followed by the visual investigation of the interaction networks among different biomolecules in Cytoscape. Then, to determine the hub genes, we used the plugin Molecular Complex Detection (MCODE) with the method of MCC performs. Finally, the top 25 genes with the highest degree values were regarded as hub genes.

To determine the target genes of 33 DE miRNAs screened from GSE10694, we used the miRNet WebServer (https://www.mirnet.ca/) [24] were utilized. The intersection of the target genes and the 25 hub genes were the final identified key genes.

Verification of expression levels based on parallel databases of key genes

The mRNA expression levels of key genes in HCC tissues vs. ANTTs were verified among the five selected datasets (Chen Liver with 197 samples [25], Mas Liver with 115 samples [26], Roessler Liver with 43 samples [27], Roessler Liver 2 with 445 samples [27] and Wurmbach Liver with 75 samples [28]) based on the ONCOMINE database (https://www.oncomine.org) [29].

The mRNA expression level of key genes in liver hepatocellular carcinoma (LIHC) and normal liver tissues were visualized by the GEPIA database (http://gepia.cancer-pku.cn/) [30].

The protein expression levels of the key genes in liver cancer tissues and non-tumor liver tissues were presented by the Human Protein Atlas (HPA, https://www.proteinatlas.org/) [31].

The prognostic and diagnostic value analysis of key genes

The prognostic significance of key genes for HCC were analyzed by Kaplan Meier survival analysis (https://kmplot.com/analysis/ index.php?p=service&cancer=liver_rnaseq). The correlation between Overall Survival (OS)/Progression-Free Survival (PFS) and the expression level of key genes was evaluated in 364/370 patients with liver cancer. Presently, patients with liver cancer were categorized into two groups according to the median expression levels of key genes. P<0.05 was considered as statistically significant.

In addition, the genomic alterations of the key genes were investigated by cBioPortal [32] with a LIHC dataset consisting of 442 samples. We assessed the association of alterations of the key genes with the Disease-Free Survival (DFS) and OS of patients with LIHC.

We downloaded the mRNA-sequencing data and clinical signatures of LIHC and non-cancerous tissues from The Cancer Genome Atlas (TCGA) (https://portal.gdc.cancer.gov/) (N=423 samples). The expression values were normalized by log2 transformation and imported into medcalc software. Receiver Operating Characteristic (ROC) curve analyses were performed to estimate the diagnostic value of key genes. The Area Under the Curve (AUC) represented the diagnostic efficacy of key genes for HCC. P<0.05 was regarded as statistically significant.

Relationship between the expression of key genes and immune infiltration in HCC

TIMER (Tumor IMmune Estimation Resource, https://cistrome. shinyapps.io/timer/) [33] is used to analyze the immune infiltration of various types of cancers in TCGA to explore the influence of key genes on the tumor immune infiltration. Herein, we investigated the association of key genes with the infiltrating levels of immune cells in HCC, including B cells, CD8+ T cells, CD4+ T cells, macrophages, neutrophils, and dendritic cells (DCs). Furthermore, we investigated the relationship between the expression of genetic markers (PD-1, PDL-1, CTLA4, TIM-3, and LAG3) of T cell exhaustion with the expression of key genes.

Association of key genes with clinicopathological profiles of patients with HCC

We determined the correlation of different clinical stages, alcohol consumption and hepatitis virus status of patients with HCC by Kaplan–Meier analysis to explore the association between key genes and clinicopathological features of HCC. In addition, we downloaded the mRNA-sequencing data and clinical signatures of patients with HCC from TCGA. After normalizing the expression values of key genes by log2 transformation, we compared the expression of key genes in different tumor stages (S1, S2 and S3-4) and differentiation grades (G1, G2 and G3-4) with the one-way ANOVA analysis by GraphPad primer. P<0.05 was considered statistically significant.

Results

DEGs and DE miRNAs

In this study, four independent datasets (GSE54236, GSE14520, GSE76427, and GSE10694), including a total of 498 HCC tissues and 357 ANTTs were selected from the GEO database (Table 1). Consequently, 23 overlapping upregulated genes (360 in GSE54236, 572 in GSE14520, and 115 in GSE76427) and 114 overlapping downregulated genes (668 in GSE54236, 755 in GSE14520, and 450 in GSE76427) were identified using a Venn diagram (Figure 1A,1B). In addition, after analyzing the GSE10694 dataset, we identified 33 (16 upregulated and 17 downregulated) DE miRNAs, as displayed in Figure 1C.