Introduction

The most frequent primary liver cancer is HCC, which also happens to be the sixth most common tumor. This malignancy has been gradually on the rise in terms of incidence and mortality in recent years. Numerous people continue to receive diagnoses of advanced-stage cancer with a wide range of symptoms, despite the significant improvements in cancer screening techniques.1 However, there are not many HCC patients who can benefit from successful therapy options. In order to improve all therapy approaches, it is vital that we explore more diagnostic biomarkers and potential therapeutic targets.

The water-soluble vitamin B3 (niacin) has an amide form called NAM, which is more frequently found in foods like meat, fish, beans, mushrooms, nuts, and cereals. By use of the enzyme nicotinamide phosphate ribosyl transferase (NAMPT), NAM is directly transformed in living cells into nicotinamide mononucleotide (NMN), which subsequently binds to ATP to create nicotinamide adenine dinucleotide (NAD+). NAM, the fundamental component of cellular energy metabolism, is crucial for the formation of NAD+.2 Cancer is characterized by metabolic dysregulation, which promotes unregulated cancer growth.3 NAM is therefore thought to impact tumor growth through modifying cell metabolism.

Based on the TCGA and GeneCards datasets, this study intended to locate nicotinamide-related genes in HCC and then confirmed the expression of these genes in multiple datasets acquired from the GEO database. Using the CIBERSORT approach, the link between niacinamide expression level and immune cell infiltration was investigated. The drug sensitivity investigation of NAM major genes related with HCC used the Genomics of Cancer Cell Line Encyclopedia (CCLE), Drug Sensitivity in Cancer (GDSC), and CellMiner databases. Our research identifies a trustworthy and feasible treatment target for HCC.

Materials and MethodsData Acquisition

The Liver Hepatocellular Carcinoma (LIHC) dataset was downloaded from TCGA database (https://portal.gdc.cancer.gov/).4 374 cases of liver cancer samples (grouping: LIHC) and 50 cases of paracancerous samples (grouping: Normal) were got from the UCSC Xena database (http://genome.ucsc.edu).5 The limma package was performed to normalize the count sequencing data of the TCGA-LIHC data set.6 Using tool GEOquery, the LIHC-related datasets GSE25097, GSE46408 and GSE844027 were obtained.7–11 The probe names of the datasets were annotated using the pertinent GPL platform files (supplementary Table S1).

Identification of Nicotinamide Metabolism-Related Genes (NMRGs)

Forty-two NMRGs were acquired from related references, and eight NMRGs were obtained from the GeneCards database.12 Additionally, the single nucleotide polymorphism (SNP) data from the TCGA-LIHC dataset’s somatic mutations were retrieved and displayed using the R package maftools.13 The TCGA-LIHC dataset’s “Copy Number Variation” data were obtained by R package TCGAbiolinks for GISTIC 2.0 analysis.14 The differentially expressed genes (DEGs) were intersected with NMRGs to produce the NMRDEGs linked to LIHC disease.

Functional Enrichment Analysis

GO and KEGG enrichment analysis on NMRDEGs was carried out using the clusterProfiler package.1516 The ‘h.all.v7.4.symbols.gmt’ genome was used to determine whether there is significant enrichment which was got from the Molecular Signatures Database (MSigDB). The expressed genes were then subjected to Gene Set Enrichment Analysis(GSEA) and Gene Set Variation Analysis (GSVA).1718

LASSO Regression Analysis to Screen NMRDEGs

The LASSO regression model’s risk score, or riskScore, was subsequently calculated based on the NMRDEGs after the penalty coefficient (lambda) of the obtained NMRDEGs was extracted.19

Prognostic Clinical Analysis

A Cox regression model was created using single/multivariate Cox regression analysis on the expression of key genes (mRNA) set to examine the clinical prognostic value of the target gene on HCC. The results were shown by a nomogram diagram which was evaluated by Calibration curves and R software package ggDCA. To figure out related genes with significant differences, K-M curves were created for NMRDEGs.

Immune Infiltration Assay (CIBERSORT)

Each TCGA-LIHC sample’s NM scores were determined using the ssGSEA method and the GSVA package of the R package. The TCGA-LIHC data set’s immune cell infiltration state was assessed using the CIBERSORT algorithm.20 First, the LIHC group’s infiltration differences of 22 different types of immune cells were examined and displayed in group comparison charts. The relevance of immune cells was then depicted using a correlation heat map.

Drug Sensitivity Analysis of Key Genes

Based on the expression levels of key genes and the drug data in GDSC,21 CCLE,22 and CellMiner databases,23 drug sensitivity analysis was conducted on key genes.

Receiver Operating Characteristic (ROC) Curve Analysis

To assess the diagnostic impact of key genes’ expression on diseases, the area under the ROC curve (AUC) was figure The diagnostic impact increases as the AUC gets closer to 1.

Western Blot

Using lysis buffer (KeyGEN, China) to extract proteins from HCC tissues, relative quantities were measured using Bradford (Thermo, USA). The protein was supplied at roughly the same amount and concentration. The following commercial antibodies were bought: NT5C, β-actin, and NADSYNI were acquired from Affinity Biosciences (Beijing, China), Cell Signaling Technology (USA), and Abcam (USA), respectively.

EdU Assay

NAM solution (10, 20, and 50 mM) were used to treat HCC cells. Cell proliferation was measured by EdU assay.

Wound-Healing Assay

The HCC cells were grown until a confluent monolayer developed. A 10-μL pipette tip was used to scrape the monolayer after that. Then the cells were treated with NAM solution (10, 20, and 50 mM). Every 24 hours, the migratory distance was shown and captured on microscope.

Migration and Invasion Assay

The transwell top chamber (pore size 12μm, BD, USA) received HCC cells, and the transwell lower chamber received 600μL of DMEM with 15% FBS. The cells were fixed, stained, and counted following a 24-hour incubation period.

Statistical Analysis

Data analysis for this research was executed using R software (version 4.1.2). For continuous data between two groups, the Mann–Whitney U-test were used for non-normally distributed variables and the independent Student’s t-test were used for normally distributed ones. The chi-square test or Fisher’s exact test were applied for comparing categorical data. The Spearman’s rank correlation analysis was used to determine the correlations between molecules. All statistical analyses were conducted with a two-sided significance level, and P-values surpassing the predetermined threshold were considered to indicate statistical significance. GraphPad Prism 10.0 were used to plot the results of the experiments.

ResultsIdentification of NMRDEGs in LIHC

Figure S1 depicts the research flow chart for our investigation. The data for LIHC datasets GSE25097, GSE46408, and GSE84402 were downloaded, and the removal of batch effects was evaluated through distribution boxplots and PCA graphs (Figure S2A–D).

The TCGA database was used to download the information for 374 LIHC samples. We utilized the limma software to identify the genes that were differentially expressed on the TCGA-LIHC data set. In the result, there are 13,305 genes that meet the criteria. The result was shown by a volcano map (Figure 1A). Finally, thirty-three NMRDEGs were obtained after intersecting DEGs and NMRGs, and the result shown by a Venn diagram (Figure 1B).

Figure 1 Analysis of differential genes in the liver cancer dataset. (A) Differential gene volcano map of dataset TCGA-LIHC. (B) Venn diagram of DEGs and NMRGs in dataset TCGA-LIHC. DEGs: differentially expressed genes. NMRGs: Nicotinamide metabolism-related genes.

NMRDEGs Functional Enrichment Analysis (GO) and Pathway Enrichment (KEGG) Analysis

The 33 NMRDEGs were analyzed using the GO (supplementary Table S2) and KEGG (supplementary Table S3) enrichment analysis, and the results were shown in the diagrams below (Figure 2A–D). The results showed that 33 NMRDEGs were most heavily involved in negative regulation of NF-kappaB transcription factor activity (GO: 0032088), reaction to hypoxia (GO: 0001666), circadian rhythm (GO: 0007623), response to oxidative stress (GO: 0006979) in BP, nucleotide activity (GO: 0008252), glycosyltransferase activity (GO: 0016757), phosphatase activity (GO: 0016791) and nucleotide-triphosphate diphosphate activity (GO: 0047429) in MF. The pathways of nucleotide metabolism (hsa01232), pyrimidine metabolism (hsa00240), purine metabolism (hsa00230), and cofactor biosynthesis (hsa01240) are primarily affected by the results of the KEGG enrichment analysis.

Figure 2 NMRDEGs functional enrichment analysis (GO) and pathway enrichment (KEGG) analysis. The results of GO functional enrichment analysis and KEGG pathway enrichment analysis of NMRDEGs are displayed in bubble diagrams (A), circular network diagram (B), Histogram display (C), circle diagrams (D). NMRDEGs: Nicotinamide metabolism-related differentially expressed genes. GO: Gene Ontology. BP: biological process. MF: molecular function. KEGG: Kyoto Encyclopedia of Genes and Genomes. The screening criteria were P.value < 0.05 and FDR value (q.value) < 0.05.

GSEA and GSVA Enrichment Analysis of LIHC Dataset

We examined gene expression and involved biological processes among various groups (Normal/LIHC) in HCC patient samples in the data set TCGA-LIHC to ascertain the effect of gene expression levels on HCC through GSEA enrichment analysis (supplementary Table S4). The findings demonstrated that all TCGA-LIHC genes were considerably enriched in five key pathways (Figure 3A), including the TP53 pathway (Figure 3B), Notch pathway (Figure 3C), Wnt pathway (Figure 3D), Jak-Stat pathway (Figure 3E), Pi3k Akt pathway (Figure 3F), etc. We conducted GSVA analysis on the data set TCGA-LIHC (Figure 3G) to calculate the functional enrichment difference between liver cancer samples and corresponding normal samples. When LIHC samples were compared to matched normal samples in the dataset TCGA-LIHC, the results revealed variations in the WNT_BETA pathway and other gene sets (supplementary Table S5).

Figure 3 GSEA enrichment analysis of dataset TCGA-LIHC. (A) Five main biological characteristics of the GSEA enrichment analysis. The differential genes were significantly prominent in TP53 pathway (B), Notch pathway (C), Wnt pathway (D), Jak Stat pathway (E), Pi3k Akt pathway (F) etc. (G) . GSVA analysis in dataset TCGA-LIHC. GSEA: Gene Set Enrichment Analysis. The significant screening criteria were P.value < 0.05 and FDR value (q.value) < 0.25.

LASSO Regression Analysis to Screen NMRDEGs

In order to identify six prognostic key genes (NAXE, NADSYN1, NT5C, NT5C3A, PNP, and NT5E), we screened NMRDEGs using LASSO regression analysis (Figure S3A). The LASSO variogram (Figure S3B) revealed that the number of genes with coefficients of 0 grew gradually as decreased and that the genes changed with the coefficient (logarithmic posterior) of the LASSO penalty term. Using a risk factor map (Figure S3C) and the ROC curve of the LASSO regression model riskScore for various outcome groups (Figure S3D), respectively, we simultaneously showed the risk factor grouping of NMRDEGs. The ROC curve demonstrated that the riskScore component of the LASSO regression model performed poorly (0.5 < AUC < 0.7) among various outcome group. Finally, we created a riskScore correlation lollipop diagram (Figure S3E) between the level of NMRDEG expression and the LASSO regression model.

Expression Analysis of Key Genes in the Data Set TCGA-LIHC and GEO Datasets

According to the findings of the differential analysis, the expression levels of the six essential genes (NAXE, NADSYN1, NT5C, NT5C3A, PNP, and NT5E) in the two groups (Normal/LIHC) were all highly statistically significant (Figure 4A). Correlation analysis was conducted using the entire expression matrix of the six important genes, and a correlation heat map was created (Figure 4B). It was discovered that the genes NAXE, PNP, and NT5E have a significant negative link, but NAXE, NADSYN1, NT5C, and NT5C3A have a significant positive connection. According to the ROC curve, there is a considerable degree of accuracy (0.7 < AUC < 0.9) in the expression levels of NAXE, NADSYN1, NT5C and NT5C3A. The expression levels of PNP and NT5E revealed poor distinction between groups (0.5 < AUC < 0.7) (Figure 4C-H). Violin plot comparisons were used to show the results of the differential analysis of the expression levels of NAXE, NADSYN1, NT5C, NT5C3A, PNP, and NT5E in the merged dataset’s different groups (Normal/LIHC) (Figure S4A). The results revealed that, after removing the missing genes from the merged dataset, the expression of the key gene NT5E varied significantly between the Normal and LIHC groups. According to the correlation heat map based on the complete expression matrix of the key genes in the merged datasets (Figure S4B), NT5E and NT5C are positively correlated, but negatively correlated with NADSYN1. The expression levels of PNP and NT5E in the various groups (Normal/LIHC) in the datasets GSE25097 and GSE84402 were highly statistically significant (Figure S4C-D).

Figure 4 The expression of key genes in the data set TCGA-LIHC. (A) The group comparison chart results of key genes in Normal/LIHC groups. (B)The correlation heat map results of key genes. ROC curve analysis of key genes NAXE (C), NADSYN1 (D), NT5C (E), NT5C3A (F), PNP (G), NT5E (H). *P < 0.05, ***P < 0.001. LIHC: Liver hepatocellular carcinoma. ROC: Receiver operating characteristic curve.

Prognostic Clinical Manifestations and Performance of Key Genes

This demonstrated a relationship between the clinical prognosis and the degree of key gene expression by univariate and multivariate Cox regression analysis (supplementary Table S6, Figure 5A). The predictive performance of the Cox regression model was subsequently evaluated through a nomogram analysis (Figure 5B), revealing that, within the diagnostic model for LIHC, NAXE exhibited the highest efficacy, whereas NT5C demonstrated the lowest effectiveness. Additionally, our multivariate Cox regression model ran calibration analyses for the 1-year, 3-year, and 5-year prognoses and created calibration curves (Figure 5C-E). We discovered that the model’s predicted survival for the patients was largely in line with their actual survival. Then, we assessed the therapeutic value of the created Cox regression prognostic model by decision curve analysis (DCA) (Figure 5F-H). The outcome demonstrated that the 5-year prognosis model has the best clinical value.

Figure 5 Prognostic performance of key genes. (A and B). Univariate regression analysis forest plot (A) and nomogram of key genes (B). (C-E). Calibration curves of the Cox regression prognostic model for 1-year (C), 3-year (D), and 5-year (E) survival probabilities. The x-axis represents the predicted survival probability, while the y-axis indicates the actual survival probability. Various colored lines and points denote model predictions at different time points. Lines closer to the gray ideal-case line suggest improved prediction accuracy. (F-H). DCA plots of the Cox regression prognostic model for 1-year (F), 3-year (G), and 5-year (H) survival probabilities. The x-axis of the DCA diagram displays the threshold probability, and the y-axis shows the net benefit. The model’s effectiveness is evident when its line consistently surpasses the all-positive and all-negative lines over an extensive x–value range.

Prognostic Performance of Key Genes

In the TCGA-LIHC datasets, Kaplan-Meier (K-M) survival and prognostic curves were generated for key genes. Our findings reveal that NAXE, NADSYN1, NT5C, and NT5C3A demonstrate statistical significance within the TCGA-LIHC datasets (Figure S5A-D).

Mutation Analysis of Key Genes in LIHC Patients

We counted the six key genes in the LIHC patient samples in order to investigate the somatic mutations of the key genes (NAXE, NADSYN1, NT5C, NT5C3A, PNP, and NT5E) through the “maftools” R package. The findings indicated that the TCGA-LIHC data set contained two primary categories of somatic mutations: missense mutations and frameshift Deletion Mutations (Frame Shift Del), with missense mutations making up the majority of them (Figure 6A). In addition, single nucleotide polymorphisms (SNPs) and a limited number of deletions (DELs) predominated among the mutation types of the six key genes in LIHC patients. Additionally, C > A is the most common single nucleotide variant (SNV). (Figure 6B). We also looked at the copy number variation (CNV) of six key genes for patients with LIHC. The findings demonstrated that there are four key genes that have numerous amplifications and deletions in LIHC patient samples, with the top three genes for amplification being NT5C, NT5E, and NADSYN, and the top three for deletion being NT5E, PNP, and NADSYN (Figure 6C).

Figure 6 Mutation analysis of key genes in LIHC patients. (A) The somatic mutation status of key genes. (B) Proportion results of key genes. (C) The copy number variation of key genes.

Drug Sensitivity Analysis of Key Genes

We also used the GDSC, CCLE and the CellMiner database to construct a ridge regression model using the pRRophytic algorithm. We predicted the sensitivity of key genes to common anticancer drugs using IC50 values. The results showed that key genes can find multiple drugs with interaction relationships in GDSC, CCLE, and CellMiner databases (Figure 7A-C).

Figure 7 Drug sensitivity analysis. The drug susceptibility analysis results of key genes in GDSC database (A), CCLE database (B) and CellMiner database (C) are displayed.

Immune Infiltration Analysis of High and Low Nicotinamide Metabolic Score Groups (CIBERSORT)

We used the ssGSEA algorithm to determine the Nicotinamide Metabolism (NMs) in the TCGA-LIHC data set based on the expression of six key genes (NAXE, NADSYN1, NT5C, NT5C3A, PNP, NT5E). The results revealed that the expression levels of the six key genes in the TCGA-LIHC patient samples were highly statistically significant between the high and low levels of NMs score (Figure S6A). Additionally, we added six key genes that are co-expressed in heatmaps (Figure S6B). In our subsequent steps, we plotted the ROC curve for the six pivotal genes, categorizing them into High/Low groups based on the NMs score from TCGA-LIHC patient samples. Notably, the ROC curve revealed that the key genes, NT5E and NT5C, effectively predicted samples within the High/Low NMs score groups. Conversely, within the same dataset, the genes NT5C3A, PNP, NADSYN1, and NAXE were less adept at differentiating the High/Low NMs score groups (Figure S6C-H).

We used the CIBERSORT method to determine the association between the sample expression profile data of 22 distinct types of immune cells in the TCGA-LIHC data set in order to further investigate the relationship between the nicotinamide metabolism score and antitumor immunity in LIHC patients (Figure 8A). The findings demonstrated that nine immune cells had statistically significant variations in the NMs score (High/Low) groupings expression levels (Figure 8B). Additionally, we proved the relationship between the degree of immune cell infiltration with statistically significant changes. The results suggest a mild positive association between the abundance of resting dendritic cells in the infiltrate and M1 macrophages. Conversely, there’s a mild negative association between the infiltrating dendritic cells and M2 macrophages, monocytes, and resting NK cells (Figure 8C). We found the interplay between the abundance of immune cell infiltration and six pivotal genes (NAXE, NADSYN1, NT5C, NT5C3A, PNP, NT5E), noting statistically significant differences. Specifically, a positive correlation emerged between macrophage M0 and the key gene NT5C (Figure 8D).

Figure 8 Immune infiltration analysis (CIBERSORT) between high and low nicotinamide metabolism scores (High/Low) groups. (A) The immune infiltration results of 22 kinds of immune cells between High/Low groups of the nicotinamide metabolism score. (B) Group comparison chart of immune cells in High/Low groups of the nicotinamide metabolism score. (C) The correlation analysis results between immune cells. (D) The correlation analysis results between key genes and immune cell infiltration abundance.

Verification of Selected Key Genes by Western Blot

Western blot was utilized to validate NADSYN1 and NT5C, which had not received much prior investigation in HCC, in tissue samples obtained from patients. Compared with paracancerous tissue, the protein levels of NADSYN1 and NT5C are high expressed in HCC tissue which were obtained from HCC patients (Figure 9A).

Figure 9 The relative expression of NAM-related genes in tissues and the effect of nicotinamide on cell migration and proliferation in HCC. (A) protein levels of NADSYN1 and NT5C in HCC tissues and paired paracarcinoma tissues. NL normal liver, HCC hepatocellular carcinoma. A CCK-8 assay (B and C) and an EdU assay (D) showed the cell proliferation ability of Huh7 and HepG2 cultured by 0, 30, or 50 mM of NAM for 48 h. Wound-healing assays (E) and Transwell assays (F) showed the cell proliferation ability of HCC cells cultured by 0, 10, 30, or 50 mM of NAM for 48 h. The data are expressed as the mean±SD; * P< 0.05, ** P < 0.01, ****P< 0.0001.

NAM Inhibits Cell Proliferation and Migration of HCC Cells

The CCK-8 assay (Figure 9B-C) and the EdU assay (Figure 9D) revealed that the experimental groups proliferated less than the control group did. Furthermore, Huh7 and HepG2 cell migration was greatly reduced by overexpressing NAM, as demonstrated by wound-healing assays (Figure 9E) and transwell assays without Matrigel (Figure 9F).

Discussion

Uncontrolled cellular energy metabolism has been linked in several studies to the onset and growth of tumors.24 NAM is well known for being a crucial regulator of REDOX processes and mitochondrial metabolism, which controls cellular energy metabolism. NAM’s contribution to cancer treatment and prevention is being supported by more and more research. Numerous researches have recently concentrated on the function of nicotinamide in the detection and treatment of tumors.25 Studies on nicotinamide’s function in HCC are still scarce, nevertheless.

In this work, we chose 48 differentially expressed NMRDEGs after a thorough analysis of the most recent online GeneCards database and the NMRGs provided by the literature. 33 NMRDEGS associated with HCC were found after looking at the expression profiles of these NMRDEGS and the differentially expressed genes of HCC patients in the TCGA database. Then, we conducted GO and KEGG analyses on 33 NMRDEGs in order to look into the potential molecular mechanism of NMRDEGS. The results revealed that NMRDEGs are primarily enriched and controlled in the following ways: the inhibition of the NF-kappaB transcription factor, the hypoxia response, the circadian rhythm, the inhibition of the oxidative stress response, the nucleotide metabolism, the pyrimidine metabolism, the purine metabolism, the co-factor biosynthesis pathway, and other pathways. Energy metabolism and REDOX processes are connected to these important pathways. It has long been understood that cancer cells behave metabolically differently from the tissues from which they originated.26 As cancer cells undergo carcinogenesis and metastasis, their metabolic activity surges, leading to an uptick in ROS production. This in turn triggers signaling pathways bolstering cancer cell survival, proliferation, and metabolic flexibility. In response to elevated ROS levels, tumor cells amplify their antioxidant defenses, further fueling cancer growth.27 Therefore, we speculate that niacinamide may influence the activity of cancer cells by influencing their capacity for metabolism and REDOX reactions, ultimately influencing the development of tumors.

By using LASSO regression analysis, we discovered six key genes—NAXE, NADSYN1, NT5C, NT5C3A, PNP, and NT5E—associated with prognosis, allowing us to investigate the connection between nicotinamide and HCC. Previous studies have demonstrated a connection between NAXE and the development of HCC.2829 In HCC tissues and cell lines, NAXE expression is markedly downregulated, and this enhances HCC invasion and metastasis both in vitro and in vivo.30

As an anti-inflammatory mediator, NT5C3A prevents the IFN response’s chromatin modifications and NF-B-mediated gene expression, which prevents the release of cytokines. Inhibiting NT5C3A boosts the production of IL-8 in response to TNF stimulation, activates an anti-inflammatory pathway, and acts as a counterbalance to inflammatory cytokine signaling.31 Cancer is known to be characterized by inflammation, which greatly encourages the emergence and growth of malignant tumors.32 The function of NT5C3A in HCC and other cancers has not been fully understood, though. The increased expression of NT5C3A may inhibited the occurrence and growth of HCC through restraining inflammatory reactions.

Research on PNP’s connection to diseases and our understanding of it are still relatively young. Nicotinamide riboside (NR), which is delivered into cells, is broken down by PNP, causing a buildup of NAM. When PNP activity is inhibited, NR can be used to synthesize NAD, and cancer is linked to dysregulation of NAD+-dependent metabolism and signaling.33 As a result, it is thought that PNP may slow the development of cancer via influencing NAD+ production. In this work, we discovered several previously unidentified genes in HCC, indicating their potential as biomarkers for the disease. We found the levels of the key genes NADSYN1 and NT5C are up-regulated in HCC tissues, and the roles of NADSYN1 and NT5C in solid tumors are still poorly studied. To solidify the link between the expression levels of these crucial genes and HCC onset, we examined the correlation between the expression of these six key genes and clinical prognosis. We observed that higher expression of these genes correlated significantly with reduced survival chances, potentially aiding in predicting HCC patient prognosis. Thus, these six genes hold both diagnostic and prognostic value for HCC.

We employed the CIBERSORT algorithm to discern statistically significant differences in the expression levels of nine immune cell types between high and low NAM score groups. This underscores that the NM degree is intricately linked to the immune infiltration in HCC cells, emphasizing the pivotal role of NMRDEGS in regulating the recruitment and invasion of HCC immune cells. Further preclinical and clinical studies are necessary before determining whether it could serve as a therapeutic target. The six key genes may have an impact on the immune milieu in the HCC microenvironment and offer points for determining which patients would respond favorably to HCC immunotherapy. The balance between intracellular activities (internal metabolites, ROS, reducing/oxidizing substrates) and external signals (growth factors, nutritional availability) determines the metabolic pathways that control immune cell maintenance and activation. Growth factors and nutrient availability are environmental factors, whereas internal metabolites, ROS, and reduction/oxidation fifth factors are intracellular activities.34 Hepatocellular carcinoma (HCC) is often associated with chronic inflammatory processes and is considered as a model immunogenic cancer. Dysregulation in the immune microenvironment is considered a key feature of HCC.35 This imbalance involves a variety of immune regulatory mechanisms, which play a central role in tumor formation, proliferation and progression, such as genetic modification, promoting cell proliferation, enhancing drug resistance, inhibiting apoptosis, interfering with genome stability, and shaping the microenvironment. According to the analysis of our data, the infiltration of M1 and M2 type macrophages was reduced in the NAM high expression group compared to the NAM low expression group, while the infiltration of M0 type macrophages was increased.36 M0 macrophages can transform into M1 macrophages in the presence of certain immune activating factors such as LPS, IFNγ and TNFα.34 However, it can be transformed into M2 type under the regulation of IL-4, IL-13, IL-10, IL-33 and TGF-β.37 M1-type macrophages will release proinflammatory cytokines, maintain chronic inflammation, and activate T cells in the initial stages of cancer, and their antitumor effect is partially dependent on the activation of Th1-type immune responses. In the tumor microenvironment, M0 macrophages tend to transition into M2 macrophages. It’s believed that M2 macrophages typically promote tumor progression, possibly through the secretion of immunosuppressive cytokines like IL-10 and TGF-β, and by influencing the recruitment of Treg cells.38 According to our study, the high expression of NAM seems to be able to inhibit the transformation of M0 macrophages to M1 and M2 macrophages, thereby reducing their inflammatory and tumorigenic effects in the tumor immune microenvironment, which may inhibit the development of cancer.27 Moreover, our findings highlighted a pronounced increase in Treg cell infiltration coupled with a marked decrease in monocyte infiltration within the high-score NAM group. Such patterns underscore NAM’s potential critical influence on tumor immune regulation. Significantly, mast cells, particularly when in their resting phase, are recognized for their strong association with the onset and advancement of HCC.39 From our observations, elevated NAM expression seems to decrease the infiltration abundance of mast cells in their resting state, potentially enhancing the prognosis for HCC patients. Dendritic cells (DCs) are pivotal in orchestrating T-cell immune reactions, with the antigen presentation by resting DCs to CD8+ T cells potentially elevating tolerance.40 Furthermore, metabolism stands as a cornerstone in the dynamic between immune and cancer cells within the tumor milieu, influencing both immune surveillance and evasion mechanisms.41 In sum, our research underscores NAM’s paramount influence on the immune landscape of HCC.

Through our experiments, we have confirmed NAM’s suppressive impact on HCC cell proliferation and migration. The metabolic interplay between cancer cells and adjacent immune cells can shape the intensity and nature of the immune response, underscoring the potential role of metabolic interactions in immune surveillance and evasion. As tumors evolve, cancer cells persistently take in nutrients to fuel their rapid growth. Simultaneously, the generation of immunosuppressive metabolites can modulate the survival and functionality of immune cells, tilting the balance towards immune evasion and further tumor advancement. NAM and ROS and reduction/oxidation go hand in hand. Therefore, alterations in immune metabolism may also be affected by aberrant expression of NAM-related genes, and the NAM metabolism may dynamically control variations in immune cell function. Using the GDSC, CCLE, and CellMiner databases, we discovered medications that interact with these six key genes. Future research must be done on the functions of these medications in HCC.

Our research offers novel perspectives into the metabolic pathways underlying HCC by presenting the first thorough examination of the association between NM and HCC. The use of data from public sources without any kind of internal validation is one of the study’s weaknesses which could result in biased interpretations of the findings. Additionally, only a small number of clinical samples were included in the validation experiments, necessitating a larger sample size in order to better examine the connection between key genes and disease diagnosis and prognosis. Furthermore, nicotinamide is not only related to the six genes we are studying (NAXE, NADSYN1, NT5C, NT5C3A, PNP, NT5E), but may also be related to other nicotinamide related genes. We are also aware that niacinamide may participate in cellular metabolism and function through multiple pathways. The effects of other nicotinamide related genes may also have an impact on HCC cells. There is still much to learn about the key genes revealed in this study’s methods of action, necessitating additional study. Finally, a thorough examination of the clinical traits reflected by the various subtypes and phases is not possible due to a lack of clinical data on HCC patients in the database.

Conclusion

In conclusion, six NMRDEGs associated with HCC were discovered to be highly associated with the prognosis of HCC. This signature gene model improves our understanding of HCC at the molecular level and may be very useful in determining the diagnosis, course of treatment, and prognosis of HCC patients.

Acknowledgments

This work was supported by grant from the National Natural Science Foundation of China (81672420, 81702406, U21A20419 and 81372563), the Natural Science Foundation of Guangdong Province of China (2016A030310207) and by grant from the Guangdong Provincial Key Laboratory of Construction Foundation (2017B030314030). This paper has been uploaded to ResearchSquare as a preprint: https://www.researchsquare.com/article/rs-3487654/v1.

Ethics Approval and Consent to Participate

Tissue samples were collected from patients after obtaining informed consent in accordance with a protocol approved by the Ethics Committee of Sun Yat-sen Memorial Hospital (Guangzhou, China; approval number 2023001498). The study was carried out in accordance with the Declaration of Helsinki.

Author Contributions

All authors made a significant contribution to the work reported, whether that is in the conception, study design, execution, acquisition of data, analysis and interpretation, or in all these areas; took part in drafting, revising or critically reviewing the article; gave final approval of the version to be published; have agreed on the journal to which the article has been submitted; and agree to be accountable for all aspects of the work.

Disclosure

The authors report no conflicts of interest in this work.

References

1. Yang JD, Hainaut P, Gores GJ, Amadou A, Plymoth A, Roberts LR. A global view of hepatocellular carcinoma: trends, risk, prevention and management. Nat Rev Gastroenterol Hepatol. 2019;16(10):589–604. doi:10.1038/s41575-019-0186-y

2. Jung M, Lee K-M, Im Y, et al. Nicotinamide (niacin) supplement increases lipid metabolism and ROS-induced energy disruption in triple-negative breast cancer: potential for drug repositioning as an anti-tumor agent. Mol Oncol. 2022;16(9):1795–1815. doi:10.1002/1878-0261.13209

3. Hanahan D. Hallmarks of Cancer: new Dimensions. Cancer Discov. 2022;12(1):31–46. doi:10.1158/2159-8290.CD-21-1059

4. Colaprico A, Silva TC, Olsen C, et al. TCGAbiolinks: an R/Bioconductor package for integrative analysis of TCGA data. Nucleic Acids Res. 2016;44(8):e71–e71. doi:10.1093/nar/gkv1507

5. Zhang Y, Guo L, Dai Q, et al. A signature for pan-cancer prognosis based on neutrophil extracellular traps. J Immunother Cancer. 2022;10(6):e004210. doi:10.1136/jitc-2021-004210

6. Ritchie ME, Phipson B, Wu D, et al. limma powers differential expression analyses for RNA-sequencing and microarray studies. Nucleic Acids Res. 2015;43(7):e47–e47. doi:10.1093/nar/gkv007

7. Wang H, Huo X, Yang X-R, et al. STAT3-mediated upregulation of lncRNA HOXD-AS1 as a ceRNA facilitates liver cancer metastasis by regulating SOX4. Mol Cancer. 2017;16(1):136. doi:10.1186/s12943-017-0680-1

8. Davis S, Meltzer PS. GEOquery: a bridge between the Gene Expression Omnibus (GEO) and bioConductor. Bioinformatics. 2007;23(14):1846–1847. doi:10.1093/bioinformatics/btm254

9. Ivanovska I, Zhang C, Liu AM, et al. Gene signatures derived from a c-MET-driven liver cancer mouse model predict survival of patients with hepatocellular carcinoma. PLoS One. 2011;6(9):e24582. doi:10.1371/journal.pone.0024582

10. Chen Y-L, Wang T-H, Hsu H-C, Yuan R-H, Jeng Y-M. Overexpression of CTHRC1 in hepatocellular carcinoma promotes tumor invasion and predicts poor prognosis. PLoS One. 2013;8(7):1.

11. Barrett T, Wilhite SE, Ledoux P, et al. NCBI GEO: archive for functional genomics data sets-update. Nucleic Acids Res. 2013;41(Database issue):D991–D995. doi:10.1093/nar/gks1193

12. Fishilevich S, Nudel R, Rappaport N, et al. GeneHancer: genome-wide integration of enhancers and target genes in GeneCards. Database. 2017;2017. doi: 10.1093/database/bax028

13. Mayakonda A, Lin D-C, Assenov Y, Plass C, Koeffler HP. Maftools: efficient and comprehensive analysis of somatic variants in cancer. Genome Res. 2018;28(11):1747–1756. doi:10.1101/gr.239244.118

14. Mermel CH, Schumacher SE, Hill B, Meyerson ML, Beroukhim R, Getz G. GISTIC2.0 facilitates sensitive and confident localization of the targets of focal somatic copy-number alteration in human cancers. Genome Biol. 2011;12(4):R41. doi:10.1186/gb-2011-12-4-r41

15. Harris MA, Clark J, Ireland A, et al. The Gene Ontology (GO) database and informatics resource. Nucleic Acids Res. 2004;32(Database issue):D258–D261. doi:10.1093/nar/gkh036

16. Kanehisa M, Goto S. KEGG: Kyoto encyclopedia of genes and genomes. Nucleic Res. 2000;28(1):27–30. doi:10.1093/nar/28.1.27

17. Hänzelmann S, Castelo R, Guinney J. GSVA: gene set variation analysis for microarray and RNA-seq data. BMC Bioinf. 2013;14(1):7. doi:10.1186/1471-2105-14-7

18. Subramanian A, Tamayo P, Mootha VK, et al. Gene set enrichment analysis: a knowledge-based approach for interpreting genome-wide expression profiles. Proc Natl Acad Sci U S A. 2005;102(43):15545–15550. doi:10.1073/pnas.0506580102

19. Kang J, Choi YJ, Kim I-K, et al. LASSO-based machine learning algorithm for prediction of lymph node metastasis in t1 colorectal cancer. Cancer Res Treat. 2021;53(3):773–783. doi:10.4143/crt.2020.974

20. Chen B, Khodadoust MS, Liu CL, Newman AM, Alizadeh AA. Profiling tumor infiltrating immune cells with CIBERSORT. Methods Mol Biol. 2018;1711:243–259.

21. Yang W, Soares J, Greninger P, et al. Genomics of Drug Sensitivity in Cancer (GDSC): a resource for therapeutic biomarker discovery in cancer cells. Nucleic Acids Res. 2013;41(Database issue):D955–D961. doi:10.1093/nar/gks1111

22. Nusinow DP, Szpyt J, Ghandi M, et al. Quantitative proteomics of the cancer cell line encyclopedia. Cell. 2020;180(2):387–402.e16. doi:10.1016/j.cell.2019.12.023

23. Tlemsani C, Pongor L, Elloumi F, et al. SCLC-cellminer: a resource for small cell lung cancer cell line genomics and pharmacology based on genomic signatures. Cell Rep. 2020;33(3):108296. doi:10.1016/j.celrep.2020.108296

24. Park SY, Lee KB, Lee M-J, Bae S-C, Jang -J-J. Nicotinamide inhibits the early stage of carcinogen-induced hepatocarcinogenesis in mice and suppresses human hepatocellular carcinoma cell growth. J Cell Physiol. 2012;227(3):899–908. doi:10.1002/jcp.22799

25. Chen AC, Martin AJ, Choy B, et al. A Phase 3 randomized trial of nicotinamide for skin-cancer chemoprevention. N Engl J Med. 2015;373(17):1618–1626. doi:10.1056/NEJMoa1506197

26. Nault JC, Ningarhari M, Rebouissou S, Zucman-Rossi J. The role of telomeres and telomerase in cirrhosis and liver cancer. Nat Rev Gastroenterol Hepatol. 2019;16(9):544–558. doi:10.1038/s41575-019-0165-3

27. Hsu PP, Sabatini DM. Cancer cell metabolism: Warburg and beyond. Cell. 2008;134(5):703–707. doi:10.1016/j.cell.2008.08.021

28. Kremer LS, Danhauser K, Herebian D, et al. NAXE mutations disrupt the cellular NAD(P)HX repair system and cause a lethal neurometabolic disorder of early, childhood. Am J Hum Genet. 2016;99(4):894–902. doi:10.1016/j.ajhg.2016.07.018

29. Van Bergen NJ, Walvekar AS, Patraskaki M, Sikora T, Linster CL. Christodoulou J.Clinical and biochemical distinctions for a metabolite repair disorder caused by NAXD or NAXE deficiency. J Inherit Metab Dis. 2022;45(6):1028–1038. doi:10.1002/jimd.12541

30. Trinh J, Imhoff S, Dulovic-Mahlow M, et al. Novel NAXE variants as a cause for neurometabolic disorder: implications for treatment. J Neurol. 2020;267(3):770–782. doi:10.1007/s00415-019-09640-2

31. Bogusławska DM, Skulski M, Bartoszewski R, et al. A rare mutation (p.F149del) of the NT5C3A gene is associated with pyrimidine 5’-nucleotidase deficiency. Cell Mol Biol Lett. 2022;27(1):104. doi:10.1186/s11658-022-00405-w

32. Parrinello S, Coppe JP, Krtolica A, Campisi J. Stromal-epithelial interactions in aging and cancer: senescent fibroblasts alter epithelial cell differentiation. J Cell Sci. 2005;118(Pt 3):485–496. doi:10.1242/jcs.01635

33. Zhang B, Chen M, Zhang Y, Chen W, Zhang L, Chen L. An ultrasonic nanobubble- mediated PNP/fludarabine suicide gene system: a new approach for the treatment of hepatocellular carcinoma. PLoS One. 2018;13(5).1.

34. Llovet JM, Castet F, Heikenwalder M, et al. Immunotherapies for hepatocellular carcinoma. Nat Rev Clin Oncol. 2022;19(3):151–172. doi:10.1038/s41571-021-00573-2

35. Donne R, Lujambio A. The liver cancer immune microenvironment: therapeutic implications for hepatocellular carcinoma. Hepatology. 2023;77(5):1773–1796. doi:10.1002/hep.32740

36. Scatozza F, Moschella F, D’Arcangelo D, et al. Nicotinamide inhibits melanoma in vitro and in vivo. J Exp Clin Cancer Res. 2020;39(1):211. doi:10.1186/s13046-020-01719-3

37. Zapata-Pérez R, Wanders RJA, van Karnebeek CDM, Houtkooper RH. NAD+ homeostasis in human health and disease. EMBO Mol Med. 2021;13(7):e13943. doi:10.15252/emmm.202113943

38. DeBerardinis RJ, Chandel NS. Fundamentals of cancer metabolism. Sci Adv. 2016;2(5):e1600200. doi:10.1126/sciadv.1600200

39. Zhang H, Sun L, Hu X. Mast cells resting-related prognostic signature in hepatocellular carcinoma. J Oncol. 2021;2021:4614257. doi:10.1155/2021/4614257

40. Probst HC, McCoy K, Okazaki T, Honjo T, van den Broek M. Resting dendritic cells induce peripheral CD8+ T cell tolerance through PD-1 and CTLA-4. Nat Immunol. 2005;6(3):280–286. doi:10.1038/ni1165

41. Xie F, Bai Y, Yang X, et al. Comprehensive analysis of tumour mutation burden and the immune microenvironment in hepatocellular carcinoma. Int Immunopharmacol. 2020;89(Pt A):107135. doi:10.1016/j.intimp.2020.107135