Multiomics-based analyses of KPNA2 highlight its multiple potentials in hepatocellular carcinoma

Available data and data processing

The clinical features of The Cancer Genome Atlas (TCGA) HCC patients (TCGA-LIHC dataset, named TCGA-HCC in this study, n = 377) and their transcriptome profiling gene expression RNA-seq data (371 primary tumors, 2 recurrent tumors, and 50 paired normal liver specimens from 371 primary HCC patients) were downloaded from the genomic data common (GDC) portal (https://portal.gdc.cancer.gov/). The clinical features of the 377 patients were shown in Table 1. The expressional counts data TMM (trimmed mean of M-values) were normalized for further analyses and the relative expression data of KPNA2 were extracted. The copy number segments after removing germline copy number variation (370 primary tumors from the 377 patients) and the methylation illuminaMethy450 data (377 primary tumors and 50 paired normal liver specimens from the 377 patients) from TCGA-HCC dataset were downloaded from Xena UCSC browser (https://xenabrowser.net/), from which KPNA2′data were further extracted. The splice variants of KPNA2 were searched through the Ensembl genome database (http://mar2017.archive.ensembl.org/). The normalized RNA-seq data of the KPNA2 transcript variants in TCGA-HCC dataset were extracted from the TCGA Pan-Cancer gene expression data which was also downloaded from UCSC Xena browser. The mRNA profiles of the KPNA2 transcripts were analyzed between the paired normal and tumor samples of 50 HCC patients. To obtain the information of KPNA2 protein variants in HCC, Human Protein Atlas (HPA), a comprehensive antibody-based protein atlas that displays expression and localization patterns of proteins in a large portion of human tissues and organs (Uhlen et al., 2005), was searched.

Insights into the copy number variation, mutation, and methylation status of KPNA2 gene in HCC

The correlations of KPNA2 expression with the gene copy number variation (CNV) and methylation status in HCC were investigated by spearman correlation analyses. KPNA2 methylation values between HCC tumors and normal liver tissues were compared by Wilcoxon test. The age-, gender-, and tumor stage- corrected prognostic effects of KPNA2 CNV and KPNA2 methylation level were investigated through ezcox package (Wang et al., 2019) in R3.6.1 software. KPNA2 mutations in TCGA-HCC samples were investigated through the GDC portal (https://portal.gdc.cancer.gov/). Besides, the genes whose mutations might result in KPNA2 expression alteration were identified by the cancer biomarker/target discovery tool muTarget (https://www.mutarget.com/).

The variables (CNV and methylation) with significant prognostic effects as well as age, gender, tumor stage, and KPNA2 expression entered into multivariable cox regression analyses to construct the risk models for HCC OS and DFS. The risk model was set as follows: ${\displaystyle \mathrm{risk score}=\sum _{i=1}^n\beta \left(i\right)\ast exp\left(i\right)}$ where n is the number of selected variables, β (i) is the coefficient of the variable i in the multivariable Cox regression analysis, and exp (i) is the gene expression, tumor stage, and KPNA2 CNV or methylation value of the CpG sites. According to the risk scores of the patients, high- and low- risk groups were distinguished according to the median risk score as the threshold. The survival differences between high- and low- risk patients were visualized via Kaplan-Meier survival analysis with log rank test. In the survival analysis, only the patients with OS time or DFS time >one month (30 days) were included. The effectiveness of the prognostic models was also evaluated in male patients and female patients, younger patients (age <= 60y) and older (age >60y) patients, early-stage (stage I and stage II) patients and late-stage (stage III and stage IV) patients respectively. For all the analyses, p < 0.05 was considered statistically significant.

Exploration of KPNA2 potential functions in HCC

To uncover the potential roles of KPNA2 in HCC, based on the tissue-specific protein-protein interaction (PPI) data from the DifferentialNet database (Basha et al., 2018), its liver-specific PPI network and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment of the interacting partners were investigated through NetworkAnalyst (https://www.networkanalyst.ca/) (Zhou et al., 2019). The DifferentialNet database provides a differential view of the human interactome by integrating current data of experimentally-detected PPIs with data of gene expression across tissues. The differential tissue-specific score of an interaction was computed as the difference between its tissue-specific score and its median score across all tissues. For the network construction, only the PPIs that score in the top 15 percentiles and bottom 15 percentiles were selected. To further explore KPNA2’s relationship with HCC immunoregulation, the microenvironmental infiltration of eight kinds of immune cells including T cells (CD3⁺T cells), CD8⁺ T cells, cytotoxic lymphocytes, NK cells, B lymphocytes, the cells originating from monocytes (monocytic lineage), myeloid dendritic cells, and neutrophils, plus two kinds of stroma cells (endothelial cells and fibroblasts) in TCGA-HCC patients were evaluated with Microenvironment Cell Populations-counter (MCP-counter) method (Becht et al., 2016). Their correlations with KPNA2 expression and HCC prognosis were investigated. In addition, KPNA2 expression in the immune infiltration cells including B cells, CD4⁺ T cells, CD8⁺ T cells, macrophages, and NK cells were compared between HCC and normal liver tissues in TCGA-HCC dataset via GEPIA (http://gepia2021.cancer-pku.cn/) with the one-way analysis of variance (ANOVA) function.

Investigation of KPNA2 association with fatty acid metabolism in HCC

As the liver plays crucial roles in fatty acid metabolism and dysregulations of fatty acids are related to fatty liver diseases and liver cancer (Barbier-Torres et al., 2020; Seo et al., 2020), we speculated that there might be connection between KPNA2 expression and fatty acid metabolism in liver. Herein, the effect of low fat diet on KPNA2 expression was evaluated in liver samples from eight individuals including four interference groups (with low-fat diet) and four control groups (with their habitual diet) in the GEO dataset GSE7117 using GEO2R (https://www.ncbi.nlm.nih.gov/geo/info/geo2r.html), an interactive web tool that allows users to compare two or more groups of samples in a GEO series. Furthermore, in consideration of the relationship of fatty acid dysregulation with HCC (Che et al., 2019; Montagner, Cam & Guillou, 2019), the correlations of KPNA2 with the genes of fatty acid metabolism pathway in TCGA-HCC primary tumors in KEGG database were investigated through Spearman correlation analysis with Hmisc package (https://CRAN.R-project.org/package=Hmisc) in R software. The potential implication of KPNA2 in fatty acid biosynthesis was evaluated from the purity-corrected Spearman correlations of KPNA2 with the 11 genes whose products were annotated to fatty acid synthase activity (http://amigo.geneontology.org/amigo/term/GO:0004312) via Tumor IMmune Estimation Resource (TIMER) (https://cistrome.shinyapps.io/timer/) and p < 0.05 was considered significant.

The correlations of KPNA2 with other genes in HCC were also investigated in the TCGA-HCC database. Gene set enrichment analysis (GSEA) was performed with the top positively and negatively KPNA2-correlated genes (Spearman’s correlation coefficient >0.3 or <−0.3) through a web-based gene set analysis toolkit (http://www.webgestalt.org/) to further uncover KPNA2-associated pathways and to see if fatty acid metabolism was in the enriched pathways.

Characterization of the variants of KPNA2 transcripts and proteins, and plasma testing

The differential expression of KPNA2 transcript variants between paired HCC and normal controls from the 50 HCC patients in the TCGA-HCC dataset was analyzed through paired sample Wilcoxon rank test and p < 0.05 was considered statistically significant. To further uncover the potential roles of KPNA2 transcripts in HCC development, the differentially expressed transcripts were evaluated to investigate their correlations with alpha fetoprotein (AFP) and albumin (ALB) gene expression of HCC samples by Spearman correlation analysis. Moreover, multivariable Cox regression analysis of the differentially expressed transcripts, with gender, age, and tumor stage as the covariates was performed to investigate the transcripts’ independent prognostic values.

For KPNA2 protein variants, the immunohistochemical staining patterns with the antibodies against different immunogenic fractions of KPNA2 in HCC and normal liver tissues were extracted from the HPA database and analyzed. The protein level of KPNA2 (full-length) in plasma was examined by enzyme linked immunosorbent assay (ELISA). The work involving the plasma specimens was reviewed and approved by the Ethics Committee of Xiangya Hospital of Central South University (approval number: 201801002). All samples were collected with informed consent in accordance with the Declaration of Helsinki. A total of 92 plasma samples with clinical information (51 samples from primary HCC patients and 41 samples from healthy individuals) were enrolled in the study (Table S1). The plasma samples were obtained between July 2018 and December 2018 from Xiangya Hospital of Central South University (Changsha, China) at the time of diagnosis before any therapy. KPNA2 concentration was measured in duplicate using commercially available ELISA kits (abx250558, Abbexa Ltd, UK) according to the manufacturer’s instructions. Absorbance was read at 450 nm with wavelength correction set at 540 nm using an ELISA plate reader (BioRad, CA, USA). Wilcoxon test was used for two-group comparison. ROC analysis was performed to investigate the diagnostic power of the factor with pROC package (Robin et al., 2011) in R software. For the analyses, p < 0.05 was considered to be statistically significant.

Genetic and epigenetic regulations of KPNA2 and their prognostic value in HCC

By correlation analysis, a significant positive correlation of KPNA2 CNV with KPNA2 expression was shown in Fig. 1A (R = 0.43, p < 0.001). Four CpG sites locating in KPNA2 promoter including cg23206777, cg22429852, cg21018429, and cg21820889 presented negative correlations (p < 0.05, Figs. 1B–1E) with KPNA2 expression respectively in HCC, indicating their potential roles in regulating KPNA2 expression. All of the 13 CpG methylation sites in KPNA2 with location information, KPNA2 expression-correlation values, and p values were showed in Table S2. By comparison analysis of the four sites methylation between HCC tumors and the control groups, lower methylation of cg22429852 and higher methylation of cg23206777 was found in HCC than in the normal, while the other two sites methylation showed no significant difference (Figs. 1F–1I). It is noteworthy that the negative relationship of the four sites methylation with KPNA2 expression was conducted in the HCC data without normal liver data. Thus we deduced that certain methylation sites might be affected by negative feedback regulation under physiological and pathological processes, and cg22429852 might be the most critical methylation site for KPNA2 upregulation. For KPNA2 gene mutation inquiry in the GDC portal (https://portal.gdc.cancer.gov/), only four patients in TCGA-HCC dataset were found to have KPNA2 mutations and no significant KPNA2 expression difference was revealed between the HCC samples with and without KPNA2 mutation (p = 0.69, Fig. S1). As for other gene mutation effects on KPNA2 expression, by muTarget analysis with mutation prevalence set at 3%, five mutated genes appeared to be related to KPNA2 expression. As shown in Fig. S2, the mutations of TP53 (p < 0.01), DNAH10 (p < 0.01), TSC2 (p < 0.01), and RB1 (p < 0.01) were associated with KPNA2 upregulation while KPNA2 was more reduced in HCC samples with BAP1 mutations (p < 0.01) than in those without BAP1 mutations.

When adjusted with age, gender, and tumor stage, the methylation status of two CpG sites of cg23206777 (p = 0.026) and cg17985418 (p = 0.00258) was shown to have a favorable or an unfavorable prognostic effect on HCC OS respectively, while KPNA2 CNV and other KPNA2 CpG sites presented no significant correlation with HCC OS (Table S3). For evaluation of the CNV and CpG associations with HCC DFS (Table S4), KPNA2 CNV (p = 0.002) and the methylation value of cg14898140 (p = 0.006) were shown to have unfavorable prognostic effects on HCC DFS independent of age, gender, and tumor stage while no significant correlation of other CpG sites was shown.

Additionally, the risk models for HCC OS and DFS were constructed by multivariable Cox regression analysis individually. As shown in Fig. 2A, tumor stage, KPNA2 expression, and cg17985418 methylation were shown to be independent prognostic factors for HCC OS and the risk model was as follows: ${\displaystyle \mathrm{Risk score}\left(OS\right)=0.40271\ast \mathrm{tumor stage}+5.49197\ast cg17985418+0.55453\ast KPNA2\mathrm{expression}.}$ There was a shorter OS for the HCC patients with high-risk scores than the ones with low-risk scores (Fig. 2B). By ROC analysis, an area under the curve (AUC) of 0.814 for 1-year survival, 0.788 for 3-year, and 0.733 for 5-year survival was shown respectively in Fig. 2C, indicating the effectiveness of the risk model for predicting HCC patients OS. To compare the predictive efficiency for HCC OS between the constitutive risk model and KPNA2 expression alone, the predictive efficiency of KPNA2 for HCC OS was calculated. By comparison, although there was no significant difference in their AUCs of 1-year OS (0.814 vs.0.779, p = 0.113), the AUCs of 3-year (0.788 vs. 0.707, p = 0.003) and 5-year OS (0.733 vs. 0.665, p = 0.034) of the constitutive model were significantly higher than those of KPNA2 expression alone (Fig. S3). It indicated the superiority of a multi-variable risk score model in the prediction of HCC prognosis, especially for the long survival status of the patients.

With regards to HCC DFS (Fig. 3), tumor stage, KPNA2 CNV, and cg14898140 methylation were presented to be independent prognostic factors for HCC patients, while the prognostic value of KPNA2 expression was not significant (Fig. 3A). The risk model was established as: ${\displaystyle \mathrm{Risk score}\left(DFS\right)=0.489\ast \mathrm{tumor stage}+0.644\ast KPNA{2}_{CNV}+11.3\ast cg14898140.}$ The HCC patients with high-risk scores had a shorter disease-free time than the patients with low-risk scores (Fig. 3B). From the ROC analysis, the AUCs for 1-year, 3-year, and 5-year HCC DFS was 0.721, 0.670, and 0.666 respectively (Fig. 3C).

Furthermore, the unfavorable prognostic effects of OS high-risk scores were shown in male HCC patients (p < 0.0001, Fig. 4A), female HCC patients (p <0.0001, Fig. 4B), young HCC patients (p < 0.0001, Fig. 4C), elderly HCC patients (p <0.0001, Fig. 4D), early-stage HCC patients (p < 0.0001, Fig. 4E), and late-stage HCC patients (p = 0.00026, Fig. 4F) individually, indicating that the OS risk score was an independent prognostic factor. Likewise, the HCC patients with DFS high-risk scores had shorter disease-free time than those with the low-risk scores in different gender groups (p < 0.01, Figs. 4G–4H), different age groups (p <0.0001, Figs. 4I–4J), and the early-stage group (p <0.0001, Fig. 4K). In the late-stage group (Fig. 4L), though only one of the 77 late-stage patients was grouped into the low-risk group, it still showed a longer disease-free time than most of the patients with high-risk scores. It indicated that the prognostic effect of the risk scores for HCC DFS was independent of gender, age, and tumor stage.

Immunoregulatory roles of KPNA2 in HCC

Via NetworkAnalyst, the liver-specific PPI network of KPNA2 was investigated and 68 proteins were found to have interactions with KPNA2 (Fig. 5A). By KEGG pathway analysis, the interacting proteins were enriched into five significant terms including viral carcinogenesis, focal adhesion, leukocyte transendothelial migration, PI3K-Akt signaling pathway, and adherens junction (Fig. 5B), suggesting the association of KPNA2 with HCC immunology.

Further study revealed that KPNA2 expression was positively correlated with infiltration of five kinds of immune cells in HCC including monocytic lineage cells (R = 0.35, p = 4e −12, Fig. 5C), T cells (R = 0.23, p = 1e −05, Fig. 5D), B lineage cells (R = 0.14, p = 0.0071, Fig. 5E), CD8⁺ T cells (R = 0.14, p = 0.006, Fig. 5F), and myeloid dendritic cells (R = 0.13, p = 0.014, Fig. 5G). While no significant correlation of KPNA2 expression with cytotoxic lymphocytes (p >0.05, Fig. 5H), neutrophils (p >0.05, Fig. 5I), and NK cells (p >0.05, Fig. 5J) was indicated. As to the two kinds of stroma cells, a negative relationship between endothelial cell infiltration and KPNA2 expression was presented (R =  − 0.19, p = 0.00033, Fig. 5K), while no significant correlation was shown between KPNA2 expression and the level of fibroblasts (p > 0.05, Fig. 5L). The correlation of KPNA2 with immune cell infiltrations in early (TNM stage I) and middle-late stage (TNM stage II/III/IV) were evaluated separately (Fig. S4). Interestingly, the Spearman correlations of the five types of immune infiltration cells with KPNA2 appeared statistically significant or more significant in HCC middle-late stage than in the early stage. It indicated that KPNA2′  connection with HCC immune infiltration is becoming increasingly apparent during HCC progress. By Cox regression analysis, among the KPNA2-correlated immune cells and stroma cells, monocytic lineage cells and myeloid dendritic cells were associated with HCC OS (p < 0.05, Fig. 6A), while none presented prognostic significance for HCC DFS (p > 0.05, Fig. 6B).

GEPIA2021 analysis resulted that among the five kinds of immune cells including B cells, CD4⁺ T cells, CD8⁺ T cells, macrophages, and NK cells, KPNA2 expression was lower in CD4⁺T cells (p = 0.02) and macrophages (p < 0.01) in HCC tumors in comparison with that in normal liver tissues, as shown in Fig. 7. Considering its important functions in normal cells (Goldfarb et al., 2004), KPNA2 dysregulation might connect to the dysfunction of the immune cells in HCC.

Associations between KPNA2 expression and fatty acid metabolism

As shown in Fig. 8A, following a low-fat hypocaloric diet for 8 weeks, there was a significant decrease of KPNA2 expression (probes: 201088_at and 211762_s_at from GPL570) in the four obese liver tissues (p < 0.05), indicating a negative effect of low-fat diet on KPNA2 expression in liver.

For the association of KPNA2 expression with fatty acid metabolism related genes in HCC, among the 42 genes in KEGG fatty acid metabolism pathway, 88.1% of them (37/42) were shown to be significantly and negatively correlated with KPNA2 expression (Table S5). The top ten negatively correlated genes were shown in Fig. 8B. When the KPNA2-negatively-correlated genes in fatty acid metabolism pathway were applied to functional enrichment analyses, multiple fatty acids metabolic processes were revealed including fatty acid degradation and beta-oxidation on the top (Fig. 8C). As the down-regulation of fatty acid degradation and beta-oxidation might lead to the increase of fatty acid in the liver, a positive correlation between KPNA2 expression and fatty acid level in HCC was deducible. In fact, as shown in Table 2, nine of the eleven genes (81.8%) related to fatty acid synthase activity were shown to be significantly positively correlated with KPNA2 expression in HCC. All of these results were consistent with the association between low fat diet and KPNA2 decrease in the liver from the beginning analysis.

In addition, the correlations of KPNA2 with other genes in HCC were investigated in the TCGA-HCC database. Its top positively and negatively correlated genes with Spearman’s correlation coefficient above 0.3 or below −0.3 were extracted (File S3). In KEGG GSEA result, KPNA2-correlated genes were most positively enriched in cell cycle pathway while negatively associated with a variety of metabolic pathways including fatty acid degradation (Fig. S5), consistent with KPNA2 negative involvement in fatty acid degradation in HCC.

The landscape of KPNA2 transcripts and protein variants in HCC

KPNA2 is encoded by the KPNA2 gene on the chromosome 17 which contains 11 exons spanning approximately 11 kb in the human genome. The full-length KPNA2 protein (UniProt_P52292) consists of 529 amino acids and weighs around 58 kDa, comprising an N-terminal hydrophilic importin β-binding domain (IBB, from 2nd to 60th amino acid), a central hydrophobic region which binds the cargo’s NLS (from 142nd to 403rd amino acid), and a short acidic C-terminus with no reported function. KPNA2 has seven splice transcript variants as shown in Fig. S6A). Both of ENST00000330459 (amount percentage range: 25.87%–100%) and ENST00000537025 (amount percentage range: 0−1.16%) code for the full length of KPNA2 protein, the former of which is the predominant transcript among all the variants (Fig. S6B)). ENST00000579754 and ENST00000584026 are derived from alternative splicing discarding the rear half of exon 5 and the posterior exons, resulting in two short protein fragments of N-terminal 143 and 134 amino acids respectively, which retain only the IBB domains. ENST00000583392 and ENST00000582898 are two non-coding RNAs (ncRNAs).

Paired samples Wilcoxon tests revealed that the two long coding transcripts (ENST00000330459 and ENST00000537025) (p < 0.05, Figs. 9A–9B) and the two ncRNAs (ENST00000583392 and ENST00000582898) (p < 0.05 Figs. 9F–9G) were higher expressed in HCC than in the normal controls. No significant expression difference of the two short protein-coding mRNA (ENST00000579754 and ENST00000584026) (p > 0.05, Figs. 9C–9D) and the nonsense mediated decay ENST00000583269 (p > 0.05, Fig. 9E) was shown between the paired normal liver and HCC samples. Among the four differentially expressed transcripts (Table S6), ENST00000330459 (r = 0.223, p = 1.583E-05) and ENST00000582898 (r = 0.189, p = 2.538E-04) were shown to be positively correlated with AFP expression. In contrast, the negative correlations of ENST00000330459 (r =  − 0.280, p = 4.384E-08) and ENST00000583392 (r =  − 0.147, p = 0.005) with ALB expression were shown. Previous studies supported that serum AFP level was positively associated with HCC nuclear KPNA2 expression (Jiang et al., 2014; Yoshitake et al., 2011). It substantiated KPNA2 involvement in liver dysfunction and HCC development.

When the four differentially expressed transcripts were applied to multivariable Cox regression analysis with gender, age, and tumor stage as variates, tumor stage and the full-length KPNA2 coding transcript ENST00000330459 independently showed unfavorable prognostic significance for HCC OS (Fig. 10A) and DFS (Fig. 10B). While the other full-length KPNA2 coding transcript ENST00000537025, which made up a tiny constitution of the transcript variants, showed a significant favorable prognostic effect on HCC OS (HR = 0.91, p < 0.05, Fig. 10A).

Although KPNA2 was reported to be upregulated in HCC tissues in a recent proteomic study (Jiang et al., 2019), the distribution of KPNA2 protein variants was unclear. Herein, the immunohistochemical staining patterns of KPNA2 variants between HCC and normal liver samples were extracted from HPA database. Two types of KPNA2 antibodies against different immunogenic fragments of KPNA2 resulted in two different staining patterns. Raised against the amino acids 480-529 in KPNA2 C-terminal, CAB015460 antibody reacted solely to the full length KPNA2 (UniProt_P52292), while HPA041270 antibody which was against the immunogenic sequence in N-terminal of KPNA2 detected all the KPNA2 protein variants including the full length KPNA2 and the other two short protein fragments (UniProt_J3KS65 and J3QLL0) which keep N-terminal domain. KPNA2 was stained by CAB015460 moderately in the cytoplasm/membrane with no nuclear staining in all of the three normal liver samples, while showed nuclear staining in all of the six HCC samples which ranged from weak to strong with rare cytoplasmic or membranous staining observed. The representative staining images were showed in Figs. 11A and 11B respectively. With the other antibody HPA041270 which indiscriminately detected all the KPNA2 protein variants, the three normal liver samples were extensively stained in both cytoplasm/membrane and nuclei (Fig. 11C). The staining differences of KPNA2 between using the two antibodies in the normal liver samples indicated that KPNA2 in the normal nuclei might be the short KPNA2 protein fragments which could only bind with antibody HPA041270, not CAB015460. The short fragments of KPNA2 might play physiological roles in normal liver cells. When HPA041270 was applied to seven HCC specimens, KPNA2 remained positive staining in cytoplasm, membrane, and nuclei in most samples. In the three HCC samples, KPNA2 was more intensive in HCC nuclei (Fig. 11D) than in normal nuclei (Fig. 11C), indicating that at least a part of KPNA2 variants were translocated to nuclei if not all, and it was most likely the full-length KPNA2 (UniProt_P52292).

We tested the full-length KPNA2 protein level in plasma. As shown in Fig. 12A, plasma KPNA2 (UniProt_P52292) was shown to be lower in HCC patients than in normal individuals (p < 0.001), which might be partly due to the translocation of KPNA2 to nuclei in HCC as described above. By ROC analysis, the AUC for plasma KPNA2 was 0.787 in discriminating HCC from normal controls (Fig. 12B). With an optimal cutoff value of 7.923 ng/ml, 65.9% (sensitivity) of the 51 HCC patients and 98.0% (specificity) of the 41 normal controls could be accurately discriminated, indicating the diagnostic power of plasma KPNA2 in HCC.