A Pyroptosis-Related LncRNA Signature for Predicting Prognosis and Treatment in Hepatocellular Carcinoma


 Background: Hepatocellular carcinoma (HCC), which carries a very bad prognosis, is a common malignant tumor worldwide. This study aim to identify a pyroptosis-related long non-coding RNA(pyLncRNA) prognostic signature in HCC by integrated bioinformatics analysis. Methods: All expression profiles of HCC were obtained from The Cancer Genome Atlas (TCGA) and pyroptosis-related genes were from the GSEA website. After identified differentially expressed pyLncRNAs, univariate Cox regression and Lasso analysis were used to identify a pyroptosis-related LncRNAs prognositic signature(py-LPS). Internal validation was used to validate the prognostic value of the py-LPS via the Kaplan-Meier(K-M) curve and receiver operating characteristic(ROC) curve. Additional, we established the nomogram and analyzed the correlation between the signature and immune immune infiltration as well as clinical treatment. Result: 7 pyLncRNAs were established the signature for HCC prognosis. K-M curves exhibited the low risk group presented a markedly longer OS than the high. Clinical subgroups analysis based age, gender, grade and stage yielded the similar results. The signature had an independent prognostic value for HCC(p<0.001). Nomogram estimated one-, three- and five-year survival were 0.777, 0.741 and 0.709. Then, gene set enrichment analysis(GSEA) demostrated significant pathways. Futhermore, we found immune cell infiltration and immunotherapy targets was associated with the signature,which could provided clinical recommendations for chemotherapy.Conclusion: In this study, a novel pyroptosis-related LncRNAs porgnostic signature of HCC, correlated with immune infiltration, could predict the survival of HCC patients and give suggestions for clinical treatment.


Background
Hepatocellular carcinoma (HCC), which accounts for about four-fths of primary hepatic carcinoma, is the fourth most common cause of cancer-related deaths in the world with its incidence rising (1,2). The prognosis of hepatocellular carcinoma is extremely poor. According to a survey of the American Cancer Society , the ve-year survival for hepatocellular carcinoma is less than 20% (2). Therefore, there is an urgent need for predicting the prognosis of patients with HCC and directing clinical treatment by identifying a robust biomarker.
Pyroptosis, a form of programmed cell death, is encoded by genes and associated with in ammatory cells (3). Pyroptosis usually occurs in infected cells and it not only can promote in ammation but plays an essential role in driving innate immune cells to in ltrate injury and infection sites (4). There are three main ways to cause pyroptosis, including classical pathways(mediated by caspase-1 and Cysteine aspartase -1), non-classical pathways(dependent on caspase-4/5/11) and pyroptosis which medicated caspase-3 dependent on tumor drug (5). Among them, the classical pathway cleaves Gasdermin D (GSDMD) into N-terminal fragments through activated caspase-1 and transports them to the cell membrane, mediating cell perforation, and leading to extracellular in ltration, cell swelling and pyroptosis (4,5). The non-classical pathway can not only cause pyroptosis by cleaving GSDMD by activated caspase-4/-5/-11, but also cause cell perforation and pyroptosis by activating the pAnnexin-1-ATP-P2X7 pathway (5). Caspase-3 can cleave the related protein DFNA5 of GSDMD to produce necrotic Nterminal fragments of DFNA5, and the N-terminal fragment of DFNA5 has similar functional activities as the N-terminal fragment of GSDMD, thereby causing cell perforation and pyroptosis (6). In recent years, pyroptosis, which has an association with cancer prognosis, has been found to in uence the occurrence and development of some cancers (7)(8)(9).
Long non-coding RNAs (LncRNAs) are regulatory non-protein coding RNAs longer than 200 nt, which are considered to play an critical role in some cellular processes, such as cell cycle, differentiation, metabolism, disease progression, infection, etc (10). Existing studies have demostrated that LncRNAs interacted with pyroptosis(11).
There is a study which established a predictive signature of pyroptosis-related LncRNAs(pyLncRNA) to predict the survival of lung adenocarcinoma (12). But it has not been discovered in hepatocellular carcinoma. Here, we aimed to nd the LncRNAs related to pyroptosis in hepatocellular carcinoma through bioinformatics methods. Meanwhile, we planned to establish a pyLncRNA-related prognostic signature(py-LPS) of for the rst time, which could predict the prognosis and direct the clinical treatment of hepatocellular carcinoma.

Identi cation of prognostic pyroptosis-related LncRNAs
To analyze the differential expression of LncRNAs between HCC and control groups, we used the "limma" package in the R software to obtain differentially expressed pyroptosis-related LncRNAs (DEpyLncRNAs) with |log2FC| > 1 and false discovery rate(FDR)<0.05 as lter conditions. After merging the DEpyLncRNAs and survival data, univariate Cox regression analyses was performed on them to obtain prognosis-related pyLncRNAs of HCC. Based on the prognostic pyLncRNAs of HCC, we drew the heatmap to show their expressions.
Establishment and validation of the pyroptosis-related LncRNAs prognositic signature(py-LPS) By the ratio of 1:1, the entire cohort(n=370) of HCC were separated into two cohort at random, including training and validation cohorts. By "glmnet" package in R software, we conducted LASSO regression analysis to reduce the over tting genes and identify more meaningful prognostic variables. Meanwhile, a py-LPS was established by LASSO regression outcomes. The risk score of trainng and validation sets were calculated as follows: risks core=∑coef i *exp i (coef: coe cient index, exp: expression level of gene, i: the number of signature genes). Then, the samples was classi ed into high and low risk groups based on the median of risk scores.
Using R software, we conducted Kaplan-Meier(K-M) survival analysis between the high and the low risk group. Then we draw the ROC curves. Predictive value of the py-LPS could be veri ed by the validation and entire cohorts. To explore the applicability of the py-LPS to patients with different clinical characteristics, we conduct K-M survival analysis on different subgroups including age, gender, tumor stage and grade. Additionally, univariate and multivariate Cox analyses were used to validate independent prognosis value. Moreover, we draw a nomogram gure, ROC and calibration curve.

GSEA
According to GSEA`s guide, the gene set conforming to normal p-value<0.05 and | NES |>1 are signi cant (24). Based on the Molecular Signatures Database v7.4 from https://www.gseamsigdb.org/gsea/msigdb/, the expressed gene sets of high/low risk group were analyzed for exploring the difference in pathway enrichment between the two groups. The 5 pathways with the most obvious enrichment in the high/low risk group were exhibited.

Immune In ltration Analysis
Based on the "estimate" package in R software, we calculated the stromal score, immune score and estimate score in tummor microenvironment(TME) of HCC. Then we analyzed the differences and correlation between TME and risk scores. Next, using the "GSVA" and "GSEABase"package of R software, single-sample gene enrichment analysis (ssGSEA) was been performed to explore immune in ltration, including the expression levels of 16 in ltrating immune cells and 13 immune pathways. Additionally, the correlation between immune in ltration and the prognosis of HCC was analyzed.
Signi cance of the py-LPS in clinical treatment Molecular targeting and immune checkpoint blocking therapy have been widely used in the treatment of HCC. PD1, PD-L1, and CTLA4 have been found to affect the prognosis of HCC (41). So we analyzed the difference in expression of these genes. Moreover, based on the TCGA-LIHC RNA-Seq data, "pRRophetic" package in R was performed to calculate IC 50 of typical chemotherapy drugs for hepatocellular carcinoma. Besides, we detected the difference in the IC50 between the two groups based on Wilcoxon test.
Statistical analysis R software(4.1.0) and Strawberry-Perl(64-bit) were utilized for statistical analysis.
Survival curve was depicted by using the K-M analysis and differences were compared by using the logrank test. Additionally, univariate and multivariate Cox analyses were conducted. ROC curves were plotted for evaluating the py-LPS. A p < 0.05 were considered statistically difference.

Results
The steps of this study was presented in Figure 1.

Identi cation of the py-LPS
Base on TCGA-LIHC RNA-seq, 492 differentially expressed pyroptosis-related LncRNAs (DEpyLncRNAs) between HCC and normal tissues were identi ed by difference analysis (Figure 2A). Then, 97 pyLncRNAs were extracted by univariate Cox analyses( Figure 2B). The expressions of 97 prognostic pyLncRNAs was shown in Figure 2C.

Establishment and validation of the prognositic signature
To identify the most powerful prognostic LncRNAs, we conducted Lasso Cox regression analysis on 97 LncRNAs( Figure 2D-E). The result demonstrated that 7 LncRNAs could predict prognosis of HCC powerfully and the coe cient of each LncRNA was shown in Figure 2F. Thereafter, the risk scores were calculuated with the following equation: riskscore= AC131009.1* 0.0772+ AC099850.3* .0735+ LINC01224* 0.1774+ AL031985.3* 0.0489+ AC116025.2* 0.2581+ MKLN1-AS* 0.3898+ AC074117.1* 0.0238. According to the median risk score, all samples could be seprated into the high and low risk groups. Kaplan-Meier curve( Figure 3A) exhibited the low risk group presented a markedly longer OS than the high in the three cohorts. In addition, area under uhe receiver operating characteristic curvea (AUCs) regarded to 1-year-survival were 0.802, 0.788 and 0.789 in training, validation and entire cohort( Figure   3B). Moreover, Figure 3C-D presented the distrubtion of risk scores and survival status with HCC samples.

Clinical subgroups analysis and independent prognosis analysis
The heatmap illusrated that 7 LncRNA expression increased with the growth of risk score( Figure 3E). Meanwhile, the heatmap( Figure 4A-B) revealed that clinicopathologic characteristics(including grade and stage) were related with risk score after removed the samples with incomplete clinical features. Subgroups survival analysis was performed for validating the validation of the py-LPS in different clinical characteristics. HCC samples were separated into 8 subgroups, including age<65, age≥65, female, male, grade 1-2, grade 3-4, stage I-II and stage III-VI. Figure 4C revealed that the low risk group had a signi cantly longer OS than the high in each subgroup, except in patient with grade 3-4. For better exploring the signi cance of our pyLncRNAs-related signature in independently predicting prognosis, both univariate and multivariate analysis were conducted on entire cohort. Result of univariate together with multivariate analysis revealed risk score was an independent risk factor with HCC patients (Figure 4D-F).

Establishment and veri cation of nomogram
Based on py-LPS, a nomogram( Figure 5A) combined risk score and clinicopathologic characteristics(age, gender, grade and stage) was constructed to estimate survival with hepatocellular carcinoma. The ROCs for estimating 1-, 3-and 5-year survival were 0.777, 0.741 and 0.709( Figure 5B). Calibration curves( Figure  5C) shown veri ed the accuracy of nomogram, and indicated the stability for predicting the survival with HCC patients.

GSEA analysis
Based on Principal components analysis(PCA), high and low risk groups were classi ed apparent into 2 clusters( Figure 6A-B). Next, GSEA was implemented to explore different pyroptosis path between the two groups. The most signi cant ve pathways with enrichment in each group were shown in Figure 6C-D.

Tumor microenvironment analysis
According to the ESTIMATE algorithm, stromal score have a range from -1622.33 to 1180.26, immune score from -861.77 to 3157.28, and estimate score from -2465.59 to 3722.93. Moreover, the average stromal score and estimate score in the low risk group was both higher than the high( Figure 7A).
Additionally, gure 7B demostrated that stromal score and estimate score reduced with rising risk score.

The relationship between immune in ltration and py-LPS
The immune in ltration levels were estimated for exploring the connection of py-LPS with tumor immune microenvironment. Among 16 immune cells, aDCs(activated dendritic cells), Macrophages, Th2 Figure 7C). Likewise among 13 immune pathways, MHC class I might be activated in the high risk group while Cytolytic activity, Type I IFN Reponse , Type II IFN Reponse might be activated in the low risk group( Figure 7D). For revealing potential association of immune in ltration with prognosis, we conducted Kaplan-Meier survival analysis. Figure 7E illustrate that differences of inmmune cell in ltration levels(including B cells, Mast cells, Neutrophils, Macrophages, NK cells, TIL and Cytolytic activity were associated with prognosis, as well as differences of immune pathway activation(such as Cytolytic activity).

cells(type-2 T helper cells) and Treg(regulatory T cell) were more active expression in the high risk group while B cells, DCs(dendritic cells), Neutrophils, Mast cells, pDCs(plasmacytoid dendritic cells), NK cells(natural killer T cells) and TIL(tumor-in ltrating lymphocytes) more active expression in the low risk group(
Signi cance of py-LPS in clinical treatment Figure 8A showed 3 immunotherapy targets of HCC had a higher expression in the high risk group. Figure   8B demostrate the association between immunotherapy targets and py-LPS. Hereafter, we accessed IC50 of six drugs for chemotherapy and immunotherapy in the two groups. The results( Figure 8C) demonstrated that 3 drugs(including Doxorubicin, Gemcitabine and Mitomycin.C) had a higher IC50 in the high risk group while only Sorafenib had a higher IC50 in the low risk group.

Disscusion
Hepatocellular carcinoma is a malignancy with an undesirable prognosis. In the past, clinical indicators(including residual liver function, total bilirubin and presence of portal hypertension), biomarkers(such as alpha-fetoprotein, osteopontin) and tumor characteristics were often used to predict its prognosis. Unfortunately, there is still a need for a more powerful prognostic signature (13). Pyroptosis is a programmed cell death with lytic and in ammatory (14). After cell swelling and lysis in pyroptosis, the pro-in ammatory factors released promote the occurrence and development of tumors (14). By affecting the expression of pyrotopia related proteins, the pyrotopia of cells can be promoted and tumor growth can be restricted (7,15). In hepatocellular carcinoma, LncRNA SNHG7 knockdown reduced the expression of SIRT1 by targeting mir-34a /SIRT1 axis but improved the expression of NLRP3 as well as caspase-1 and interleukin 1β, causing the pyroptosis (11). Besides, increasing researches have proved the interactions between various LncRNAs and cell pyrophosis (16)(17)(18). Therefore, this study established a noval pyroptosis LncRNA-related signature, which could strengthen the prediction with the prognosis and guide the treatment of HCC.
In this study, 97 prognosis-related pyLncRNA were identi ed by univariate Cox analyses. Then, 370 patients with HCC were ramdomly divide into training and validation cohort. Subsequently, the py-LPS was constructed based on LASSO regression analysis and it could accurately recognize patients in the high/low risk group. Notably, the high risk group presented a markedly shorter OS compared to the low. So it was able to predict the survival with HCC in different cohorts. Among different clinical subgroups, the py-LPS was proved it a clinically useful predictor. After identi ed as an independent predictive factor for prognosis, risk score was used to establish nomogram with clinical characteristics for enhance the prediction of HCC prognosis.
GSEA interprets gene expression information via focusing on gene sets participate in common biological processes (24). GSEA revealed the 5 most signi cant pathways in the high score group, namely "oocyte meiosis", "cell cycyle", "ubiquitin mediated proteolysis", "spliceosome" and "homologous recombination". Existing studies (25) revealed that LncRNA could in uened cell cycle progression via promoting expression of genes controlling cell proliferation. In addition, the other LncRNAs can in uence the selective shearing of precursor mRNAs by regulating splicerons and participate in the progression of cell cycle (26). Ubiquitin ligases could promote the ubiquitination of key proteins to regulate cellular proliferation and differentiation (27). In pancreatic cancer, Linc01232 enhances tumor metastasis via restraining ubiquitin-mediated degrdndegradation of HNRNPA2B1 (28). Homologous recombination is a way of repairing DNA when two double helix strands are damaged at the same time by obtaining genetic information in the form of sister chromatids or homologous chromosomes (29). Miller revealed that Olaparib had different therapeutic effects in tumor patients with homologous recombination repair de ciency (30). Thus, we considered these pathways were highly related to our signature. TME is the cellular environment surrounding tumors or cancer stem cells, mainly including immune cells and stromal cells (31). ESTIMATE as a new algorithm uses the transcription data of cancer samples to calculate tumor purity. The result of our study revealed risk score had a potential relevance with stromal score while not with immune score. For further exploring the potiental relationship between py-LPS and immune in ltration, single-sample gene enrichment analysis(ssGSEA) was performed. Result illustrated that there were differences in immune in ltration degree(B cells, Mast cells, Neutrophils, Macrophages, NK cells, TIL and ) between the two groups. Based on TCGA sequencing data, Li (32) demostrated that the in ltration degree of immune cells is related to the prognosis of a variety of tumors. Subsequent our survival analysis further veri ed immune cell in ltration was associated with the OS in HCC. CD20+ B cell predicts a higher HCC survival, and activated CD8+ T cells and CD56+ NK cells in TME will affect the B cell content, resulting in enhanced local anti-tumor immune response (33). TIL played a vitial role on antitumor immunity or immune evasion,such as CCR5+ TIN ,which induced an antitumor immune response by releasing IFN-γ(interferon gamma) in bladder cancer (34). Similiarly, tumor-associated macrophages(TAMs) are the main component in TME and it can enhance tumor metastasis by inducing the epithelial-mesenchymal transition(EMT) (35). CD8+ TIL has an anti-tumor effect, and TIL therapy has been used in the treatment of various cancers (36). Besides, zhang`s review (37) concluded that cytolytic activity, depended on NK cells, was activated and inhibitted by the in ltrating immune cells and immune factor of TME, such as Tregs, TAMs , IL-10, TGF-β and etc. However, a previous study illustrated that mast cells could produce IL-17, which accelerated the development of HCC and mast cell is related to the poor prognosis of patients (38). Our result demostrated HCC patients with high mast cell in ltration had longer OS and it was consistent with Yao's study (39). He considered that mast cell might play an anti-tumor role via inducing TILs (39).
Given that immunotherapy is a very effective anti-tumor therapy and universally used for treatment of tumors (40). In this study, the high risk group had a signi cant higher expression of PD-1, PD-L1 and CTLA-4 than the low risk group(p<0.05), which indicated immune checkpoint inhibitors (ICIs) therapy might produce a better bene t in the high risk group. Furthermore, our study demostrated the py-LPS was associated wtih sensitivity of Chemotherapy drugs which commonly used in clinical treatment with HCC. The result illustrated the signature provided clinical recommendations for chemotherapy.
Several limitations appeared in this study. First, there are little studies reports about these candidate LncRNAs in cancer, especially in HCC. The second limitation is that we used internal veri cation instead of external veri cation due to lack of a vailable data about lncRNA sequencing in other external databases. Third, the precise mechanism of pyrophosis regulation of LncRNA in HCC was not explored. The above needs to further experimental demonstration.

Conclusion
In conclusion, a pyroptosis-related LncRNAs porgnostic signature of HCC was rst established in this study. The signature, correlated with immune in ltration, could predict the survival of HCC patients and give suggestions for clinical treatment.