Pan-Cancer Profiling and Digital Pathology Analysis Reveal Negative Prognostic Biomarker ZPR1 Associated with Immune Infiltration and Treatment Response in Hepatocellular Carcinoma

Introduction

Malignant tumor is the second leading cause of death in the world and poses a severe threat to the life and health of all humankind. By 2020, hepatocellular carcinoma (HCC) was the sixth most frequent cancer and the third leading cause of cancer death globally.¹ The 5-year relative survival rate of HCC in the United States is 20%.² With the rapid development of precision tumor therapy, it is crucial to develop more sensitive and effective biomarkers to guide the selection of optimal clinical treatment plans effectively. Pan-cancer research combined with multiple omics analysis in specific tumors has been used to identify tumor markers and signaling pathways and help us understand the molecular mechanism of biomarkers in cancer development.

ZPR1 zinc finger (ZPR1), also known as ZNF259, is ubiquitously expressed and evolutionarily conserved in mammalian cells. ZPR1 was first found to bind to the tyrosine kinase domain of cytoplasmic epidermal growth factor receptor (EGFR) through the zinc finger domain.³ ZPR1 has been shown to interbinds with eukaryotic translation elongation factor 1A (eEF1A) to form the ZPR1-eEF1A complex, which is essential for cell replication, cell cycle progression, and survival.^4–6 Increasing evidence has demonstrated that ZPR1 is associated with spinal muscular atrophy (SMA),^7–15 type 2 diabetes,¹⁶ non-alcoholic fatty liver disease (NAFLD),¹⁷ metabolic syndrome,¹⁸ hypertriglyceridemia,^17–23 and growth restriction.²⁴ In tumors, ZPR1 serves a specific role in carcinogenesis and development. Previous research has reported that ZPR1 is upregulated and promotes the malignant progression of breast cancer.²⁵ In contrast, Shan et al found that ZPR1 was downregulated and inhibited the proliferation and invasion of lung cancer.²⁶ To date, the expression level, prognostic value, and potential molecular functions of ZPR1 in pan-cancer and HCC are not clear.

Based on the integrated data analyzed from TCGA, GTEx, GENT2, UALCAN, and Kaplan–Meier databases, as well as the experiment validated by IHC assay, our research applied the comprehensive bioinformatics methods to explore the ZPR1 expression in pan-cancer and its correlation with prognostic significance and immune infiltrate in HCC, aiming to reveal the ZPR1 biological function in malignant progression and its potential as a novel independent prognostic predictor for HCC.

Materials and Methods

Data Source Collection

The pan-cancer gene expression profiles (HT Seq-TPM) and clinical data of 33 tumors based on Genotype-Tissue Expression (GTEx) and The Cancer Genome Atlas (TCGA) databases were acquired from UCSC Xena platform (https://xena.ucsc.edu). Detailed information on 33 cancers was presented in Supplementary Table 1.

Pan-Cancer Expression and Clinical Correlation Analysis of ZPR1

The ZPR1 pan-cancer mRNA expression was then extracted and visualized using R package limma. UALCAN-CPTAC platform (http://ualcan.path.uab.edu/) was utilized to retrieve ZPR1 protein expression difference in pan-cancer. The ZPR1 expression in different tissues/cells and subcellular localization were retrieved from Human Protein Atlas platform (https://www.proteinatlas.org/). We investigated DNA methylation profile of ZPR1 by retrieving the SMART platform (http://www.bioinfo-zs.com/smartapp/).We looked into the association between ZPR1 expression and pan-cancer clinicopathological features (including age, gender, tumor status, and pathological stage). The values of the pan-cancer stemness index,²⁷ TMB (tumor mutation burden),²⁸ and MSI (microsatellite instability)²⁹ were obtained from previous reports. Spearman correlation analysis of ZPR1 expression with pan-cancer stemness index, TMB, and MSI values was performed.

Expression Difference Analysis and Prognostic Analysis of ZPR1 in HCC

The GENT2 database (http://gent2.appex.kr/gent2/) collects the expression profiles of multiple GEO datasets based on the GPL96 and GPL570 platforms. The GENT2 tool was adopted to visualize the ZPR1 mRNA expression difference in the GEO database. The chi-square test was utilized to look into the correlation between ZPR1 and pathological parameters. Kaplan-Meier (K-M) survival analysis for progression-free survival (PFS) and overall survival (OS) were carried out on the TCGA HCC samples and verified by using Kaplan-Meier plotter database. The independent prognostic evaluation of ZPR1 expression and clinicopathological characteristics were performed using cox regression analyses. R package rms was applied to construct a nomogram to predict the 1-year and 3-year OS of HCC patients based on ZPR1 expression.

Interaction Factors and Functional Enrichment Prediction of ZPR1 in HCC

With the criteria of p < 0.001 and Pearson correlation coefficient > 0.3, co-expression analysis was adopted to determine the interaction factors associated with ZPR1. The harmonizome database (https://maayanlab.cloud/Harmonizome/), collecting the transcription factor-target gene prediction information from the CHEA, ENCODE, and JASPAR databases, was adopted to predict the target genes of the potential transcription factor ZPR1. The correlation scatter plots between ZPR1 and potential target genes were undertaken. Kyoto Encyclopedia of Genes and Genomes (KEGG), Gene Ontology (GO), and Gene set enrichment analysis (GSEA) were performed to investigate the potential molecular functions, biological processes, and tumor-related pathways of ZPR1. We screened the ZPR1-associated differentially expressed genes (DEGs) with a threshold of |log Fold Change| > 0.5 and FDR < 0.05. The KEGG and GO analyses of ZPR1-related DEGs were carried out with R package ClusterProfiler. GSEA was conducted with the GSEA tool (https://www.gsea-msigdb.org/gsea).

Assessment of Immune Characteristics in HCC

R package ESTIMATE was adopted to evaluate tumor microenvironment (TME) scores for TCGA-HCC samples. The scores of 22 immune cell infiltration in HCC samples were analyzed by the robust method CIBERSORT. Furthermore, the multiple algorithms quantified immune cell infiltration scores of the TCGA-HCC samples were acquired from TIME database (http://cistrome.shinyapps.io/timer/). Pearson correlation analysis with a significance threshold of p < 0.001 was performed between ZPR1 expression and 79 previously reported immune checkpoint genes (ICG).³⁰ Using the single-sample gene set enrichment analysis (ssGSEA) method with R package GSVA, the degree of immune cell infiltration and the activity level of immune-related functions were evaluated in HCC samples. The TISCH2 database (http://tisch.comp-genomics.org/) was utilized to investigate the correlation between the ZPR1 and TME-associated cells at the single-cell level.

Treatment Response Prediction

The significance of ZPR1 in predicting clinical response to HCC therapy was further evaluated. We then predicted the ZPR1 sensitivity to immune checkpoint blockade therapy using the TIDE tools (http://tide.dfci.harvard.edu/).³¹ TIDE score, including T cell dysfunction score and exclusion scores, was used to evaluate the ZPR1 resistance to immune checkpoint blockade (ICB) treatment. The dysfunction score predicts T cell exhaustion status and tumor immune evasion signature. T cell exclusion score predicts immune escape by evaluating immunosuppressive cells. Furthermore, we compared the ZPR1 susceptibility to common HCC therapeutic agents using R package pRRophetic (p < 0.01).

Patient Specimen Collection and Immunohistochemistry (IHC)

Two tissue microarrays (TMA) were constructed from the paraffin-embedded tissue blocks with 70 HCC tissues and 30 adjacent liver tissues collected from the Department of Pathology, Liaoning Cancer Hospital & Institute. All cases were HCC specimens confirmed by pathology. The IHC protocol was conducted as previously described.³² IHC experiment was conducted following the manufacturer’s instructions with ZSGB immunohistochemical staining kit (PV-9000, ZSGB-BIO, Beijing, China). The TMAs were incubated overnight with commercial antibodies at 4 °C. Detailed information about antibodies is provided in Table 1. PD-L1 IHC 22C3 pharmDx kit was applied to evaluate the PD-L1 expression on the DAKO Link-48 Autostainer system.

Table 1 The Detailed Information About antibodies in This Research

Digital Quantitative IHC Evaluation by QuPath Software

The TMAs were scanned by the digital whole slide scanner (Hamamatsu NanoZoomer) and quantitatively analyzed using QuPath software (version 0.4.1, https://qupath.github.io/).³³ The histochemistry scores (H-score) calculated by QuPath software were adopted for subsequent analyses.³⁴ In short, the plug-in TMA dearrayer was used to identify each tissue point, and then plug-in cell detection was applied to perform cell identification and quantification. An object classifier is trained to classify individual cells as tumor or stromal. Finally, the H-score was obtained by setting the identification intensity feature according to different antibody staining distributions (Table 1). H-score intensity thresholds were set as 0.1, 0.2, and 0.3, respectively. Combined positive score (CPS) and tumor cell proportion score (TPS) are commonly used clinical methods for evaluating CD274 (PD-L1). With 1% as the positive threshold, the correlation between ZPR1 and CPS, TPS in HCC was analyzed.³⁵

Statistical Analysis

R software (version 4.1.0) was utilized to perform all bioinformatics statistical analyses. In IHC, the scatter plot of the ZPR1 expression difference was analyzed using the Wilcox signed-rank test. The association of ZPR1 with immune cell makers was evaluated by Spearman correlation analysis. Unless otherwise indicated, statistical significance was defined as p < 0.05.

Results

The ZPR1 Was Differentially Expressed in Pan-Cancer

As a first step, we investigated the ZPR1 expression difference in pan-cancer. Considering that some cancer types do not have corresponding normal tissues, we integrated the normal tissue expression profile of GTEx database on the basis of TCGA database and systematically explored the ZPR1 mRNA expression differences in 33 cancer types. We found that ZPR1 mRNA expression was up-regulated in most tumors, including BRCA, colorectal cancer (CRC), KIRC, GBM, HNSC, lung cancer (LC), PAAD, STAD, and HCC (Figure 1A and B). Furthermore, the tumor tissues paired analysis of TCGA samples suggested that ZPR1 mRNA expression was also significantly increased in CRC, HNSC, KIRC, LC, STAD, PAAD, and HCC tissues (Figure 1C). In the UCLACN-CPTAC database, a significant up-regulation of ZPR1 protein expression was observed in CRC, KIRC, LC, HNSC, GBM, PAAD, and HCC tissues (Figure 1D). The results retrieved from the HPA website demonstrated that ZPR1 was widely expressed in various human tissues and cell types (Figure 1E and F). ZPR1 was highly expressed in the liver, testis, skeletal muscle and bone marrow tissues, and myeloid and lymphoid cell lines. Immunofluorescence studies showed that ZPR1 was scattered localized in the nucleus and cytoplasm, proposing that ZPR1 may be localized to the vesicles (Figure 1G).

Figure 1 The ZPR1 expression in pan-cancer. (A) ZPR1 expression in pan-cancer tissues from TCGA database. (B) ZPR1 expression in pan-cancer tissues based on TCGA and GTEx databases. (C) ZPR1 expression in paired tumor and adjacent normal tissues from TCGA database. (D) ZPR1 protein expression in pan-cancer tissues from the UALCAN-CPTAC database. ZPR1 mRNA expression in various tissues (E) and cancer cell lines (F) from the HPA database. (G) Cellular localization of ZPR1 from the HPA database. *p < 0.05; **p < 0.01; and ***p < 0.001. Abbreviation: Ns, no significance.

The DNA Methylation Profile of ZPR1 in Pan-Cancer

Further, we investigated the ZPR1 methylation level in attempt to elucidate the reasons for the high ZPR1 expression in most tumors. As shown in Figure 2A and B, 14 methylation probes associated with ZPR1 were identified on chromosome 11, including seven island regions (cg01218534, cg23953395, cg22654495, cg06541888, cg27192902, cg19345149, cg17490921), three S Shore region (cg12395125, cg20570279, cg05447673), one S Shelf region (cg14815609), one N Shore region (cg16408784), one N Shelf region (cg26548001), and an unidentified region (cg06595719). We then explored the ZPR1 methylation level in pan-cancer. In BLCA, BRCA, SARC, LC, PRAD, and HCC, the ZPR1 methylation level was significantly lower than in normal tissues (Figure 2C). The above results indicate the up-regulation of ZPR1 in various tumors may be due to its hypo-methylation status.

Figure 2 DNA methylation analysis of ZPR1. (A) Chromosomal distribution of the 14 methylation probes associated with ZPR1 on chromosome 11. (B)Detailed location of 14 probes in CpG Island. (C) ZPR1 methylation levels in tumor and normal samples cross different types of cancers. p > 0.05; *p <= 0.05; **p <= 0.01; ***p <= 0.001; ****p <= 0.0001. Abbreviation: Ns, no significance.

ZPR1 Correlated with Pan-Cancer Clinicopathological Characteristics, TMB, MSI, and Stemness Values

We then evaluated the clinical correlation of ZPR1 in patients with different tumor types. We observed that ZPR1 expression was closely related to age in BRCA, LUAD, KIRP, STAD, READ, THCA, TGCT, and UVM (Figure 3A). In ACC, LGG, PRAD, and THYM, ZPR1 expression was closely correlated with tumor status (Figure 3B). Furthermore, we found that the ZPR1 expression in HCC was differentially expressed in gender and pathological stage (p < 0.05, Figure 3C and D). In addition, ZPR1 expression was closely correlated with TMB, MSI, and stemness values in HCC (Figure 3E–G). The pan-cancer clinical correlation results indicated that ZPR1 might serve a critical role in HCC development. Thus, the relationship between ZPR1 expression and prognosis of HCC was further investigated.

Figure 3 Clinical correlation analysis between the ZPR1 expression and patient’s age (A), gender (B), tumor status (C), pathological stage (D), TMB (E), MSI (F), and stemness (G) in pan-cancer. *p < 0.05; **p < 0.01; ***p < 0.001.

High ZPR1 Expression Was Strongly Correlated with the Unfavorable Outcome of HCC

The up-regulation of ZPR1 in HCC retrieved from the GENT2 database (Figure 4A) was in accordance with the above TCGA differential expression analysis (Figure 1A and B). Then, we looked into the correlation between the ZPR1 expression and HCC clinicopathological features using chi-square test. As demonstrated in Figure 4B, ZPR1 expression was remarkably correlated with pathological stage, T stage, and tumor grade of HCC. K-M survival analyses indicated that high ZPR1 expression group had shorter OS (p < 0.001, Figure 4C). The median OS of the ZPR1 high expression group was 2.2932 years, compared to 5.8384 years for the low expression group. Consistent results were observed for PFS (1.3836 vs 2.4438 years, p = 0.019; Supplementary Figures 1A). Additionally, the negative connection between high ZPR1 expression and favorable prognosis in HCC patients was verified using the Kaplan-Meier plotter database (Supplementary Figures 1B and C). Furthermore, subgroup survival analysis stratified by different clinicopathological features indicated that high expression of ZPR1 was closely associated with unfavorable outcomes (Table 2). The above results show that the ZPR1 expression is negatively related to the survival of HCC patients, and therefore inhibiting ZPR1 expression may be a useful intervention strategy.

Table 2 The Prognostic Values of ZPR1 (RNAseq ID: 8882) Expression Stratified by Different Clinicopathologic Features from Kaplan-Meier Plotter Databases

Figure 4 ZPR1 could be considered as an independent predictor for HCC patient survival. (A) The difference of ZPR1 expression in HCC retrieved from the GENT2 database. (B) The heatmap showing the clinical correlation analysis of ZPR1 expression in HCC by the Chi-square test. (C) The overall survival analysis of ZPR1 in TCGA HCC samples. Univariate Cox analysis (D) and multivariate Cox analysis (E) of ZPR1 expression and clinicopathological characteristics. (F) The ROC curves of 1-year and 3-year OS of ZPR1. (G) Nomogram of ZPR1 expression and clinicopathological features. (H) The calibration curves of the nomogram for 1-year and 3-year OS.

ZPR1 Could Be Considered as an Independent Predictor for HCC Patients

We next adopted Cox regression analyses to assess the potential prognostic variables impacting the prognosis of HCC patients. Univariate Cox regression analysis showed that the TNM stage (HR = 1.680) and ZPR1 expression (HR = 1.745) could be candidate predictors of OS for HCC patients’ survival (Figure 4D). Furtherly, multivariate Cox regression analysis indicated that high expression of ZPR1 could be an independent predictor of adverse prognosis in HCC patients (HR = 1.530, Figure 4E). The ROC curve analysis of 1-year (AUC = 0.673) and 3-year (AUC = 0.622) showed the high accuracy of ZPR1 prognostic efficiency (Figure 4F). Furthermore, we establish a nomogram based on ZPR1 expression and multiple clinicopathological characteristics to predict 1-year and 3-year OS incidence of HCC patients, with a higher total score indicating worse survival (Figure 4G). The calibration curves demonstrated that the nomogram-predicted OS outcomes based on the ZPR1 expression were in high agreement with the observed probabilities (Figure 4H).

Potential Target Prediction of ZPR1 in HCC

The co-expression circos plot presented the 11 candidate genes associated with ZPR1 expression (Figure 5A). The top 6 positively correlated genes are DDX10, BUD13, ZW10, DPAGT1, HMBS, and DCUN1D5, and the top 5 negatively correlated genes are LCAT, ANG, MT-RNR1, RNASE4, and PGLYRP2. Among them, the correlation of ZPR1 with DDX10 (R = 0.742) and BUD13 (R = 0.729) was over 0.7. Accumulating studies have demonstrated that ZPR1, containing two zinc-finger domains, may act as a potential transcription factor to regulate downstream target gene transcription.^10,14,36 The prediction information of CHEA, ENCODE, and JASPAR databases indicated that the target gene of ZPR1 might be CREB1 (R = 0.44, Figure 5B and C). Additionally, the CHEA and ENCODE databases also suggest that other genes (CTCF, KDM5A, KDM5B, MYC, HNF4A, FLI1) were the candidate target genes of ZPR1.

Figure 5 The functional prediction of the ZPR1 gene in HCC. (A) Co-expression analysis of ZPR1. (B) Target gene prediction through CHEA, ENCODE, and JASPAR databases. (C) Pearson correlation analysis between ZPR1 and predicted target gene. (D) Differential expression genes in high- and low-ZPR1 expression groups. (E and F) GO and KEGG analyses showing biological processes and molecular functions involved in ZPR1-related DEGs. (G) GSEA analysis showing cancer- or metabolism-related pathways enriched in high-ZPR1 groups.

Functional Enrichment Analysis of ZPR1 in HCC

The top 50 DEGs associated with ZPR1 were then identified (Figure 5D). Among them, the significantly elevated expression level of SEZ6L, IL11, CDH17, and CD177 was observed in the high ZPR1 expression group. The GO and KEGG analysis revealed that ZPR1-associated DEGs were localized in microtubules, spindle, and nuclear chromatin regions and implicated in various biological processes associated with cell cycle and DNA replication segregation, nuclear division, and tubulin binding (Figure 5E and F). Interestingly, these ZPR1-related DEGs were significantly enriched in transcriptional co-regulatory activity molecular function and nucleocytoplasmic transport pathway. GSEA analysis proposed that ZPR1 high-expression group was involved in neurological disorders (like Alzheimer’s disease and Huntington’s disease), lipid metabolism (including glycerophospholipid metabolism, sphingolipid metabolism, and glycerolipid metabolism), DNA damage repair pathways³⁷ (including homologous recombination, nucleotide excision repair, base excision repair, and mismatch repair), and several tumor development-related pathways (ERBB, MAPK, MTOR, NOTCH, and WNT signaling) (Figure 5G).

Immune-Related Analysis of ZPR1 in HCC

We next investigated the correlation between ZPR1 and TME. As illustrated in Figure 6A, high expression of ZPR1 was negatively associated with the TME estimation and immune scores, indicating that ZPR1 high expression group had low tumor purity. The CIBERSORT method suggested that ZPR1 expression was positively correlated with the degree of M0 macrophage infiltration (p = 0.008, Figure 6B). The infiltration fraction of immune cells quantified by CIBERSORT, CIBERSORT-ABS, MCPCOUNTER, and QUANTISEQ algorithms showed that ZPR1 expression was closely related to macrophage infiltration (Figure 6C). Additionally, we screened 15 ICGs significantly associated with ZPR1, including effective immunotherapy targets CD47, CD274 (PD-L1), and CD276 (B7-H3) (Figure 6D). Then, we evaluated the immune-related functions in TCGA-HCC samples using ssGSEA method. As shown in Figure 6E, the immune infiltration scores of CD8 T cell, B cell, neutrophils, T helper cell, NK cell, and tumor infiltrating lymphocytes (TILs) were markedly decreased in the group with high ZPR1 expression. The high ZPR1 expression group had enhanced activity in MHC class I and attenuated activity in several immune functions, including type II IFN response, T cell co-stimulation, T cell co-suppression, inflammation promotion, and cytolytic activity (Figure 6F).

Figure 6 Immune landscape analysis of the ZPR1 gene. (A) TME score comparison in low- and high-ZPR1 groups. (B) The immune cell infiltration correlation analysis of the ZPR1 gene using the CIRBSORT algorithm. (C) The immune cell infiltration correlation analysis of the ZPR1 gene by multiple methods. (D) Correlation analysis between ZPR1 gene and immune checkpoint genes. (E and F) The correlation between the ZPR1 gene and immune cell infiltration scores or immune-related functions. *p < 0.05; **p < 0.01; ***p < 0.001.

Single-Cell Analysis of ZPR1 in HCC

In the next step, we applied TISCH database to looked into the relationship between ZPR1 and TME-associated cells in HCC. In GSE125449 dataset, ZPR1 was mainly expressed in malignant cells, stromal cells (endothelial, fibroblasts), and immune cells (B, CD8Tex, Mono/Macro, and Plasma). In contrast, ZPR1 expression was lower in the hepatic progenitor (Figure 7A and B). Notably, the ZPR1 expression was significantly up-regulated in the PD-L1/CTLA4 treatment group (Figure 7C). In GSE140228 dataset, ZPR1 was expressed in various immune cells (CD4 Tconv, CD8 T, Mono/Macro, and NK) (Figure 7D–F). In GSE146409 dataset, ZPR1 was highly expressed in endothelial, fibroblasts, Mono/Macro, and malignant cells (Figure 7G–I).

Figure 7 ZPR1 expression in TME-associated cells. Annotation of all cell types, ZPR1 expression percentages and differential expression in GSE125449 (A-C), GSE140228 (D-F), and GSE145409 (G-I) datasets. *p < 0.05; **p < 0.01; and ***p < 0.001. Abbreviation: Ns, no significance.

Validation of ZPR1 Up-Regulation in HCC with IHC

Subsequently, the IHC assay of HCC TMAs was applied to perform experimental validation. Compared to the weak scattered expression in the normal liver hepatocytes, ZPR1 staining showed a strong diffusely positive distribution in the HCC tissues (Figure 8A). In addition, we observed that ZPR1 was aggregated in the nucleus of some HCC cells. Using quantitative IHC analysis, ZPR1 expression was remarkably elevated in HCC tissues compared to adjacent normal hepatocytes (p < 0.001, Figure 8B). Surprisingly, we found that ZPR1 expression was negatively associated with CD4 (R= −0.28) and CD56 (R= −0.25), and positively correlated with CD68 (R= 0.26) and CD274 (R=0.24) (Figure 8C). Furthermore, the ZPR1 expression was considerably increased in PD-L1 TPS-positive (TPS>1%, p=0.0047) and CPS-positive (CPS>1%, p=0.032) groups (Figure 8D). The AUC value of ZPR1 for tumor diagnosis was 0.785, indicating that ZPR1 is robust in pathologically identifying tumors (Figure 8E). As demonstrated in Figure 8F, high expression of ZPR1 positively correlated with advanced T stage, pathological stage, and tumor grade in HCC.

Figure 8 ZPR1 immunostaining in HCC tissue samples. ZPR1 immunostaining showing low expression in the cytoplasm of normal hepatocytes, and prominent positive expression in the cytoplasm of HCC cells. (A) The ZPR1 expression between HCC cells and adjacent normal hepatocytes in the same visual field (200 ×). (B) The H-score difference of ZPR1 expression in HCC tissues. (C and D) Correlation between ZPR1 expression and common immune cell markers, as well as PD-L1. (E) The AUC value for tumor diagnosis of ZPR1. (F) Correlation analysis of ZPR1 expression with clinicopathological characteristics.

Treatment Response Prediction of ZPR1 in HCC

Finally, we applied TIDE algorithm to evaluate the ZPR1 response to ICB immunotherapy. As shown in Figure 9A, no difference in TIDE scores was observed between ZPR1 low and high expression groups. Furtherly, we found that the T cell exclusion score was markedly up-regulated in the ZPR1 high expression group, while the T cell dysfunction score was the opposite. Among the predictors, the immunosuppressive cells MDSC (myeloid-derived suppressor cells) and TAM-M2 (M2 subtype of tumor-associated macrophages) prediction scores were significantly elevated in ZPR1 high expression group. The predicted scores of IFNG and CD8 were negatively related to ZPR1. Furthermore, we investigated the impact of ZPR1 expression on response to chemotherapy and targeted therapy (Figure 9B). The comparative drug susceptibility analysis suggested that the IC50 values of camptothecin (DNA topoisomerase I inhibitor) and elesclomol (cuproptosis inducer) were significantly increased in the ZPR1 high expression group, suggesting that ZPR1 high expression is relatively resistant to elesclomol and camptothecin. In contrast, the IC50 values of paclitaxel (microtubule-associated inhibitor), doxorubicin (DNA topoisomerase inhibitor), rapamycin (mTOR inhibitor), gemcitabine (DNA synthesis inhibitor), sunitinib (multitargeted receptor tyrosine kinase inhibitor), and sorafenib (Raf inhibitor, VEGFR inhibitor, and ferroptosis activator) were decreased in the ZPR1 high expression group.

Figure 9 Treatment response prediction of ZPR1. (A) TIDE score prediction (B) The IC50 difference of common chemotherapeutic or targeted drugs in low- and high-ZPR1 expression groups. *p < 0.05; **p < 0.01; ***p < 0.001.

Discussion

ZPR1 has been previously proposed to play an essential function in SMA, hyperlipidemia, type 2 diabetes, and NAFLD. However, there is no systematic bioinformatics method to evaluate the role of ZPR1 in cancer through pan-cancer profiling and quantitative IHC analysis. This study aims to systematically evaluate the ZPR1 expression in different cancer tissue types and reveal the ZPR1 biological function in malignant progression and its potential as a novel independent prognostic predictor for HCC.

Firstly, we demonstrated that ZPR1 was upregulated in most cancer tissues. ZPR1 mRNA and protein expression were highly expressed in CRC, KIRC, LC, HNSC, GBM, PAAD, and HCC tissues. In addition, the hypo-methylation status of ZPR1 may be a potential cause of high expression of ZPR1 in tumors such as BRCA, HCC, and LC. Liu et al demonstrated that ZPR1 was over-expressed in BRCA and promoted breast cancer invasion and migration through ERK signaling.²⁵ Shan et al found that ZPR1 was lowly expressed in the LC cytoplasm and affected LC proliferation and invasion through the FAK-AKT axis.²⁶ Previous tumor studies of ZPR1 have been involved in LC and breast cancer, and our pan-cancer research provides directions for future research on ZPR1 in other malignancies.

In the HPA database and IHC experiment, we observed that ZPR1 was localized to vesicles and present in both the nucleus and cytoplasm of cancer cells. Liu et al found that ZPR1 was mainly strongly positively expressed in the cytoplasm and observed scattered in the nucleus.²⁵ Previous studies also reported that ZPR1 diffusely localized in cytoplasm and nucleus, and cycle-dependently co-localizes with SMN1 in the gems and Cajal bodies of the nucleus.¹⁴ After mitogen treatment induced quiescent cells in the proliferating state, ZPR1 was translocated from the cytoplasm to the nucleus.³⁸ Interestingly, emerging studies suggest that ZPR1 may be a putative transcription factor. Previously reported that ZPR1 deficiency caused defects in transcription and S-phase progression.¹⁴ In plant cells, ZPR1 was identified as a novel clock-associated nuclear DNA-binding regulator involved in the transcriptional regulation of the circadian rhythm regulatory gene BBX24.³⁶ The latest research demonstrated that over-expression of ZPR1 could promote SMN2 gene transcription in SMA patient cells.¹⁰ Here, we proposed that CREB1 may be a potential target gene of ZPR1. CREB1 has been demonstrated to be linked to hepatocarcinogenesis.^39,40 Chittilla et al have demonstrated the down-regulation of ZPR1 expression in high-fat diet, and hypothesized that ZPR1 might lead to cognitive decline by reducing PPAR-γ-CREB-BDNF expressions through transcriptional binding.⁴¹ These findings provide theoretical clues for our future cell biology experiments of ZPR1.

Furthermore, our research identified that BUD13 was significantly co-expressed with ZPR1. ZNF259-BUD13 is a gene cluster localized at 11q23.3. Increasing epidemiological studies demonstrated that the ZNF259-BUD13 gene cluster polymorphisms were strongly associated with plasma triglyceride levels in many regions, including China, Korea, India, Iran, European, Mexico, and the United States.^19,20,42,43 The ZNF259-BUD13 gene cluster is also related to liver metabolic processes.⁴⁴ The clinical pan-cancer correlation analysis indicated that ZPR1 was strongly correlated with gender, pathological stage, TMB, MSI, and stemness index of HCC. TMB and MSI were predictive indicators to assess tumor response to immunotherapy, therefore we explored the relevance of ZPR1 to immunotherapy in HCC. The stemness index was strongly correlated with tumor proliferation, metastasis, and drug resistance. Kannan et al found that ZPR1 was closely associated with R-loop accumulation and DNA damage.¹⁰ Our GSEA analysis also suggests that ZPR1 is closely associated with neurodegenerative diseases, DNA damage repair, lipid metabolism, and tumorigenesis-related pathways. DNA damage repair was reported as correlated with neuronal death⁴⁵ and tumor response to radiotherapy, chemotherapy, and immunotherapy.⁴⁶ Lipid metabolism is closely related to tumor development, even TME and immunotherapy.⁴⁷

TME is crucial for carcinogenesis, immune evasion, and immunotherapy resistance. The TME results demonstrated that ZPR1 was negatively correlated with the immune and estimated scores in HCC. Furthermore, a significant positive correlation of ZPR1 with immune infiltration of macrophages was observed. Single-cell analysis revealed that ZPR1 was closely associated with the immune infiltration of macrophages, T cells, and NK cells. The ZPR1 high expression group had lower CD8+ T cell and NK cell infiltration, as well as attenuated cytolytic activity and type II IFN response, showing the microenvironment characteristics of the cold tumor in HCC. In addition, we found that ZPR1 was positively correlated with CD47 (macrophage-associated ICG) and PD-L1 (T cell-associated ICG), and MHC1 immune function enhancement. The well-known “don’t eat me” signal MHC I-LILRB1⁴⁸ and CD47-SIRPα⁴⁹ axes play essential roles in inhibiting macrophage phagocytosis from clearing tumor cells and exerting immune escape function. The PD1-PD-L1 axis could also act as a phagocytic checkpoint involved in the inhibition of phagocytosis and tumor immunity.⁵⁰ These results suggest that the tumor cells with high ZPR1 expression attenuate the infiltration of immune cells, and potentially exert immune escape function by inhibiting macrophage phagocytosis through CD47, PD-L1 signaling. By Qupath quantitative analysis, we found that ZPR1 was highly expressed in HCC samples and closely correlated with the expression of CD4 (helper T cell phenotypic marker), CD56 (NK cell phenotypic marker), CD68 (macrophage phenotypic marker), and PD-L1. ZPR1 was closely correlated with TPS and CPS positive of PD-L1. The TIDE prediction results suggest that ZPR1 is associated with tumor immune escape and can be considered a regulator of ICB resistance. The results of ZPR1 response to immunotherapy, chemotherapy and targeted therapy can provide a novel theoretical foundation for the precise treatment and combination therapy in HCC.

Although pan-cancer and specific tumor bioinformatics analysis provides critical insights into our systematic understanding of the ZPR1 role in malignancy and TME, our research also has some limitations. For example, treatment strategy and response information for HCC patients are not available. The specific molecular mechanism of ZPR1 affecting immune infiltration needs to be elucidated by forthcoming in vitro (including flow cytometry analysis and cell biology experiments) and in vivo experiments. Nevertheless, our study showed that ZPR1 was up-regulated in HCC and closely associated with immune cell infiltration. ZPR1 may be a promising biomarker to predict patient prognosis and efficacy of anticancer therapy. These findings provide a theoretical basis for precise and individualized anticancer treatment strategies.

Conclusion

The pan-cancer profiling suggests that ZPR1 was differentially expressed in various cancers and is closely related to methylation status, clinicopathologic characteristics, TMB, MSI, and stemness index. Digital pathology analysis indicates that ZPR1 was remarkably upregulated in HCC and is closely associated with immune infiltration, which may be a potential biomarker for treatment response and therapeutic target for HCC.