Comprehensive Analysis of NKX3.2 in Liver Hepatocellular Carcinoma by Bigdata

Background and Objectives: The gene NKX3.2 plays a role in determining cell fate during development, and mutations of NKX3.2 have been studied in relation to human skeletal diseases. However, due to the lack of studies on the link between NKX3.2 and cancer, we aimed to provide insights into NKX3.2 as a new prognostic biomarker for liver hepatocellular carcinoma (LIHC). Materials and Methods: The clinical significance of LIHC was investigated using open gene expression databases. We comprehensively analyzed NKX3.2 expression in LIHC using Gene Expression Profiling Interactive Analysis 2, Tumor Immune Estimation Resource (TIMER), and Kaplan–Meier plotter databases. Then, we investigated the association between NKX3.2 expression and tumor-infiltrating immune cells (TIICs). Results: NKX3.2 expression was higher in the primary tumor group compared to the normal group, and expression was higher in fibrolamellar carcinoma (FLC) compared to other subtypes. When the prognostic value of NKX3.2 was evaluated, highly expressed NKX3.2 significantly improved the overall survival and had an unfavorable prognosis. In addition, NKX3.2 expression was associated with immune cell infiltration. Patients with low gene expression and high macrophage expression had a poorer survival rate than those with low NKX3.2 and low macrophage expression (p = 0.0309). Conclusions: High NKX3.2 expression may induce poorer prognosis in LIHC. In addition, these findings can be used as basic data due to the lack of available related research. However, further in vivo studies are essential to gain a deeper understanding of the biological role of NKX3.2 in LIHC and its potential implications for cancer development and progression.


Introduction
Liver hepatocellular carcinoma (LIHC) is the most common cancer worldwide with a poor prognosis [1,2].The most common type of LIHC is hepatocellular carcinoma (HCC), which begins in the main type of liver cell, the hepatocyte.Its other types, such as intrahepatic cholangiocarcinoma and hepatoblastoma, are much less common.There are limited data on the treatment of advanced liver cancer, even though LIHC is the common cancer among solid tumors [3].Moreover, advanced LIHC has high recurrence and metastasis rates, and a low survival rate [4].To overcome this challenge, finding potential prognostic biomarkers through big data research is urgently needed.
The first finding of NKX3.2 (Bapx1) was in Drosophila as a transcription factor [5,6].It is a member of the NK family of developmental genes and has a central role in determining cell fate during the development of various organs and visceral mesoderm during midgut myogenesis [5].NKX3.2 is involved in axial and limb skeleton formation and is known to protect chondrocytes from programmed apoptosis [7]; it is highly expressed in prochondrocytes and proliferating chondrocytes, but downregulated in mature chondrocytes, and plays Medicina 2023, 59, 1782 2 of 11 a role in the hypertrophic maturation of articular chondrocytes in osteoarthritis [5,7,8].NKX3.2 mutation is known as one of the genetic causes of spondylo-megaepiphysealmetaphyseal dysplasia (SMMD), a human skeletal disease [9].Interestingly, despite its roles in developmental and skeletal processes, research on NKX3.2 in the context of cancer has not been well-studied.Given its diverse functions and implications in developmental and skeletal processes, further exploration of NKX3.2 s potential involvement in cancer is warranted.
Herein, we demonstrated that NKX3.2 could be a target for the prediction of several cancers, for the first time.As a preliminary study on NKX3.2 in cancers, publicly available gene expression databases such as Total Cancer Genome Atlas (TCGA) and Tumor Immune Estimation Resource (TIMER) were used as recent studies [10,11]; NKX3.2 was highly expressed in most cancers compared to normal tissues.NKX3.2 expression was related to both cancer progression and rate of metastasis especially, in LIHC.Furthermore, we demonstrated that NKX3.2 expression was associated with the survival rate and immune cell infiltration in liver cancer patients.We conclude that NKX3.2 may be a novel biomarker for predicting prognosis and provides insight into the tumor immunogenicity of LIHC.

TIMER Database Analysis
TIMER is a TCGA database that includes immune infiltration and clinical data (age, sex, race, and stage) across 32 cancer types, for the NKX3.2 gene.We analyzed NKX3.2 gene expression and survival value, as well as the association between immune cell infiltration and cancer prognosis in LIHC, using correlation analysis among immune cells.In addition, the correlation between NKX3.2 gene expression and TIICs, such as with CD4+ and CD8+ T cells, B cells, macrophages, neutrophils, and dendritic cells, was analyzed through the gene module, and gene expression levels were plotted against tumor purity.In addition, we analyzed the association between immune cell infiltration in LIHC.

UALCAN Database Analysis
The UALCAN database (http://ualcan.path.uab.edu)uses TCGA level 3 RNA-seq and clinical data from LIHC (accessed on 1 February 2023).And then, the relative expression of genes across tumor and non-tumor tissues was compared according to cancer stages, tumor grade, or other clinicopathological features.

Immunohistochemistry (IHC) Staining Analysis
HPA (https://www.proteinatlas.org/) is a human proteome atlas database that contains information on protein distribution in human tissues and cells (accessed on 1 February 2023).To analyze the NKX3.2 proteomic expression levels, we obtained IHC images of NKX3.2 proteomic expression from the HPA for LIHC.The proteomic expression level of NKX3.2 was graded as not detected, low, medium, or high, based on the intensity of staining and the fraction of stained cells.

KM Database Analysis
KM is based on an online database, and includes survival rates, such as overall survival (OS), relapse-free survival (RFS), progression-free survival (PFS), disease-specific survival (DSS), and clinical data in patients with LIHC.Survival data of NKX3.2 were investigated using a log-rank test.

OSlihc Database Analysis
OSlihc (http://bioinfo.henu.edu.cn/DatabaseList.jsp)provides a platform for researchers to discover new prognostic biomarkers and may provide the opportunity to create novel targeted therapies for LIHC (accessed on 22 February 2023).To clarify the predictive value of genes in OSlihc, survival data such as OS, disease-free interval (DFI), progression-free interval (PFI), and DSS were generated, and OS was measured in all cohorts and combined cohorts, while DFI, PFI, and DSS were analyzed using TCGA.

LinkedOmics Database Analysis
The LinkedOmics database (http://www.linkedomics.org/admin.php) is a web-based platform for analyzing 32 TCGA cancer-associated multidimensional datasets (accessed on 26 February 2023).Correlation in NKX3.2 associated genes was analyzed by Pearson's correlation coefficient, and presented in volcano plots, heat maps, or scatter plots.

GEPIA2 Database Analysis
GEPIA2 (http://gepia2.cancer-pku.cn/#index)was used to integrate the clinical data acquired from 9736 tumors and 8587 normal tissues from TCGA and GTEx projects (accessed on 9 March 2023).GEPIA2 performs survival analyses according to gene or isoform expression levels.Based on gene expression in 33 different types of cancer, OS and RFS were generated by GEPIA2.Heatmaps of OS and DFS based on NKX3.2 gene expression across all TCGA tumor types were obtained in the "Survival Map" module.GEPIA2 showed KM curves based on NKX3.2 gene expression with the log-rank test and the Mantel-Cox test.

Statistical Analysis
Gene expression data were analyzed using online tools and the TIMER database.Survival analyses were performed by KM plots, TIMER, and OSlihc.Survival results are presented as HR and p values by a log-rank test.The correlation between gene expression and immune signature score was performed by Spearman's correlation.All data in this study were extracted from public open databases, and we conducted statistical analysis by web tools.All results are presented as p-values from the log-rank test.Statistical significance was set at p < 0.05.
We specifically investigated the clinical characteristics of NKX3.2 expression in LIHC using the UALCAN database to validate the results from TIMER.The NKX3.2 gene showed significant results in Asian and Caucasian patients.In addition, it was overexpressed in the primary tumor group compared to the normal group (Figure 1B), and its expression level was significantly different according race (Figure 1C).NKX3.2 was more significantly expressed with increasing tumor histological and clinical stages (Figure 1D,E).Furthermore, lymph node metastasis was associated with NKX3.2 expression (Figure 1F).When compared to other histological subtypes, NKX3.2 showed higher expression in FLC (Figure 1G).

Protein Expression of NKX3.2 in LICH
To assess the expression of NKX3.2 at the protein level, the protein expression levels of NKX3.2 were obtained from the HPA database.IHC images from the HPA revealed that the protein expression of NKX3.2 was undetected in normal liver tissue, whereas protein expression of NKX3.2 was significantly higher in LICH tissues than in normal tissues (Figure 2).Our results revealed overexpression of NKX3.2 at transcriptional and translational levels in patients with LIHC.

Protein Expression of NKX3.2 in LICH
To assess the expression of NKX3.2 at the protein level, the protein expression levels of NKX3.2 were obtained from the HPA database.IHC images from the HPA revealed that the protein expression of NKX3.2 was undetected in normal liver tissue, whereas protein expression of NKX3.2 was significantly higher in LICH tissues than in normal tissues (Figure 2).Our results revealed overexpression of NKX3.2 at transcriptional and translational levels in patients with LIHC.

Multifaceted Prognostic Value of NKX3.2 in LIHC
The current study used three databases to evaluate NKX3.2 in LIHC.To date, the prognostic value of NKX3.2 has not been reported, and the correlation between the expression level of NKX3.2 and prognosis in LIHC were evaluated using the GEPIA2 database.These results demonstrated the effect of the expression level of NKX3.2 on cancer survival and recurrence rates.

Multifaceted Prognostic Value of NKX3.2 in LIHC
The current study used three databases to evaluate NKX3.2 in LIHC.To date, the prognostic value of NKX3.2 has not been reported, and the correlation between the expression level of NKX3.2 and prognosis in LIHC were evaluated using the GEPIA2 database.These results demonstrated the effect of the expression level of NKX3.2 on cancer survival and recurrence rates.

Association of NKX3.2 Expression with Immune Cell Infiltration in LIHC
Immune cell infiltration usually accelerates cancer progression and affects survival outcomes.Therefore, we examined the association between NKX3.2 and six tumor immuneinfiltrating cells using the TIMER database.Neutrophils, macrophages, and CD8+ T, CD4+ T, B, and myeloid dendritic cells all showed a positive correlation.The results were NKX3.2 and CD8+ T cells (rho = −0.022,p < 0.001), CD4+ T cells (rho = 0.34, p = < 0.001), and B cells (rho = 0.279, p < 0.001), macrophages (rho = 0.267, p < 0.001), neutrophils (rho = 0.193, p < 0.001) and dendritic cells (rho = 0.399, p < 0.001) (Figure 4A).We then investigated the correlation between the tumor immune infiltrating cells, NKX3.2 expression, and survival time in LIHC.The results showed that patients with high NKX3.2 and high neutrophil expressions had a lower survival rate than patients with high gene expression and low neutrophil expression, though it did not achieve statistical significance (p = 0.0763).The patient group with low gene expression and high macrophage expression had a poor survival rate compared to the patient group with low NKX3.2 and low macrophage expression (p = 0.0309) Also, we investigated whether high NKX3.2expression was related with prognosis and TIICs in LIHC.Patients with high NKX3.2expression and high neutrophil infiltration levels showed worse prognosis than those with low NKX3.2 expression and low neutrophil infiltration levels.Low NKX3.2 expression and high neutrophil infiltration levels were associated with a worse prognosis than low NKX3.2 expression and low neutrophil infiltration levels.High NKX3.2 expression and low CD4+T cell infiltration levels were associated with a worse prognosis than low NKX3.2 expression and high CD4+T cell infiltration levels.High NKX3.2 expression and high macrophage infiltration levels predicted a worse result than low NKX3.2 expression and low macrophage infiltration levels.Moreover, low NKX3.2 expression and high macrophage infiltration levels were associated with a worse prognosis than low NKX3.2 expression and low macrophage infiltration levels.B cells, CD8+T cells, and dendritic cells showed no significant level of prognosis (Figure 4B).
Therefore, the analysis of the cumulative curve indicated that the infiltration of immune cells was associated with the NKX3.2 gene in LIHC, and affected the prognosis.These data suggest that NKX3.2 can be the prediction gene for the immune invasion process in LIHC.

Co-Expression Analysis of the NKX3.2 Gene
To obtain biological insights into NKX3.2 in LIHC, we analyzed the co-expressed genes of NKX3.2 using the LinkedOmics database, and aimed to identify differentially expressed genes that correlated with the NKX3.2 gene in LIHC.By Pearson correlation analysis, we determined positively and negatively correlated differentially expressed genes for NKX3.2.Of the 19,914 genes analyzed, 13,131 were positively correlated with the NKX3.2 gene in LIHC, and 6782 were negatively correlated.Volcano plots showing these correlations are depicted in Figure 5A (Left).Additionally, we picked 50 genes to make lists of positive and negative correlations with the NKX3.2 gene based on our analysis (Figure 5A (Right)).The top five positively correlated genes were NKX3.2,LOC285548, RAVER1, PKM2, and C11orf84, and the top five negatively correlated genes were aQP9, tAKR, SLC27A2, HAO1, and CES2.Furthermore, GO (biological process, cellular component, and molecular function) and KEGG pathway analyses of NKX3.2-related genes were performed.GO annotations presented that NKX3.NKX3.2 and high neutrophil expressions had a lower survival rate than patients with high gene expression and low neutrophil expression, though it did not achieve statistical significance (p = 0.0763).The patient group with low gene expression and high macrophage expression had a poor survival rate compared to the patient group with low NKX3.2 and low macrophage expression (p = 0.0309) Also, we investigated whether high NKX3.2expression was related with prognosis and TIICs in LIHC.Patients with high NKX3.2expression and high neutrophil infiltration levels showed worse prognosis than those with low NKX3.2 expression and low neutrophil infiltration levels.Low NKX3.2 expression and high neutrophil infiltration levels were associated with a worse prognosis than low NKX3.2 expression and low neutrophil infiltration levels.High NKX3.2 expression and low CD4+T cell infiltration levels were associated with a worse prognosis than low NKX3.2 expression and high CD4+T cell infiltration levels.High NKX3.2 expression and high macrophage infiltration levels predicted a worse result than low NKX3.2 expression and low macrophage infiltration levels.Moreover, low NKX3.2 expression and high macrophage infiltration levels were associated with a worse prognosis than low NKX3.2 expression and low macrophage infiltration levels.B cells, CD8+T cells, and dendritic cells showed no significant level of prognosis (Figure 4B).
Therefore, the analysis of the cumulative curve indicated that the infiltration of immune cells was associated with the NKX3.2 gene in LIHC, and affected the prognosis.These data suggest that NKX3.2 can be the prediction gene for the immune invasion process in LIHC.

Prognosis Analysis of the NKX3.2-Related Genes
We investigated the prognostic value of NKX3.2-related genes in LIHC using data from the ZEPIA2 database.NKX3.2-related genes are most likely to be high-risk genes for LIHC.Among the positively related genes of NKX3.2, 35 showed a high hazard ratio (HR) for OS (Figure 6A), and 24 showed a high HR for RFS (Figure 6B).In contrast, 24 genes of the negative genes had low HR (p < 0.05) in OS (Figure 6C), and 16 genes had low HR in RFS (Figure 6D).Also, we showed that positively related genes of NKX3.2 indicated a high

Discussion
LIHC is a malignant tumor, and its incidence is continuously increasing.After primary treatment, such as surgical resection and transplantation, LIHC has a high recurrence rate and 5-year survival rate of 50-70% [12].To date, NKX3.2 has been studied extensively; however, our findings are the first to show the relation with cancer.This study will help to predict the prognosis of patients with LIHC and improve future outcomes.
NKX3.2, a well-known gene for its roles in developmental processes and human skeletal diseases, has not been previously studied in the context of cancer.Through comprehensive analyses of gene expression databases, we have identified a significant correlation between NKX3.2 expression and clinical characteristics, especially in LIHC.
The pathological types in which LIHC develops occur in chronic hepatitis and fibrosis [13].In particular, viral hepatitis infection is a potential risk factor for LIHC.Asian and Caucasian patients showed a higher expression of NKX3.2 in LIHC; however, only Asian patients had significant OS results.According to previous studies, 60% of Asian patients with LIHC are HBV-positive, whereas only 25% of European patients are carriers [14,15].Similarly, OS results from men showed poor outcomes compared to results from women.NKX3.2 will be a significant target to predict other types of liver diseases, such as chronic hepatitis, cirrhosis, and non-alcoholic fatty liver disease.Therefore, NKX3.2 may have an important role in the carcinogenesis of Asian patients with LIHC.
Primary liver cancer is a frequent cancer worldwide, and hepatocellular carcinoma (HCC) is the main type of primary liver cancer, accounting for 75-85% of diagnoses [16].Our results show that NKX3.2 expression is high in primary tumors.Additionally, FLC, as a rare and unknown cancer type, is more prevalent in younger patients than other subtypes [17,18].Recent studies have shown that FLC has a better overall prognosis than other primary liver tumors such as HCC and intrahepatic cholangiocarcinoma.However, due to the lack of data on FLC, it is said that it is difficult to develop new treatments.Also, because of the rarity of FLC, it has an aggressive tumor with a low overall cure rate, and recurrence is common even after surgery [18].In this study, NKX3.2 was overexpressed in FLC.Our study is limited because it is unknown whether a high expression of NKX3.2 causes primary tumors and FLC.Additional experiments will be helpful in the development of therapeutic agents.
Regarding immunotherapy, the tumor microenvironment (TME) is known to affect the efficacy of immunotherapy [19,20].As far as is known, studies have shown that immunotherapy provides a survival advantage in LIHC treatment [21,22].In this study, the level of immune infiltration of NKX3.2 in LIHC was explored.LIHC is an inflammationrelated cancer and is known to be closely related to the immune response [23].In particular, the immune response of LIHC is regulated by various activating and inhibitory signaling pathways.Notably, NKX3.2 plays an important role in regulating LIHC immunology.Immunopenetration analysis showed that NKX3.2 was significantly associated with various immune cells in LIHC.Analyses using the TIMER database showed that NKX3.2 was positively correlated with the infiltration status of immune cells including neutrophils, macrophages, and CD4+ T, CD8+ T, B, and dendritic cells.Therefore, TIICs are closely related to tumorigenesis and its development.In our study, high NKX3.2expression and cell infiltration of neutrophils showed a worse prognosis than low NKX3.2 expression.
According to previous studies, inflammation is essential for cancer initiation by damaging tissues, and neutrophils play an essential role in this process [24].In addition, neutrophils have recently received a lot of attention because they are associated with cancer growth promotion and metastasis [24].Future studies are needed to determine the role of the NKX3.2 gene in neutrophil expression.
It is also important to note that although some associations between NKX3.2 expression, immune cell infiltration and survival outcomes were observed, not all reached statistical significance.Further studies and studies with larger sample sizes are needed to validate these results and understand the precise mechanism for the role of NKX3.2 in immune invasion and cancer progression in LIHC.

Figure 2 .
Figure 2. Protein expression levels of NKX3.2 in LIHC.Protein expression of NKX3.2 was analyzed using the Human Protein Atlas (HPA) database.

Figure 3 .
Figure 3. Prognostic significance of high expression of NKX3.2 in LIHC.(A) The prognostic value of NKX3.2 was analyzed using KM.(B) Evaluation of the prognostic value of NKX3.2 in Oslihc.p-value, confidence interval (95%CI) and number at risk are as shown.The y-axis represents survival rate and the x-axis represents survival time (months).

Figure 4 .
Figure 4. Correlation between NKX3.2 expression and infiltrating immune cells in LIHC.(A) The correlation between NKX3.2 and infiltrating immune cells was analyzed using the TIMER database.(B) Survival analysis according to NKX3.2 gene expression and expression of TIICs.
2 co-expressed genes were predominantly co-expressed during cell cycle G2/M phase transition, chromosome segregation, DNA conformation and change, mitotic cell cycle phase transition, cell cycle checkpoint, etc. (Figure 5B, Left).KEGG pathway assays demonstrated enrichment of the cell cycle, DNA replication, spliceosome, and microRNAs in cancer (Figure 5B, Right).

Figure 5 .
Figure 5. Correlation of NKX3.2 and related genes in LIHC.(A) The volcanic map showed correlations between NKX3.2 and genes differentially expressed in LIHC.Top 50 genes positively associated with NKX3.2 in LIHC and top 50 genes negatively associated with NKX3.2 in LIHC (B) GO and KEGG analyses performin LIHC.

Figure 5 .
Figure 5. Correlation of NKX3.2 and related genes in LIHC.(A) The volcanic map showed correlations between NKX3.2 and genes differentially expressed in LIHC.Top 50 genes positively associated with NKX3.2 in LIHC and top 50 genes negatively associated with NKX3.2 in LIHC (B) GO and KEGG analyses performin LIHC.

Figure 6 .
Figure 6.Prognostic value of NKX3.2-related genes in LIHC.Survival map of the positively related genes of NKX3.2 in OS (A) and RFS (B).Survival map of the negatively related genes of NKX3.2 in OS (C) and RFS (D).Survival map presenting the log10 (HR) of the genes in LIHC.A square with bold border represents a p value < 0.05 in the survival analysis.