A newly identified pyroptosis-related gene signature for predicting prognosis of patients with hepatocellular cancer.

Background: Hepatocellular carcinoma (HCC) is a very heterogeneous illness, making prognosis prediction a huge problem. Pyroptosis, which has recently been shown to be an inflammatory type of programmed cell death, is involved in HCC. Nevertheless, the role of pyroptosis-related genes in HCC has not been fully elucidated. Thus, this study aimed to construct a prognostic signature based on pyroptosis-related genes for HCC.
Methods: The messenger RNA expression patterns of HCC patients, as well as the accompanying clinical information, were retrieved from The Cancer Genome Atlas (TCGA) database for this research. After differentially expressed pyroptosis-related Gene in tumor and normal groups were identified, Cox regression analyses were performed to construct a prognostic signature which was then assessed through independent prognostic analysis.
Results: A signature consisting of four genes (CASP8, GSDME, NOD2, and PLCG1) was constructed to predict overall survival (OS) for HCC. The signature was identified to be independent by the cox regression analysis and obtained the largest area under the receiver operating characteristic (ROC) curve (AUC) was 0.691, 0.628, and 0.632 for survival at 1, 2, and 3 years, respectively. Conclusions: We discovered that the levels of pyroptosis-related genes expression differed across HCC patients and were associated with both survival and prognosis. This suggested that targeting pyroptosis as a treatment strategy for HCC may be a viable option. 2022 Translational Cancer Research. All rights reserved.

Introduction

On the global scale, primary liver cancer is the sixth commonest form of cancer with a high incidence rate and has been reported to be the fourth major contributor to cancer-related deaths (1,2). Hepatocellular carcinoma (HCC) is the commonest kind of liver cancer (affecting 75–85% of all cases) (3) and is associated with various well-known etiologies, such as nonalcoholic steatohepatitis, types 1 and 2 diabetes mellitus, hepatitis B virus (HBV) infection, obesity, non-alcoholic fatty liver disease, alcohol abuse, and hepatitis C virus (HCV) infections (4). Owing to the multiple etiologic variables associated with HCC, prognosis may be difficult to be predicted. Additionally, HCC has an unfavorable prognosis, with just an 18% 5-year survival probability in the United States (5). Surgical therapies, such as percutaneous ablation, liver transplantation, and liver resection, are often considered to be the most efficacious treatments options for HCC. However, only around 20% of patients with liver cancer are surgical candidates owing to the presence of multiple lesions and extrahepatic metastases. Given the limits of current HCC treatment techniques, novel therapeutic targets are required to enhance HCC clinical outcomes. Hence, credible innovative prognosis models are needed promptly to increase the feasibility of targeted treatments, and there is an added necessity for the creation of new prognostic models.

Pyroptosis is an inflammatory type of programmed cell death characterized by the rupture of the membrane; cytoplasmic swelling; formation of cell membrane pores; and the discharge of cytosolic components including interleukin (IL)-1β into the extracellular environment, which attenuates the systemic or local inflammatory impacts (6). Pyroptosis might have a dual function in the pathogenic mechanisms of cancers. On one hand, as a kind of cell death, it has the potential to inhibit the onset and progression of cancer (7), Additionally, several studies have shown that chemotherapeutic medications, natural products, reagents, and targeted treatment pharmaceuticals may induce or trigger pyroptosis in many types of cancer (8-11). As a consequence, pyroptotic death might be a novel therapeutic target in treating cancer. On the other hand, the many signaling pathways and inflammatory mediators that are activated during pyroptosis have been shown to be intimately associated with tumor progression (12). Although the function of pyroptosis in tumor progression has been noted as research has progressed, the precise activities of pyroptosis in HCC have received less attention. Therefore, we conducted thorough research to evaluate the levels of pyroptosis-related genes expression in HCC patients who were drawn from publicly available databases in order to discover new prognostic biomarkers and targeted therapies for HCC. We present the following article in accordance with the TRIPOD reporting checklist (available at https://tcr.amegroups.com/article/view/10.21037/tcr-22-366/rc).

Methods

Data collection

The Cancer Genome Atlas (TCGA) database (https://portal.gdc.cancer.gov/) was used to retrieve HCC transcriptome data as well as relevant clinical information on the disease. HCC patients’ messenger RNA expression data from 374 tumors together with 50 normal samples and clinical variables [which included tumor-node-metastasis (TNM) stage, grade, gender, clinical stage, and age] were obtained for this study. The study was conducted in accordance with the Declaration of Helsinki (as revised in 2013). The Medical Ethics Committee of Dongyang Hospital affiliated with Wenzhou Medical University exempt our study due to the data of our study came from the public and open TCGA database.

Detection of pyroptosis-related genes with differential expression

provides a list of 33 pyroptosis-related genes acquired from prior reviews. The “limma” R package was employed to analyze the differences in the 33 pyroptosis-related genes expression levels between HCC and normal control samples [false discovery rate (FDR) <0.05, |log2fold change (FC)| ≥0.5]. The differentially expressed genes (DEGs) were notated as illustrated below: * if P<0.05, ** if P<0.01, and *** if P<0.001. To detect gene-to-gene relationships, Pearson correlation analyses were carried out. The heatmap of these genes was plotted using the R software program (ver. 3.8; R Foundation for Statistical Computing, Vienna, Austria).

Table 1

The full names of 33 pyroptosis-related genes in our study

Genes	Full-names
AIM2	Absent in melanoma 2
CASP1	Cysteine-aspartic acid protease-1
CASP3	Cysteine-aspartic acid protease-3
CASP4	Cysteine-aspartic acid protease-4
CASP5	Cysteine-aspartic acid protease-5
CASP6	Cysteine-aspartic acid protease-6
CASP8	Cysteine-aspartic acid protease-8
CASP9	Cysteine-aspartic acid protease-9
ELANE	Elastase, neutrophil expressed
GPX4	Glutathione peroxidase 4
GSDMA	Gasdermin A
GSDMB	Gasdermin B
GSDMC	Gasdermin C
GSDMD	Gasdermin D
GSDME	Gasdermin E
IL-18	Interleukin-18
IL-1β	Interleukin-1 beta
IL-6	Interleukin-6
NLRC4	NLR family CARD domain containing 4
NLRP1	NLR family pyrin domain containing 1
NLRP2	NLR family pyrin domain containing 2
NLRP3	NLR family pyrin domain containing 3
NLRP6	NLR family pyrin domain containing 6
NLRP7	NLR family pyrin domain containing 7
NOD1	Nucleotide binding oligomerization domain containing 1
NOD2	Nucleotide binding oligomerization domain containing 2
PJVK	Pejvakin/deafness, autosomal recessive 59
PLCG1	Phospholipase C gamma 1
PRKACA	Protein kinase cAMP-activated catalytic subunit alpha
PYCARD	PYD and CARD domain containing
SCAF11	SR-related CTD associated factor 11
TIRAP	TIR domain containing adaptor protein
TNF	Tumor necrosis factor

Construction and verification of a prognostic pyroptosis-related gene signature

In order to investigate the predictive significance of the pyroptosis-related genes in the TCGA, we performed a Cox regression analysis of the associations between each gene and survival state. We set 0.05 as the threshold P value and six genes associated with survival were selected for additional investigation. The least absolute shrinkage and selection operator (LASSO)-penalized Cox regression analysis (R package “glmnet”) was employed to build a prognostic model to reduce the risk of overfitting. Finally, the four candidate genes and the coefficients for the four genes were identified, and the penalty parameter (λ) was determined based on the minimal requirements. The computation of the risk score was done after the data from the TCGA expression analysis had been centralized and standardized (using the “scale” package in R) and the risk score equation was as described below: , (X: coefficients, Y: gene expression level).

Following stratification according to the median risk score value, the HCC patients from TCGA were divided into two groups: high- and low-risk groups. Principal component analysis (PCA) was conducted using the “prcomp” utility of the “stats” R package based on the genes expressed in the signature. The “survival”, “survminer”, and “timeROC” R packages were employed to conduct a 3-year receiver operating characteristic (ROC) curve analysis.

Analysis of the prognosis significance of the risk score independently

The clinical data of patients in the TCGA cohort were obtained. Our regression model was constructed using clinical data in conjunction with the risk score. The analyses were conducted using multivariate and univariate Cox regression models.

Functional enrichment analysis

The “clusterProfiler” R package was utilized to conduct Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) analyses based on the DEGs. The “gsva” package was utilized to conduct the single-sample gene set enrichment analysis (ssGSEA) to calculate the infiltrating score of 16 immune cells and the activity of 13 immune-related pathways.

Statistical analysis

In order to evaluate the levels of gene expression between HCC and adjoining normal tissues, a single-factor analysis of variance was utilized, and Pearson’s chi-square analysis was employed to contrast categorical variables. The Kaplan-Meier analysis and the log-rank test were utilized to examine the overall survival (OS) across various groups. Multivariate and univariate Cox regression models were employed to investigate the independent prognostic significance of risk scores. The R software (version 4.0.2) was employed to carry out all of the statistical analyses. As previously stated, unless otherwise stated, a P value <0.05 was considered significant, and all P values were two-tailed.

Results

The TCGA dataset and the characteristics of the patients

The investigation in the present study included 374 samples of HCC tissue and 50 samples of adjoining normal tissue from the TCGA. A sum of 365 HCC samples was included in the study, each having complete survival data with an average age of patients was 61 years. Clinical information encompassed the patient’s clinical stage, gender, grade, age, and TNM stage.

Detection of DEGs between tumor and normal tissues

The levels of 33 pyroptosis-related gene expressions were compared utilizing TCGA data of 374 HCC tissue and 50 adjoining normal tissue samples, and we identified 18 DEGs (all P<0.01). Specifically, two genes (IL-1β and IL-6) were downregulated, and 16 other genes (AIM2, CASP3, CASP8, TIRAP, PYCARD, PLCG1, NOD2, PJVK, NOD1, NLRP7, NLRP1, GSDME, GSDMD, GSDMC, GSDMB, and GPX4) were upregulated in the tumor tissues (). The correlation between these genes is presented in .

Figure 1

The expression of 33 pyroptosis-related genes and the relationships between DEGs. (A) A heatmap of pyroptosis-related genes compared between the normal (brilliant blue) and the tumor tissues (red) is shown. P values are presented as **P<0.01, ***P<0.001. (B) The correlation network of DEGs is represented by a red line for a positive correlation or a blue line for a negative correlation. The darker the color, the stronger the correlation. N, normal tissue; T, tumor tissue; DEGs, differentially expressed genes.

Construction of a gene model for prognosis

A sum of 365 HCC samples was matched with relevant patients whose survival data was complete. We employed univariate Cox regression to examine 18 DEGs, and six genes (GSDME, PLCG1, CASP8, NOD2, NOD1, and CASP3) were selected with P<0.05 as a selection condition (). A LASSO Cox regression model was utilized to determine the genes that were the most highly predictive as prognostic markers for HCC. The λ value was chosen once the median of the total of squared residuals was the least (). Four potential predictors (CASP8, GSDME, NOD2, and PLCG1) were found to be prognostic markers for HCC. The following is the equation for calculating the risk score: risk score = (0.093 × CASP8 expression) + (0.126 × GSDME expression) + (0.311 × NOD2 expression) + (0.056 × PLCG1 expression). According to the median value generated by the risk score equation, 365 patients were categorized into two comparable groups: low- and high-risk groups (). The results of the PCA illustrated that patients in the various risk groups were spread in two distinct directions (). When compared with individuals in the low-risk group, those in the high-risk group exhibited a greater likelihood of dying sooner than those in the low-risk group, as presented in . Similarly, the Kaplan-Meier curve illustrated that high-risk group patients displayed a considerable dismal OS as opposed to the ones in the low-risk group (; P<0.001). Time-dependent ROC analysis was conducted to examine the specificity and sensitivity of the prognostic model, and we discovered that the area under the ROC curve (AUC) was 0.691, 0.628, and 0.632 for survival at 1, 2, and 3 years, respectively ().

Figure 2

Survival analysis and gene selection in the anticipating of prognosis of patients with HCC. (A) Hazard ratios for survival-associated pyroptosis-related genes in HCC were determined using forest plots. (B) A LASSO Cox regression model was utilized to plot partial likelihood deviance against log (λ). (C) The λ parameter represents the coefficients of special characteristics. Each line represents a gene (line 1: NOD1, line 2: CASP8, line 3: GSDME, line 4: CASP3, line 5: NOD2, line 6: PLCG1). HCC, hepatocellular carcinoma; LASSO, least absolute shrinkage and selection operator.

Figure 3

Prognostic evaluation of the four-gene signature model in the TCGA cohort. (A) The risk score distribution and median values in the TCGA cohort (the dashed lines on the left portion are the populations at low risk while those on the right portion are the populations at high risk). (B) PCA plot for HCC on the basis of the risk score. (C) Each patient’s survival status (the dashed lines on the left portion are the populations at low risk while those on the right portion are the populations at high risk). (D) Kaplan-Meier graphs illustrating the OS of patients classified as high- or low-risk groups. (E) ROC curves validated the predictive effectiveness of the risk score. PCA, principal component analysis; AUC, area under the ROC curve; ROC, receiver operating characteristic; TCGA, The Cancer Genome Atlas; HCC, hepatocellular carcinoma; OS, overall survival.

The independent predictive significance of the risk model

The risk score, clinical stage, M stage, and T stage of pyroptosis-related genes all had an impact on the patients’ prognosis, according to univariate analysis (P<0.05). N stage, sex, grade, and (P>0.05) (). The risk score still independently predicted OS in the multivariate Cox regression analysis even after adjusting for additional confounders ().

Figure 4

Analyses of the risk score utilizing multivariate and univariate Cox regression. (A) The forest plot shows the results of the univariate Cox regression analysis in HCC. (B) The forest plot shows the results of a multivariate Cox regression evaluation in HCC (P<0.05). T, tumor; M, metastasis; N, node; HCC, hepatocellular carcinoma.

Functional analyses based on the DEGs

In order to further explore the possible pathways related with pyroptosis in HCC, GO enrichment analysis and KEGG pathway analysis were then performed based on DEGs. The results indicated that the DEGs were mainly correlated with the IL-1β production, inflammasome complex, cytokine receptor binding (). KEGG pathway analysis revealed that these DEGs were mainly enriched in the NOD-like receptor signaling pathway, TNF signaling pathway, NF-kappa B signaling pathway ().

Figure 5

Functional analysis based on the DEGs. (A) Bubble graph for GO enrichment. (B) Bubble graph for KEGG pathways. (The bigger bubble means the more genes enriched, and the increasing depth of red means the differences were more obvious; BP, biological process; CC, cellular component; MF, molecular function; DEGs, differentially expressed genes; GO, Gene Ontology; KEGG, Kyoto Encyclopedia of Genes and Genomes.

Comparison of the immune activity between high-risk group and low-risk group

In order to explore the relationship between pyroptosis-related gene and immunity in HCC, we further compared the enrichment scores of 16 types of immune cells and the activity of 13 immune-related pathways between the low- and high-risk groups by employing the ssGSEA. The high-risk group had decreased infiltrations of the mast cells, natural killer (NK) cells, and plasmacytoid dendritic cells (pDCs) and increased infiltrations of the activated dendritic cells (aDCs), macrophage and Treg compared to the low-risk group (P<0.05, ). Besides, the high-risk group exhibited suppression in the immune pathways, including the cytolytic activity, type I interferon (IFN) response and type II IFN response (P<0.05, ).

Figure 6

Comparison of the ssGSEA scores for immune cells and immune pathways. (A) Comparison of the enrichment scores of 16 types of immune cells between low- (green box) and high-risk (red box) group. (B) Comparison of the enrichment scores of 13 immune-related pathways between low- (blue box) and high-risk (red box) group. P values were showed as: ns, not significant; *P<0.05; **P<0.01; ***P<0.001. aDCs, activated dendritic cells; DCs, dendritic cells; iDCs, immature dendritic cells; NK, natural killer; pDCs, plasmacytoid dendritic cells; Tfh, T follicular helper cell; Th1, T-helper 1; Th2, T-helper 2; TIL, tumor-infiltrating lymphocyte; APC, antigen presenting cell; CCR, C-C chemokine receptor; HLA, human leukocyte antigen; MHC, major histocompatibility complex; IFN, interferon; ssGSEA, single-sample gene set enrichment analysis.

Discussion

We initially examined the messenger RNA levels of 33 known pyroptosis-related genes in HCC and normal tissues and discovered that 18 of them exhibited differential expression. Our research indicated that pyroptosis-related genes were extensively expressed in tumor tissue, suggesting that they perform an instrumental function in anticipating the prognosis of HCC patients. Additionally, pyroptosis-related genes were closely connected with each other in modulatory networks, indicating that they cooperate in the onset and progression of cancer. The LASSO algorithm examined all independent variables and identified the most significant ones (13), making it far more accurate than the conventional regression methods. Four of the 33 genes (CASP8, GSDME, NOD2, and PLCG1) were found to be prognostic variables for HCC utilizing LASSO Cox analysis. ROC curve analysis was conducted to examine the predictive capacity of pyroptosis-related genes in determining the prognosis of HCC patients, and the findings indicated that the four abovementioned possible predictors could have detrimental impacts in HCC patients owing to their correlation with poor dismal survival. We utilized multivariate and univariate Cox analyses to determine risk factors in order to develop a predictive model on the basis of these four genes. Eventually, the signature composed of four genes was validated as an independent predictive factor for HCC.

Pyroptosis, a novel kind of programmed cell death, has emerged as an innovative target of investigation in cancer. Numerous research reports have illustrated a strong correlation between pyroptosis and the onset and progression of many malignancies (14-16). Nonetheless, the extent to which pyroptosis-related genes are associated with patient survival and their involvement in HCC has not been elucidated. In this study, we developed a predictive model integrating four pyroptosis-related genes (CASP8, GSDME, NOD2, and PLCG1) and discovered that it was capable of predicting OS in HCC patients.

Caspase-8 is a critical component of the apoptotic process and is encoded by the CASP8 gene. Fritsch et al. established that caspase-8 functions as the pyroptotic molecular switch (17). The CASP8 gene is a critical cancer susceptible gene, as shown by numerous researches examining a variety of tumor forms, such as ovarian cancer, lung, prostate, and breast (18-20). Soung et al. discovered that the CASP8 gene is commonly altered in HCC and that this mutation could result in the loss of the cell death role of CASP8 gene, thereby leading to the occurrence and progression of HCC (21).

GSDME/DFNA5 (deafness, autosomal dominant 5) is a gasdermin superfamily gene correlated with autosomal dominant nonsyndromic hearing impairment (22). GSDME was recently discovered as a pore-forming protein that is triggered by caspase-3 cleavage, leading to secondary necrosis upon cell apoptosis or initial necrosis without the need for an apoptotic stage, dubbed pyroptosis-like necrosis (23). GSDME has been shown to operate as a potential tumor suppressor gene in a significant proportion of gastric, colorectal, and breast cancers (24,25). Compared to the GSDME low-expression group, the 5-year survival rate of the high-expression group was considerably elevated, indicating that GSDME may serve as a positive prognostic marker in squamous esophageal cancer (26). GSDME was shown to be overexpressed in HCC tissues and to have a negative connection to survival rates, which led us to infer that it was a cancer-promoting gene. In view of the limited available data of TCGA and often contradictory results in different tumors, our results about GSDME provide some insights for further research.

NOD2, an affiliate of the NOD-like receptor family, has been known as an innate immune sensor capable of inducing powerful immune responses against infections (27). A large number of innate immune sensors have been identified to have a significant role in the progression of cancer. NOD2 polymorphisms have been discovered to be correlated with an elevated risk of a variety of malignancies, particularly laryngeal cancer, breast, endometrial, ovarian, colorectal, and gastric (28). Recent research has illustrated that the NOD2 protein engages in the process of cell pyroptosis in sepsis (29). The function of NOD2 in the progression of HCC, on the other hand, is still incompletely understood. As a result, our work serves as a starting point for a thorough investigation.

PLCG1, an affiliate of the enzymes phospholipase C family, participates in the signal transduction pathway controlled by receptor tyrosine kinases, hence regulating cell growth, differentiation, migration, and apoptosis (30). Recent research by Kang et al. found that knocking down PLCG1 decreased GSDMD-N-stimulated cell death and suggested that the function of GSDMD and pyroptosis might be mediated by PLCG1 (31). We discovered that the elevated level of PLCG1 expression was correlated with unfavorable survival outcomes, which could be the result of the negative modulation of pyroptosis that this gene exerts. In summary, of the other geneses in the prognostic model, four genes (CASP8, GSDME, NOD2, and PLCG1) were shown to be mediators of pyroptosis. In the present investigation, these genes were all shown to be elevated in HCC tumor tissue and to be correlated with an unfavorable prognosis. Nonetheless, additional research is needed to determine how these genes interface with one another throughout the process of pyroptosis.

In conclusion, the present study indicated that pyroptosis is significantly associated with HCC since the expression levels of major pyroptosis-related genes were shown to have a differential expression of HCC and normal tissues. The risk scores, which were centered on four pyroptosis-related genes, were shown to be an independent risk factor indicator for anticipating HCC. Using these findings, we have identified a unique gene signature that may be used to anticipate the prognosis of patients with HCC. In clinical practice, we can screen high-risk patients by a unique gene signature, so as to provide personalized treatment plans and improve the prognosis of patients. At the same time, through the research on the cells and animal models of these prognostic markers, the targets with clinical therapeutic effects are screened out, which provides the basis for the targeted therapy of liver cancer. Nevertheless, one limitation of this study is that our study was based on individual sources of TCGA and was not validated from an independent cohort. More research is needed to further clarify these findings.

The article’s supplementary files as