Data acquisition

Gene expression data for HCC and non-tumoral samples, along with clinical information for patients of Caucasian and Asian descent, were obtained from the TCGA-LIHC cohort (https://portal.gdc.cancer.gov/) [10, 11]. The dataset was obtained through open-access data, with demographic filters applied to select patients based on their self-reported ethnicity, as recorded in the clinical metadata, specifically focusing on Caucasian and Asian populations. All data were downloaded and merged into a raw count table in CSV format.

Differential gene expression analysis

To identify differentially expressed genes (DEGs) between Caucasian and Asian HCC, the PyDESeq2 implementation of DESeq2 was utilized (https://pydeseq2.readthedocs.io/en/latest/). A DeseqDataSet object was created using the filtered count matrix and corresponding metadata. The design matrix was set to compare the “Caucasian” and “Asian” groups. For normalization, DESeq2 uses a method called size factor normalization to account for differences in sequencing depth between samples [12]. This step adjusts raw counts by dividing each sampleʼs counts by a sample-specific size factor, which ensures that comparisons across samples are not biased by library size.

Dispersion estimation, differential expression thresholding, and gene analysis

Dispersion estimates were calculated for the data to model the gene-specific variance. The results were extracted using the DeseqStats function, which provided fold changes (FCs) and adjusted p-values of false discovery rate (FDR) for each gene. Genes with FDR < 0.05 and |log₂FC| ≥ 1.0 were considered significantly differentially expressed. Genes with log₂FC ≥ 1.0 were labelled upregulated, and those with log₂FC ≤ −1.0 were labelled downregulated. These genes were retained for further analysis. Log₂FC values were computed, and genes with log₂FC ≥ 1.0 or ≤ –1.0 were classified as significantly upregulated or downregulated, respectively, and included in downstream analyses.

The analysis pipeline was implemented in a Jupyter Notebook environment (https://jupyter.org/) using Python packages including pandas  (https://pandas.pydata.org/), NumPy  (http://numpy.org/), Matplotlib  (https://matplotlib.org/), and seaborn  (https://seaborn.pydata.org/). All procedures were conducted with an emphasis on reproducibility and transparency in gene selection and data visualization.

Gene annotation and visualization

To facilitate biological interpretation, gene symbols were mapped using the Sanbomics package (https://github.com/mousepixels/sanbomics). The Ensembl gene identifiers in the DESeq2 results were converted to human gene symbols based on mappings provided by the package using id_map tools.

Visualization of expression trends was carried out using the seaborn library (https://seaborn.pydata.org/). Heatmaps were generated to represent the log₂FC values of the conserved genes. Annotated heatmaps were produced using a diverging color scale centered at zero to distinguish upregulated and downregulated patterns clearly.

Survival analysis

The data of the TCGA-LIHC and the Genotype-Tissue Expression (GTEx) (https://www.gtexportal.org/home/) portals [10, 11, 13] were also analyzed and visualized by the GEPIA online tool [14]. Several biomarkers previously reported to play roles in HCC, including glypican-3 [15], Plasmalemmal Vesicle Associated Protein (PLVAP) [16], ficolin-3 [17], and OIT3 [18], were also considered in the interpretation of expression results. The prognostic relevance of DEGs on overall survival was evaluated using validated data from the Human Protein Atlas (HPA) (https://www.proteinatlas.org/) that uses the TCGA-LIHC cohort [19]. Patients were classified into two groups (high and low expression) based on the best expression cut-off, which refers to the fragments per kilobase of transcript per million mapped reads (FPKM) value that yields the maximal difference with regard to survival at the lowest log-rank p-value.

Gene set enrichment analysis (GSEA)

To gain insights into biological processes (BPs) enriched in the DEGs, GSEA was performed using the GSEApy package (https://gseapy.readthedocs.io/en/latest/introduction.html). The gene rankings based on DESeq2 statistics were used as input. The analysis focused on Gene Ontology (GO) terms, including BPs, cellular components (CCs), and molecular functions (MFs). The enrichment results were visualized using bubble plots.

Samples characteristic

The HCC cohort from the TCGA-LIHC database included 184 Caucasian and 157 Asian subjects, who were predominantly male. Notably, the Asian group had a markedly higher proportion of male subjects, approximately three times more than females. The mean age was significantly different between the two groups, 63.5 ± 13.9 years for Caucasians and 55 ± 13.9 years for Asians (p = 0.0001). Additionally, most of the Caucasian subjects were elderly patients (71.19%), whereas the Asian cohort had a younger age distribution. As controls, 34 non-tumoral samples were included for the Caucasian cohort and 6 for the Asian cohort. Both cohorts included patients with HCC ranging from stage I to stage IV, with the majority diagnosed at an early stage (stage I).

Common DEGs in Asian and Caucasian LIHC populations compared to non-tumoral tissues

Compared to non-tumoral tissues, a total of 503 genes were upregulated in the Asian LIHC and 880 genes in the Caucasian LIHC populations, of which 387 genes were upregulated in both populations (Figure 1A, Table S1). Among them, 116 genes were uniquely upregulated in Asians, while 493 were exclusive to Caucasians. In both groups, 250 genes were commonly downregulated; additionally, 126 genes were uniquely downregulated in Asians, and 288 in Caucasians (Figure 1A, Table S2). The full set of the DEGs in both populations—potentially representing shared molecular markers of HCC—is provided in Table S3.

Several genes were consistently dysregulated in both LIHC populations compared to the non-tumoral samples (Figure 1B). Notably, GPC3 and PLVAP were upregulated, while FCN3 and OIT3 were downregulated in both the Asian and Caucasian populations. In contrast, distinct sets of genes exhibited population-specific expression patterns. In the Asian population, genes such as AKR1B10, UBE2C, and S100P were significantly upregulated, while CSRNP1, RND3, and RCAN1 were significantly downregulated. Meanwhile, in the Caucasian cohort, MDK, LCN2, and NQO1 were significantly upregulated, whereas ECM1, TRIB1, and CYP2C19 were significantly downregulated.

Display full size

Figure 1. Differentially expressed genes (DEGs) in Asian and Caucasian LIHC cohorts compared to non-tumoral tissues. (A) Venn diagrams of upregulated and downregulated genes compared to non-tumoral tissues. (B) Volcano plots of downregulated and upregulated genes in hepatocellular carcinoma (HCC) in Asian [left panel] and in Caucasian [right panel] populations, compared to non-tumoral samples. FDR: false discovery rate; DE: differentially expressed.

Prognostic relevance of common DEGs in both populations

Using GEPIA2 (http://gepia2.cancer-pku.cn/), the expression levels of these genes were obtained (Figure 2). Based on the HPA, the 5-year survival of patients with high GPC3 expression was 46% against 60% in low-expressing patients, whereas patients with high PLVAP expression showed 54% versus 46% survival probability. Poorer prognosis was seen in patients with low FCN3 expression (37% vs. 55% in high-expressing patients) and with low OIT3 expression (42% vs. 63% in high-expressing patients).

Display full size

Figure 2. Differential gene expression in the TCGA cohort and its association with overall survival in Asian and Caucasian cohorts. Boxplot comparing the expression levels of consistently dysregulated genes in both cohorts. Red box: tumor tissue, grey box: normal tissue. *: p < 0.05.

The expression of these genes in the TCGA-LIHC cohort, as well as the survival curves, is presented in Figure 3. Of particular interest is AKR1B10, which is distinctly upregulated in the Asian cohort, in which the 5-year survival probability in high-expressing patients was only 43% versus 58% in low-expressing patients. Overall, the high expression of genes in the Asian cohort resulted in a lower 5-year survival probability. In the Caucasian cohort, high MDK showed a poor 5-year survival rate of 38% versus 52% in low MDK patients. The same trend was observed for NQO1 (39% in high-expressing versus 52% in low-expressing) and for LCN2 (43% versus 51% in low-expressing).

Display full size

Figure 3. Differential gene expressions in the TCGA cohort and its association with overall survival. Boxplot comparing the expression levels of Asian-specific upregulated genes (AKR1B10, S100P, and UBE2C) and Caucasian-specific upregulated genes (MDK, NQO1, and LCN2). TPM: transcripts per million; LIHC: Liver Hepatocellular Carcinoma. Red box: tumor tissue, grey box: normal tissue. *: p < 0.05.

To further explore the difference between the two groups, gene expression profiles between Caucasian and Asian cancer samples were compared using each population alternately as a reference. This analysis revealed several genes that were differentially expressed depending on the direction of comparison.

To complement the transcriptomic findings from the gene expression analysis, Figure 4 presents a list of genes that are differentially expressed in Asian and Caucasian cohorts. Heatmap analysis of DEGs revealed distinct transcriptional profiles in each of the Asian and Caucasian HCC populations. In the Asian population, 54 genes were upregulated compared to Caucasians, with IFITM3P2 being the most highly upregulated gene with a log₂FC of around 2.7, followed by CYP2D6 and PAGE5. In contrast, 25 downregulated genes included UGT2B17 (−3.80 log₂FC), LINC01595, and AVPR1A, suggesting suppressed metabolic and regulatory pathways in the Asian group.

When comparing gene expression between the Caucasian and Asian cohorts, 100 genes showed higher expression in the Caucasian group, with UGT2B17 (log₂FC = 3.80), LINC01595, and AVPR1A standing out as the most upregulated. These genes are associated with increased xenobiotic metabolism and immune regulatory activity. Meanwhile, 16 genes, such as IFITM3P2, were significantly downregulated. The consistent q-value distribution across both populations supports the reliability of these findings.

Display full size

Figure 4. Differentially expressed genes between Asian and Caucasian HCC patients. Heatmaps displaying significantly upregulated and downregulated genes in Asian (top panel) and Caucasian (bottom panel) cohorts, based on RNA-seq data from the TCGA-LIHC dataset. The color scale represents the magnitude of log₂FC. FC: fold change; HCC: hepatocellular carcinoma; TCGA: The Cancer Genome Atlas; LIHC: Liver Hepatocellular Carcinoma.

To simplify Figure 4, differential gene expression analysis between Asian and Caucasian HCC samples identified a distinct subset of genes with opposing expression patterns across the two populations, as shown in Figure 5. A total of 16 genes were significantly upregulated in the Asian cohort but downregulated in the Caucasian group, while 25 genes showed the reverse trend—upregulated in Caucasians but downregulated in Asians. The bar plot illustrates the log₂FC of these genes, with positive values indicating higher expression in Asians and negative values indicating higher expression in Caucasians. Among the 16 genes upregulated in Asians, IFITM3P2 had more than 2.5-fold higher expression in Asians compared to Caucasians. Conversely, among the 25 genes downregulated in Asians but upregulated in Caucasians, UGT2B17 and LINC01595 showed over 2.5-fold higher expression in Caucasians.

Display full size

Figure 5. Genes with opposite expression patterns in Asian and Caucasian HCC patients. Bar plot showing log₂ fold changes of genes that are upregulated in one population but downregulated in the other, based on TCGA-LIHC data. Positive values indicate higher expression in Asians, and negative values indicate higher expression in Caucasians. DEGs: differentially expressed genes; FC: fold change; HCC: hepatocellular carcinoma; TCGA: The Cancer Genome Atlas; LIHC: Liver Hepatocellular Carcinoma.

Enriched BPs associated with the alteration of gene expression

Both groups underwent GO BP enrichment analysis to gain a better understanding of the shared functional impact of upregulated and downregulated genes between the two populations. The identification of the biological pathways and functional roles associated with the DEGs highlighted key BPs in each group. The enriched BP, MF, and CC keywords for commonly upregulated and downregulated genes are shown in Tables S4 and S5, respectively.

Functional enrichment analysis of the DEGs

To identify population-specific molecular signatures in HCC, functional enrichment analysis was performed on DEGs between the Asian and Caucasian populations. In the Asian population, upregulated genes were significantly enriched in BP, such as positive regulation of ubiquitin protein ligase activity, mitotic metaphase/anaphase transition, and anaphase-promoting complex-dependent catabolic processes (Figure 6A). Enrichment in CC included blood microparticles, extracellular exosomes, and the anaphase-promoting complex, while MF terms involved antigen binding, structural constituents of ribosome, and protein tyrosine/serine/threonine phosphatase activity. Conversely, downregulated genes in the Asian population were associated with BP, including dibenzo-p-dioxin metabolic process, porphyrin-containing compound metabolism, steroid catabolic process, and various chemotactic responses such as lymphocyte and macrophage chemotaxis. These genes were also enriched in CC, such as high-density lipoprotein particles, endocytic vesicles, and cytoplasmic microtubules. In terms of MF, significant reductions were observed in activities related to estrogen hydroxylase, glutathione transferase, and enzyme binding.

In the analysis of upregulated genes in the Caucasian population (Figure 6B), several BP were significantly enriched, including positive regulation of interleukin-1β secretion, macrophage and neutrophil migration, cytokine production, dibenzo-p-dioxin and porphyrin-containing compound metabolism, peroxide catabolism, and positive regulation of cell adhesion molecules. At the cellular level, these genes were predominantly associated with the extracellular space, exosomes, and the cell surface. Additional enrichment in mitochondrial and plasma membrane components highlights the interplay between energy metabolism, vesicle trafficking, and membrane-associated signaling in HCC progression. MF analysis revealed dysregulation in redox homeostasis and hormone metabolism, with significant enrichment in oxidoreductase and hydroxylase activities, particularly those involved in estrogen and fatty acid metabolism. Binding activities related to heme, iron ions, and cofactors were also prominently enriched.

Conversely, downregulated genes in the Caucasian population were enriched in pathways related to xenobiotic metabolism, xenobiotic catabolism, and arachidonate metabolism. At the cellular level, these genes were associated with blood microparticles, the immunoglobulin complex, and mitochondria. MF analysis indicated suppression of key activities, such as antigen binding, aromatase activity, oxidoreductase activity, and heme/iron ion binding.

Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analyses revealed distinct pathway signatures between ethnic cohorts and controls (Figures S1–4). In the Asian cohort compared with controls, the cell cycle pathway was the most significantly enriched (FDR = 7.2 × 10^–5), characterized by a general upregulation of cyclin-cyclin-dependent kinase (cyclin-CDK) complexes and checkpoint regulators, including the anaphase-promoting complex/cyclosome (APC/C) ubiquitin-mediated proteolysis system, suggesting a strong proliferative drive. In the Caucasian cohort, the cell cycle pathway was also enriched, though more modestly (FDR = 0.076), with dysregulation predominantly involving CDK1/CDK2 and APC/C.

Cross-ethnic comparisons further highlighted divergent patterns. In the Asian tumors, strong enrichment was observed in xenobiotic metabolism and chemical carcinogenesis pathways (FDR = 2.0 × 10^–6), driven by cytochrome P450 enzymes and DNA adduct repair genes, consistent with processes related to alcohol metabolism and detoxification of environmental carcinogens (Figure S4). Conversely, the Asian-up vs. Caucasian analysis revealed enrichment in endocrine resistance (FDR = 0.15), implicating estrogen receptor signaling and activation of the PI3K/AKT pathway.

Display full size

Figure 6. Gene ontology analysis comprises biological processes, cellular components, and molecular functions. (A) Unique upregulated and downregulated genes in the Asian population. (B) Unique upregulated and downregulated genes in the Caucasian population.

Conclusion

This study reveals significant ethnic-specific differences in gene expression, immune responses, and metabolic pathways between Asian and Caucasian liver cancer patients, highlighting distinct molecular mechanisms underlying HCC pathogenesis. These findings underscore the need for population-specific therapeutic strategies and support the integration of molecular profiling into precision medicine to improve treatment outcomes across diverse patient groups. Future research should focus on translating these molecular insights into targeted immunotherapies and personalized treatment approaches to enhance clinical efficacy and reduce health disparities in liver cancer care.