A Novel Nine-Gene Signature Associated With Immune Infiltration for Predicting Prognosis in Hepatocellular Carcinoma

Background: Hepatocellular carcinoma (HCC) is one of the most common malignant tumors worldwide, and its prognosis remains unsatisfactory. The identification of new and effective markers is helpful for better predicting the prognosis of patients with HCC and for conducting individualized management. The oncogene Aurora kinase A (AURKA) is involved in a variety of tumors; however, its role in liver cancer is poorly understood. The aim of this study was to establish AURKA-related gene signatures for predicting the prognosis of patients with HCC. Methods: We first analyzed the expression of AURKA in liver cancer and its prognostic significance in different data sets. Subsequently, we selected genes with prognostic value related to AURKA and constructed a gene signature based on them. The predictive ability of the gene signature was tested using the HCC cohort development and verification data sets. A nomogram was constructed by integrating the risk score and clinicopathological characteristics. Finally, the influence of the gene signature on the immune microenvironment in HCC was comprehensively analyzed. Results: We found that AURKA was highly expressed in HCC, and it exhibited prognostic value. We selected eight AURKA-related genes with prognostic value through the protein-protein interaction network and successfully constructed a gene signature. The nine-gene signature could effectively stratify the risk of patients with HCC and demonstrated a good ability in predicting survival. The nomogram showed good discrimination and consistency of risk scores. In addition, the high-risk group showed a higher percentage of immune cell infiltration (i.e., macrophages, myeloid dendritic cells, neutrophils, and CD4+T cells). Moreover, the immune checkpoints SIGLEC15, TIGIT, CD274, HAVCR2, and PDCD1LG2 were also higher in the high-risk group versus the low-risk group. Conclusions: This gene signature may be useful prognostic markers and therapeutic targets in patients with HCC.


INTRODUCTION
According to global data, hepatocellular carcinoma (HCC) is the most common primary liver tumor and the third most common risk factor for cancer-related deaths worldwide (Sung et al., 2021). In 2020, there were approximately 906,000 newly diagnosed patients with liver cancer worldwide and approximately 830,000 liver cancer-related deaths (Sung et al., 2021). The main causes of liver cancer include hepatitis virus infection, smoking, alcoholic cirrhosis, chemical drugs, and aflatoxin infection (Forner et al., 2018). Approximately 400,000 people die annually in China, accounting for >50% of liver cancerrelated deaths globally (Chen et al., 2016). Early-stage HCC is insidious and difficult to detect; consequently, a large number of patients already have advanced-stage disease at the time of diagnosis. At present, liver cancer is mainly treated by surgical resection, supplemented by other methods, such as ablation therapy, targeted therapy, and immunotherapy (Hartke et al., 2017). However, even with timely intervention, the recurrence and mortality rates remain high due to the high degree of malignancy in liver cancer, rapid disease progression, and poor prognosis  There are numerous markers used for predicting the prognosis of patients with liver cancer; nevertheless, their effectiveness is currently limited. Therefore, there is an urgent need to identify more effective biomarkers for predicting the prognosis of patients with liver cancer.
Aurora kinase A (AURKA) is a member of the Aurora kinase family, which consists of AURKA A, B, and C (Carmena and Earnshaw, 2003). Human AURKA is located on chromosome 20q13 and encodes a protein of 403 amino acids. It mainly regulates mitotic spindle formation, stability, and chromosome segregation, and plays an important role in cell cycle regulation (Lindon et al., 2016). The abnormal expression of AURKA can lead to chromosomal abnormalities and instability of the cell genome, which is a risk factor for tumor formation (Wu et al., 2018). AURKA is abnormally expressed in a variety of tumors and regulates tumor proliferation, migration, invasion, and metastasis (Yan et al., 2016) In addition, it is involved in multiple signaling pathways, such as the TP53 pathway, Ras/ mitogen-activated protein kinase (MAPK) pathway and NFKB pathway (Katayama et al., 2004;Briassouli et al., 2007;Umstead et al., 2017). Previous studies have confirmed that AURKA is related to the prognosis of a variety of cancers (breast, colorectal, pancreatic, gastric, and head and neck) and may be a therapeutic target (Jeng et al., 2004;Reiter et al., 2006;Zhang et al., 2015).
Investigations showed that AURKA played an important role in liver cancer progression. Jeng et al. confirmed that AURKA was overexpressed frequently and correlated with high grade and high stage in HCC (Jeng et al., 2004). Lu et al. reported that AURKA mediated c-Myc's oncogenic effects in HCC (Lu et al., 2015). Zhang et al. revealed AURKA promoted chemoresistance through targeting NF-kappaB/microRNA-21/PTEN signaling pathway in HCC (Zhang et al., 2014). Chen et al. suggested that AURKA promoteed cancer metastasis through regulating epithelial-mesenchymal transition and cancer stem cell properties in HCC . However, the specific mechanism of AURKA in HCC is still not very clear and needs to be further explored.
Thus far, no study investigated the role of AURKA gene and AURKA-related prognostic genes in liver cancer. It is well established that the tumor immune microenvironment influences tumor progression (Locy et al., 2018). Currently, the immunological value of AURKA in liver cancer has not been reported. In this study, we first analyzed its clinical value in HCC and selected prognostic genes associated with AURKA. Furthermore, we developed an AURKA-related gene signature in HCC. Next, we constructed a nomogram by combining risk scores and clinical characteristics. Finally, we evaluated the relationship between the gene signature and tumor immunity in HCC.
Differential gene expression analysis of mRNAs was performed based on TCGA database. Limma package of R software (R version 3.6.2) was used to conduct differential gene expression analysis. The adjusted p-value was analyzed to correct for false positive results in TCGA or GTEx. "Adjusted p < 0.05 and Log (Fold Change) >1 or Log (Fold Change) <−1" were regarded as the thresholds of differential expression of mRNAs. The study flowchart was presented in Figure 1.

Protein-Protein Interaction Network and Gene Enrichment Analysis
STRING (version 11.0; http://string.embl.de/) is an open-access biological database that predicts comprehensive interactions of genes at the protein level from multiple organisms including Homo sapiens (Mering et al., 2003). To screen AURKA-related genes, we used STRING to explore AURKA-related genes as well as conduct PPI network analysis on AURKA-related genes. The protein-protein interactions (PPI) with medium confidence >0.4 were regarded as significant.
To further confirm the underlying function of potential targets, the data were analyzed by functional enrichment. Gene Ontology (GO) is a widely-used tool for annotating genes with functions, especially molecular function (MF), biological pathways (BP), and cellular components (CC). Kyoto Encyclopedia of Genes and Genomes (KEGG) Enrichment Analysis is a practical resource for analytical study of gene functions and associated high-level genome functional information. ClusterProfiler package (version: 3.18.0) in R was employed to analyze the GO function of potential targets and enrich the KEGG pathway. p < 0.05 was set as the cut-off criterion.

Identification of Prognostic Genes
Gene Expression Profiling Interactive Analysis (GEPIA) is an interactive web application based on The Cancer Genome Atlas (TCGA) and Genotype-Tissue Expression databases (Tang et al., 2017). We used the Gene Expression Profiling Interactive Analysis to screen the prognostic value of AURKA and AURKA-related genes. The identification of prognostic genes was based on the following criteria: 1) significant differences in gene expression levels between normal liver samples and liver tumor samples and 2) significant association of the gene with both overall survival (OS) and disease-free survival (DFS) of patients with liver hepatocellular carcinoma (LIHC). p-values <0.05 denoted statistically significant differences.

Development and Validation of the Gene Signature
Primary screening of the LIHC data from TCGA database was performed for missing data. After deleting the samples with missing data, The available LIHC data from the TCGA database (340 samples) was divided into two subsets: a training set (n 240) and a validation set (n 100) randomly using the random function of Microsoft Excel program. The training set was used to train a predictive model and the test set was applied for validation. Another database International Cancer Genome Consortium (ICGC) (dcc.icgc.org) was also used as validation set.
Identified prognostic genes were submitted to the multivariate Cox regression model to calculate each prognostic gene's coefficient and risk scores. We used X-tile plot (version 3.6.1, http://www.tissuearray.org/rimmlab) to determine the optimum cutoff of AURKA-related gene signature risk score. X-tile which could calculate the best cut-point of sub-populations, is a software developing by team Rimm Laboratory from Yale University (Camp et al., 2004). A train set and two validation sets were divided into three sub-groups using the same best cut-off value. Risk score analysis, including risk score distribution, survival status, and gene expression heatmap among sub-groups were performed. Kaplan-Meier curves were analyzed using Kruskal-Wallis test and visualized through GraphPad Prism (version 8.0). Receiver operating characteristic (ROC) curve analysis was performed to assess the predictive value of AURKA-related gene signature. The clinical characteristics of three cohorts were shown in Table 1. The correlation among clinicopathological characteristics and risk groups were analyzed by the chisquare test. For sample sizes of less than 40 or theoretical frequencies (T) of less than 1, Fisher exact probability method was used.

Distribution of Somatic Mutations
To identify the somatic mutations of patients with LIHC, we downloaded single-nucleotide polymorphism (SNV) data and clinical follow-up information from TCGA database. The downloaded single-nucleotide polymorphism data were organized in the multiple alignment (MAF) format and visualized using the "maftools" package in R software. The horizontal histogram showed the genes with the highest frequency of mutation.

Immune Infiltration Analysis
To explore the associations between different subgroups and immune cells infiltration, we employed Tumor Immune Estimation Resource (TIMER), which is a useful resource for comprehensive analysis of tumor-infiltrating immune cells . The infiltration of six type of immune cell, including B cell, Macrophage, Myeloid dendritic cell, Neutrophil, T cell CD4 + and T cell CD8 + , were calculated. SIGLEC15, IDO1, CD274, HAVCR2, PDCD1, CTLA4, LAG3 and PDCD1LG2 were selected to be immune checkpoints and the expression values of these eight immune checkpoints among sub-groups were explored. Differences between the three groups were assessed using the Kruskal-Wallis test. p < 0.05 was considered statistically significant. All the above analysis methods and R package were implemented by R foundation for statistical computing (2020) version 4.0.3 and software packages ggplot2 and pheatmap.

Construction and Validation of a Prognostic Model
A nomogram model was constructed to predict the probability of survival at 3 and 5 years for liver cancer patients. Briefly, the prognostic value of the clinicopathological characteristics for OS was estimated through univariate and multivariate Cox regression analyses in both training set and validation set. The performance of the risk model was validated by internal validation and external validation. Internal validation was performed by bootstrap Cox proportional regression analysis based on 1,000 bootstrap samples. Validation set was conducted based on another HCC patients from the TCGA database. Those parameters with p-values <0.05 in both training set and validation set were identified as potential prognostic factors, which were were included in multivariate Cox regression model and visualized using R package "rms." The Calibration curves were plotted to analyze the diagnostic performance of the nomogram. The ROC curve were conducted to determine the clinical value of the nomogram.

AURKA mRNA Level in HCC Samples
The results of the DEG analysis are shown in Figure 2. In TCGA database, a total of 421 HCC samples were selected for the DEG analysis. The analysis showed that the AURKA gene was upregulated (log 2 fold change >1.5 and adjusted p-value <0.05), which indicated that the mRNA expression levels of AURKA differ significantly between normal liver tissue and liver cancer tissue. The mRNA expression levels of AURKA were also determined using three Gene Expression Omnibus data sets (GSE14323, n 115; GSE14520, n 445; and GSE25097, n 557). The results also showed that AURKA was significantly highly expressed in HCC.
The Gene Ontology (GO) enrichment analysis was composed of three parts: GO biological process (GO-BP), GO cellular component (GO-CC), and GO molecular function (GO-MF). In GO-CC ( Figure 3B), these genes were enriched in condensed nuclear chromosome, centromeric region, mitotic spindle, and spindle pole. In GO-BP ( Figure 3C), AURKA and its related genes were significantly enriched in mitotic nuclear envelope disassembly, positive regulation of ubiquitin protein ligase activity, and anaphase-promoting complex-dependent catabolic process. In GO-MF ( Figure 3D), genes were mainly enriched in protein serine/threonine kinase activity, protein heterodimerization activity, and protein kinase binding. For the Kyoto Encyclopedia of Genes and Genomes pathway analysis, nine pathways ( Figure 3E) were observed, namely

Development and Validation of a Nine-Gene Signature
The results of the identification of prognostic genes are shown in Figure 4. According to the screening strategy and criteria described above, two AURKA-related genes, namely CCNB2 (OS Log-rank p > 0.05) and UBE2C (OS Log-rank p > 0.05), were excluded. Eight other AURKA-related genes were significantly highly expressed in HCC tissue compared with normal liver tissue and correlated with prognosis. The AURKA and the other eight AURKA-related genes were used to construct a risk model. The risk score was calculated through multivariate regression analyses. The cutoff value of the risk score was identified using the X-tile software. According to the results, the HCC samples were divided into three subgroups. The cutoff values of the FIGURE 4 | Identification of prognostic genes. AURKA, BIRC5, CDC20, PLK1, TPX2, CDK1, CENPA, DLGAP5 and TACC3 gene expression levels were significantly higher in HCC patients, and these nine genes are related to overall survival and disease-free survival in HCC.
Frontiers in Genetics | www.frontiersin.org November 2021 | Volume 12 | Article 730732 6 risk score were 0.17 and 0.95. HCC samples with a risk score <0.17 and >0.95 classified into the low-and high-risk groups, respectively. The remaining samples were assigned to the moderate-risk group.
The details of the risk model showed in Figure 5, which revealed that higher risk scores were associated with higher expression levels of the nine-gene signature. Furthermore, higher risk scores also indicated worse OS. These results were similar in both the training and other two validation sets. Collectively, these results suggest that the risk model had potential value in predicting the prognosis of HCC. Furthermore, as the ROC curves showed, the area under the ROC curve (AUC) of the training set and two validation sets were higher than 0.5, which indicated that AURKA-related gene signature risk model had important value in predicting prognosis ( Figures 5D-F).
The clinicopathological characteristics among the risk groups in training set were shown in Table 2. Age, T stage, pathological TNM (pTNM) stage, grade, recurrence, and survival endpoint were significantly different between the three risk groups (p < 0.05). Somatic Mutation Results Figure 6 illustrates the somatic landscape of the three risk subgroups. Information on the mutation status of each gene in each sample was shown in the waterfall plot, where different colors with specific annotations at the bottom indicated the various types of mutation. The barplot above the legend exhibited the number of mutations. The results showed that CTNNB1 was the most commonly mutated gene in the lowrisk group. Tumor protein p53 (TP53) was the most frequently mutated gene in both the moderate-and high-risk groups. Hence, we further grouped the HCC samples into two groups based on the TP53 mutation status. A total of 101 and 241 HCC samples were assigned to the TP53 mutant-and wild-type groups, respectively. The survival analyses of both the TP53 mutantand wild type cohorts yielded similar results. Higher risk scores were associated with worse prognostic outcome (p < 0.05).

Association Between Risk Score and Immune Infiltration
We further investigated differences in the degree of immune infiltration in various risk groups. Using the TIMER database, we evaluated the immune cell infiltration in samples from the three aforementioned groups. As shown in Figure 7, we found that the level of immune cell infiltration (macrophages, myeloid dendritic cells, neutrophils, and CD4T cells) was significantly different between the three groups. To further investigate the levels of immune cell infiltration on the gene level, the following eight immune-related genes were selected: sialic acid binding Ig like lectin 15 (SIGLEC15), T cell immunoreceptor with Ig and ITIM domains (TIGIT), CD274, hepatitis A virus cellular receptor 2 (HAVCR2), programmed cell death 1 (PDCD1), cytotoxic T-lymphocyte associated protein 4 (CTLA4), lymphocyte activating 3 (LAG3), and programmed cell death 1 ligand 2 (PDCD1LG2). We compared the expression levels of immune-related genes among the three risk groups. The results are shown in Figure 8. Except for CTLA4, LAG3, and PDCD1, the expression levels of the other five immune-related genes differed significantly between the three groups. Overall, the high-risk group showed significantly higher levels of immune gene expression and immune cell infiltration compared with the other groups.

Construction and Evaluation of the Nomogram Model
As shown in Table 3, we performed the univariate and multivariate analyses using the SPSS software (version 23.0; IBM Corporation, Armonk, NY, United States) to identify independent prognostic factors predicting OS in patients with HCC. In the training data set, the risk score, T stage, M stage, pTNM stage, and recurrence were identified as independent prognostic factors. The results of the testing and training data sets were similar. Overall, the risk score, T stage, and pTNM stage were identified as independent prognostic factors in both data sets. Hence, we used age, T stage, N stage, M stage, pTNM stage, recurrence, and risk group as estimated factors in the construction of our model. A nomogram was constructed to estimate the probabilities for three-and 5-year survival. Calibration curves analysis showed that the new nomogram model had good predictive accuracy ( Figure 9B). Furthermore, the new nomogram model achieved an area under the curve (AUC) of 0.78 at 5 years, which was better than that of a model without the gene signature (AUC 0.73) or pTNM stage (AUC 0.70) ( Figure 9C).

DISCUSSION
The prognosis of patients with HCC varies greatly. The 5-year survival rate after resection of early HCC can be as high as 70%, whereas that of patients with vascular invasion or advanced HCC is markedly lower (Mazzaferro et al., 2011;Grandhi et al., 2016). Therefore, the use of effective prognostic indicators can promptly identify high-risk patients and assist in implementing individualized treatment to improve the prognosis.
In the present study, we found that the AURKA gene was highly expressed in HCC and an independent prognostic risk factor. The GO analysis of AURKA and related genes indicated that these genes are involved in a variety of biological process, including nuclear chromosome condensation, centromeric regions, mitotic spindle and spindle poles, mitotic nuclear envelope disassembly, and ubiquitin protein ligase activity. Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analysis showed that AURKA and related genes regulate the progression of HCC through multiple pathways, such as the TP53 signaling pathway and FOXO signaling pathway. Numerous studies have demonstrated that the TP53 signaling pathway is involved in the development of a variety of tumors and plays a regulatory role in tumor immunity (Stegh, 2012;Blagih et al., 2020;Muñoz-Fontela et al., 2016). The role of the FOXO signaling pathway in tumors has also received extensive attention (Farhan et al., 2017;Coomans de Brachène and Demoulin, 2016). In general, the enrichment analysis revealed some potential mechanisms and possible pathways of AURKA and its related genes in HCC. In addition, it provided some new ideas for the treatment of patients with HCC. Among the nine genes associated with AURKA, eight genes (BIRC5, CDC20, PLK1, TPX2, CDK1, CENPA, DLGAP5, and TACC3) were identified as prognostic genes in HCC. Therefore, we used these nine genes to construct the AURKA-related gene signature for the prediction of prognosis of patients with liver cancer. In gastric cancer, high expression of BIRC5 promotes gastric cancer metastasis and is associated with poor prognosis (Zou et al., 2019). However, in lung cancer, high expression of BIRC5 may prolong OS and DFS (Vischioni et al., 2004). In pancreatic cancer, CDC20 can promote tumor cell proliferation and affect the progression of pancreatic cancer (Chang et al., 2012). CDK1 is considered a synthetic target for KRAS-mutated  tumors and has been identified as a prognostic marker in numerous types of cancer (Sung et al., 2014;Costa-Cabral et al., 2016). CENPA plays a key role in cell mitosis. Of note, it is abnormally expressed in tumors and regulates tumor cell activity (Valdivia et al., 2009). DLGAP5 expression is regulated by the ubiquitin-proteasome pathway and participates in tumor cell migration and invasion . PLK1 is a key regulator of mitosis and is involved in multiple stages of mitosis . Downregulation of PLK1 can inhibit the invasion and metastasis of esophageal cancer cells (Li et al., 2014). In colorectal cancer, PLK1 may promote the growth, invasion, and metastasis of colorectal cancer cells through the PDKI-PLK1-MYC signaling pathway (Tan et al., 2013). TACC3 can activate the Akt/RAS/ mitogen-activated protein kinase kinase/extracellular signal-regulated kinase (Akt/RAS/MEK/ERK) signaling pathway to promote the malignant transformation of cells (Burgess et al., 2018). High TACC3 expression has been found in a variety of tumors and is closely related to poor prognosis . TPX2 is mainly involved in centrosomal maturation and spindle formation (Gruss and Vernos, 2004). In gliomas, MiR-1294 can target TPX2 to inhibit tumor cell proliferation and enhance sensitivity to chemotherapy (Chen et al., 2018). Several previous studies have successfully constructed multi-gene signatures for risk stratification and prognosis prediction in HCC Ouyang et al., 2020). In this study, we constructed a new nine-gene signature to predict the prognosis of patients with liver cancer. This gene signature was verified using an internal verification data set. Gene signatures can effectively classify patients into high-, moderate-, and low-risk groups. Higher risk scores in the training and validation data sets indicated a poor prognosis in HCC. Finally, we constructed a personalized nomogram based on the risk scores, with a concordance index of 0.78.
As an important tumor suppressor gene, TP53 plays a vital role in cell cycle regulation. TP53 mutation is a common mutation in tumors and the most important mutation in liver cancer. This mutation can promote the proliferation, migration, and invasion of tumor cells and increase resistance to drugs (Warren et al., 2013). We divided the HCC cohort of TCGA data set into two groups (TP53 mutation and wild type) and investigated the relationship between the gene signature and these two groups. We found that the gene signature could effectively predict the risk of patients in the TP53 mutation group, thereby effectively classifying these patients into low-, moderate-, and high-risk groups. Numerous recent studies have also confirmed the close relationship between TP53 mutation and tumor immunity (Long et al., 2019;Wu et al., 2020;Sun et al., 2020). Based on the immunoprognostic model established by TP53 mutation, Long et al. found that the levels of T cell follicular helper proteins, T cell regulatory proteins, and macrophages M0 were higher in the high-risk HCC group versus the other groups (Long et al., 2019). Notably, the expression of immune checkpoints CTLA4, PDCD1, and T-cell immunoglobulin mucin family member 3 (TIM3) were also higher in the high-risk group. The investigators suggested that TP53 mutations significantly reduced the immune response in liver cancer.
Immune cell infiltration affects tumor progression. Numerous immunotherapies have been used to regulate immune cells in tumors. Therefore, we investigated the immune cell infiltration in different risk groups. We found that the numbers of macrophages, myeloid dendritic cells, neutrophils, and CD4 + T cells differed significantly in different risk score groups. Higher risk scores were linked to higher numbers of these four types of immune cells. Macrophages play a dual role in the tumor microenvironment, promoting tumor formation and development as well as inhibiting tumor growth (Kim and Bae, 2016). It has been confirmed that the degree of macrophage infiltration in the tumor microenvironment is related to prognosis (Conway et al., 2016). It is currently thought that neutrophils in the tumor microenvironment directly kill or Frontiers in Genetics | www.frontiersin.org November 2021 | Volume 12 | Article 730732 stimulate the immune system to inhibit tumor cells and can also promote immune escape of tumor cells (Kim and Bae, 2016). Myeloid dendritic cells mainly play an antigen-presenting role to activate T cells and induce immune responses (Garris and Luke, 2020). CD4 + T cells mainly support CD8 + T toxic lymphocytes and enhance their anti-tumor immune effect (Ghiringhelli et al., 2006). Our results showed differences in the distribution of immune cells in the tumor microenvironment among the different risk groups and revealed the underlying reason for the poor prognosis observed in the high-risk group. To the best of our knowledge, this is the first study to investigate the relationship between AURKA and related genes and tumor immunity. Immune checkpoint inhibitors are a new approach to the treatment of tumors. Several immune checkpoint inhibitors have been used effectively in the treatment of liver cancer (Liu and Qin, 2019). In this study, we also assessed the relationship between the risk score and immune checkpoints. The results found that the expression of the five immune checkpoints (SIGLEC15, TIGIT, CD274, HAVCR2, and PDCD1LG2) varied in different risk groups. Higher risk scores were associated with higher expression of the five immune checkpoints. These five immune checkpoints play an important role in the activation of T cells. PDCD1LG1 and PDCD1LG2 are the two ligands of PDCD1. In the tumor microenvironment, PDCD1 on the surface of immune cells binds to the PDCD1LG1 and PDCD1LG2 receptors on the surface of tumor cells to activate a series of signal factors in immune cells. This process initiates a series of signaling factors in immune cells to inhibit T cell activation and promote T cell failure, thus helping tumor cells to evade immunosurveillance (Barclay et al., 2018;Ai et al., 2020). SIGLEC15 is a newly discovered immune checkpoint. Wang et al. reported that high expression of SIGLEC15 in tumors can significantly inhibit the activity of T cells . Furthermore, the inhibition or knockout of SIGLEC15 expression can improve the anti-tumor ability of T cells in mice . HAVCR2 is thought to play a dual role, inducing immune tolerance and promoting tumor cell apoptosis (Das et al., 2017). TIGIT is mainly expressed on T cells and natural killer cells. Joller et al. observed significant T cell proliferation in TIGIT-knockout mice (Joller et al., 2011). The expression of these immune checkpoints significantly affects tumor prognosis. At present, the use of single immune checkpoint blockers or combinations of these agents has shown good efficacy in different tumors (Darvin et al., 2018).
This study is characterized by several limitations. Firstly, all analyses were based on public databases. The specific mechanisms of AURKA and related genes in liver cancer have not been thoroughly investigated. Furthermore, the gene signature was associated with immune infiltration and immune checkpoint expression in HCC and affected the prognosis of patients with this disease. Further studies are needed to examine the value of gene signatures in immune invasion and prognosis in HCC.
In summary, we found that nine AURKA-related genes with prognostic value can be used as prognostic markers for liver cancer. The gene signature based on AURKA successfully classified patients with liver cancer into high-, moderate-and low-risk groups. Hence, the gene signature can may be an effective marker for the prognosis of HCC. In addition, the risk score was related to immune cell infiltration and immune checkpoint expression in HCC.

DATA AVAILABILITY STATEMENT
The original data presented in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding author.