RNA expression data and related clinical information were obtained from TCGA (https://cancergenome.nih.gov/). A total of 424 samples in TCGA-LIHC were used in the following study as a training cohort ( Supplementary File 1 ). In addition, data from 231 HCC patients from ICGA-LIRI-JP (https://dcc.icgc.org/) were downloaded as an independent, external validation cohort ( Supplementary File 2 ). This research strictly follows TCGA and ICGC access rules and publication guidelines. Detailed information is shown in Table 1 . The starvation-related gene set was obtained from Gene Set Enrichment Analysis (GSEA) ‘GOBP RESPONSE TO STARVATION ‘ in The Molecular Signatures Database(https://broadinstitute.org/gsea/msigdb/). It contained 196 genes responsible for the changes in the state/activity of a cell/organism as a result of a starvation stimulus ( Supplementary File 3 ).

In the training set, we first identified 142 differentially expressed genes associated with starvation metabolism in 374 samples using R (P < 0.05) ( Supplementary File 4 ). Univariate Cox regression analysis (P < 0.05) was used to obtain 39 genes associated with prognosis. Finally, multivariate Cox regression analysis was used to construct a signature containing nine genes, detailed in Table 2 . The risk score for each patient was calculated using the following equation: risk score = (β1 × expression of gene1) + (β2 × expression of gene2) + … + (β9 × expression of gene9). All patients were divided into high- and low-risk groups based on the median risk score. Kaplan-Meier survival curve and two-sided log-rank test were used to compare the overall survival (OS) of the high- and low-risk group patients. Receiver-operating characteristic (ROC) curves were applied to assess the diagnostic efficacy of each clinicopathological characteristic and the prognostic signature. Stratified survival analysis was performed to examine the accuracy of the prognostic signature in predicting patient survival outcomes. Furthermore, we performed univariate and multivariate Cox regression analyses to evaluate whether the risk score was independent in determining the prognosis of the HCC patients. The M and N stages were not analyzed because data were missing for several patients. P < 0.05 was considered statistically significant.

The mRNA expression profile matrix of 231 HCC patients from ICGC was used as an external independent validation cohort to validate the accuracy of the 9 gene signature.

The biological processes, molecular functions, and cell component Gene Ontology (GO) of mRNAs associated with survival were identified using GO enrichment analysis. The main signaling pathways of mRNA regulation were identified using the Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analysis.

We constructed a nomogram by integrating clinicopathologic characteristics, such as age, stage, sex, and grade, as well as the risk score derived from the prognostic signature to analyze the probable 3-year and 5-year OS of the patients with HCC. Calibration charts were used to evaluate the performance of the Nomogram.

GSEA software 4.0.1 was used to identify starvation-related gene sets in 50 HCC tissues and their adjacent tissues. Patients with HCC were divided into low- and high-risk groups based on the median risk value. GSEA was used to further analyze gene expression differences between the high- and the low-risk resistance groups. The Hallmark gene sets (h.all.v7.4.symbols.gmt), KEGG gene sets (c2.cp.kegg.v7.4.symbols.gmt) and GO gene sets (c5.go.v7.4.symbols.gmt) were downloaded from the Molecular Signatures Database (https://www.gsea-msigdb.org/gsea/msigdb/genesets.jsp). The gene sets were filtered using the maximum and minimum gene set sizes of 500 and 15 genes, respectively. The enriched gene sets were obtained based on a P-value < 0.05 and a false discovery rate (FDR) < 0.25 after performing 1,000 permutations.

Liver cancer cells (Hep-3B and Huh-7) were obtained from the Institute of Chinese Academy of Sciences. Hep-3B and Huh-7 cell lines were cultured in DMEM (Corning Incorporated, Corning, NY, USA) supplemented with 10% fetal bovine serum (FBS, Gibco, Thermo Scientific, Waltham, MA, USA). All cells were cultured in an incubator with an atmosphere of 95% O₂ and 5% CO₂ at 37°C.

Total protein was extracted with RIPA lysis buffer (Beyotime Institute of Biotechnology, Jiangsu, China) containing protease and phosphatase inhibitors. The protein concentration was detected using a BCA protein detection kit (Jiangsu Beyotime Biotechnology Research Institute). Equivalent proteins were separated on 10% SDS-PAGE gels and then transferred to polyvinylidene fluoride (PVDF) membranes. After blocking with skimmed milk (dissolved in TBST) for 2 h, the membranes were subsequently probed using antibodies against B-actin (Cell Signaling Technology, USA, 1:10000), EIF2S1 (Cell Signaling Technology, USA, 1:1000),p-EIF2S1 (Cell Signaling Technology, USA, 1:1000), Vimentin(Cell Signaling Technology, USA, 1:1000) and E-cadherin(Cell Signaling Technology, USA, 1:1000) overnight at 4°C. The membranes were then washed with Tris-buffered saline containing Tween and incubated with an HRP‐conjugated anti-rabbit antibody at 37°C for 1 h. Finally, the protein bands on the membranes were observed with an Odyssey Scanning System.

Small interfering RNA (siRNA) was purchased from GenePharma Biological Technology (Shanghai, China). Lipofectamine 2000 was transfected according to the manufacturer’s protocol. Cells were transfected with EIF2S1 siRNAs (siRNA-1: sense, GCCAUAAUCGUCCUCACCA; siRNA-2: sense, CCAUAAUCGUCCUCACCAA, siRNA-2: sense, CCAUAAUCGUCCUCACCAA) at a concentration of 50 nM for 6 h. After 48 h, the treated cells were collected for subsequent experiments.

The cells were immobilized with 4% paraformaldehyde, planted evenly on a slide, infiltrated with TRITON, then sealed with goat serum. The primary and fluorescent secondary antibodies were incubated overnight. Finally, nuclei were counterstained with DAPI.

The ability of cells to migrate and invade was analyzed using a transwell chamber. A total of 8×104 cells were directly and uniformly distributed in the wells in serum-free medium for the migration experiment. Similarly, the invasion test was performed almost identically to the migration test, except that the upper chamber was first coatedwith the matrix according to the manufacturer’s instructions. For both tests, medium containing 15% FBS was added to the lower chamber as a chemical attractant. Cells were incubated in 5% CO₂ for 8 h (migration) or 12 h (invasion). The membrane was wiped with a cotton swab; cells were removed from the upper surface of the cavity, fixed with 4% formaldehyde, and stained with 0.5% crystal violet.

Statistical analyses were conducted using GraphPad Prism 5.0 (San Diego, CA). The data were processed using the PERL programming language (version 5.30.2, http://www.perl.org). All statistical analyses were performed using R (version 3.6.2, https://www.r-project.org/). P < 0.05 was regarded as statistically significant.

GSEA was used to determine whether there were significant differences in the starvation-related gene set between HCC samples and paired adjacent normal samples. The results suggested that the starvation-related gene set was significantly enriched in HCC samples (NES = 1.64, nominal P < 0.001, FDR < 0.001) ( Figure 1A ). A total of 196 starvation-related genes were used in the following study ( Figure 1B ).

In TCGA database, we first identified differentially expressed genes related to hunger (P < 0.05). As shown in Figures 2A, B , 32 of the 142 differentially expressed genes were downregulated, and 110 were upregulated. Then, we identified 39 prognostic differential genes using the univariate Cox regression analysis, among which 37 had a positive and 2 had a negative correlation with risk ( Figure 2C ). It can be seen from the hunger-related prognostic gene protein interaction network in Figure 2D that EIF2S1 is at the core-site. The correlation of hunger-related prognostic genes is shown in Figure 2E .

GO analysis and KEGG pathway enrichment analysis were used to verify whether the genes screened were involved in hunger-related energy metabolism. As shown in Figures 3A, B , the most notable correlation in GO is related to starvation metabolism. The same conclusion was obtained by KEGG pathway enrichment analysis ( Figures 3C, D ). GO analysis and KEGG enrichment analysis further verified that our candidate genes were closely related to hunger metabolism.

We carried out a multivariate analysis of the 39 genes obtained above and obtained 9 genes: EHMT2, HNRNPL, EIF2S1, PPARGC1A, RRP8, FOXK1, CAD FOXK2, and MYBBP1A. Among these genes were four protective genes, those with HR < 1, and five potentially harmful genes, those with HR > 1 ( Table 2 ). We built a signature based on these nine genes and calculated the risk score for each patient based on the resulting model. Based on the median risk score, we divided the patients into low- and high-risk groups ( Figure 4A ). According to the graph, the number of patient deaths increases with increased risk values in the training and validation sets ( Figure 4B ). As can be seen from the training and validation sets, there were significant differences in OS between the high- and the low-risk groups (P < 0.001 and P < 0.001, respectively), indicating a higher mortality rate in the high-risk group ( Figure 4C ). ROC was used to validate the model; the AUC values for 5-year survival for the training and validation cohorts were 0.73 and 0.76, respectively, demonstrating the high accuracy of this model ( Figure 4D ).

Univariate and multivariate analyses were used to determine whether our hunger-related gene model could be an independent prognostic factor. In TCGA, risk scores and sex, age, grade, and stage of starvation-related gene-building models were used for univariate and multivariate analyses. In the univariate analysis, stage (HR = 2.479, 95% CI 1.698–3.619, P < 0.001) and risk score (HR = 1.243, 95% CI 1.182–1.307, P < 0.001) were associated with OS. In the multivariate analysis, stage (HR = 2.047, 95% CI 1.374–3.049, P < 0.001) and risk score (HR = 1.208, 95% CI 1.146–1.273, P < 0.001) were associated with OS ( Figure 5A ). In the ICGC, the risk scores of the starvation-related gene-building prognostic models and clinical characteristics such as age, sex, and stage were also analyzed by univariate and multivariate analyses. The results of the univariate analysis showed that gender, stage, and risk scores were associated with OS. The results of the multivariate analysis showed that gender, stage, and risk scores were also associated with OS ( Figure 5B ).

In the training and validation cohorts, it is evident that the starvation-related gene prognostic model we constructed can act as an independent prognostic factor in patients with liver cancer. areas under the curve (AUC) values in the training and validation sets were 0.745 and 0.731, respectively, indicating high accuracy of the risk score as an independent prognostic factor ( Figure 5C ).

We further analyzed the relationship between the model of genes related to starvation and clinical characteristics. It is evident that the model is not related to age and gender ( Figures 6A, B ) but is closely related to the liver cancer grade and liver cancer stage ( Figures 6C, D ). The higher the stage and grade of the patients with increased risk value, the higher the model of hunger-related gene construction and HCC progression were closely related.

We conducted a stratified analysis of age, sex, staging, and grading to verify the accuracy of our model. We divided the patients into low- and high-risk groups based on the median risk score. The results showed that our model had excellent predictive significance at ages > 65 years and < 65 years for both males and females, grade 1–2 and grade 3–4, stage 1–2 and levels 3-4 ( Figures 7A–D ).

To provide clinicians with a practical clinical tool for predicting 3-year and 5-year OS incidence in liver cancer patients, we constructed a nomogram based on clinicopathological characteristics (age, sex, grade, stage) and risk score based on the 9-mRNA signature ( Figure 8A ). The 3-year and 5-year overall survival (OS) calibration curve is a better predictor than the ideal model ( Figures 8B, C ).

The expression level of EHMT2, HNRNPL, EIF2S1, RRP8, FOXK1, CAD FOXK2, and MYBBP1A increased with an increase in the risk coefficient of the patient, while the expression level of PPARGC1A decreases with an increase in the risk coefficient ( Figure 9A ). The expression level of EHMT2, HNRNPL, EIF2S1, RRP8, FOXK1, CAD FOXK2, and MYBBP1A was significantly higher in cancer, while the expression of PPARGC1A was lower in liver cancer ( Figure 9B ). Figure 9C shows the correlation of genes in the nine models. CAD and HNRNPL had the strongest positive correlation, while PPARGC1A and EHMT2 had the strongest negative correlation. GSEA was conducted to identify starvation-related biological processes and carcinogenic signaling pathways. The results revealed that “Hallmark analysis” gene sets involving cell cycle signals, PI3K/AKT/mTOR pathway, glycolysis, and p53 pathways related to cancer biological processes were enriched in the high-risk group. In addition, several typical pathways from the GO and KEGG genomes, including cell cycle pathways, mTOR signaling pathways, and apoptotic responses, were highly enriched in high-risk phenotypes ( Figures 9D–F ).

As shown in Figures 10A, B , EIF2S1 expression level in liver cancer is increased and is closely related to the degree of malignancy and prognosis of the disease. Thehigher the expression level of EIF2S1, the worse the prognosis of patients. It can be seen from Figures 10C, D that EIF2S1 expression levels were higher in HCC patients with high stage or grade HCC. EIF2S1 expression and phosphorylation levels were higher when HCC cells 3B and Huh-7 were in the starvation state ( Figures 10E, F ). Figure 10G further demonstrated that starvation could increase EIF2S1 expression in HCC cells Huh-7. Hunger and false hunger can induce cancer metastasis (10). We knocked down the expression levels of EIF2S1 in HCC cells 3B and Huh-7 ( Figure 10H ). EIF2S1 knockdown can reduce the invasive and metastatic ability of HCC cells 3B and Huh-7, both under starvation and normal conditions ( Figure 10I ). In addition, upon EIF2S1 knockdown in Huh-7 cells, the protein expression of E-cadherin increased while that of Vimentin decreased, suggesting that EIF2S1 may affect the invasion and metastasis ability of HCC cells through EMT ( Figure 10J ).