A Pyroptosis-Related Gene Panel in Prognosis Prediction and Immune Microenvironment of Human Endometrial Cancer

As the second common diagnosed cancer among gynecological tumors, endometrial cancer (EC) has heterogeneous pathogenesis and clinical manifestations. Therefore, prognosis prediction that considers gene expression value and clinical characteristics, is helpful to patients with EC. We downloaded RNA expression and clinical data from the TCGA database. We achieved 4 DEPRGs and constructed the PRG panel by univariate, lasso and multivariate Cox analysis. Based on the median value of the risk score, patients were divided into two groups. The Kaplan–Meier curve suggested that the patients with lower risk scores had better clinical outcomes of EC. AUC of ROC curves suggested the panel can be used as an independent predictor. Future analysis indicated the positive correlations between risk score and clinical characteristics. What’s more, we performed GO and KEGG functional analysis and immune environment exploration to get an understanding of the potential molecular mechanism and immunotherapeutic target. To future validate the panel, we found that the relapse-free and overall survival probability of 4 prognostic DEPRGs between high-expression group and low-expression group were different through the Kaplan–Meier plotter in UCEC. In addition, GEPIA database and RT-PCR experiment indicated GPX4 and GSDMD were highly expressed in UCEC compared to normal endometrial tissue, and TIRAP and ELANE were downregulated. This study identified a PRG panel to predict the prognosis immune microenvironment in human EC. Then, Kaplan–Meier analysis and AUC below the ROC curves was used to validate the panel. In addition, Chi-square was used to show the clinical significance. GO, KEGG and GSEA were used to show the functional differences. Different immune-related databases were used to analyze the immune characteristics. The Kaplan–Meier plotter website was used to assess the effect of genes on survival. GEPIA and RT-PCR were used to analyze the expression level. In summary, we identified 4 prognosis-associated pyroptosis-related genes (ELANE, GPX4, GSDMD, and TIRAP). The panel can also predict prognosis prediction and immune microenvironment in human endometrial cancer.


INTRODUCTION
There were 417267 new cases and 97370 new cancer deaths reported in 2020, endometrial cancer (EC) was ranked in the top 10 most common cancers in female patients across the world, and is the second most commonly diagnosed cancer among gynecological tumors (Sung et al., 2021). Most cases are diagnosed in the late stage of life, but there are increasing numbers being diagnosed at an early stage, meaning the onset of EC is younger (Amant et al., 2005;Moore and Brewer, 2017). Though the mortality is relatively lower than other gynecological tumors, the pathogenesis and clinical manifestations of EC are heterogeneous (Gupta, 2017). Therefore, prognosis prediction considering the gene expression value and clinical characteristics is helpful for patients with EC.
Pyroptosis mediated by gasdermin is a kind of programmed cell death, which is characterized by cell swelling until cell membrane rupture, leading to the release of cell contents and then activating a strong inflammatory response (Bergsbaken et al., 2009;He et al., 2015;Vande Walle and Lamkanfi, 2016;Kovacs and Miao, 2017;Man et al., 2017;Shi et al., 2017). For morphology, pyroptosis cells showed swelling under a light microscope, and there were many bubble-like protrusions. Compared with necrotic cells, the swelling degree of focal death is lower. Under the electron microscope, it can be seen that before the rupture of the membrane, the pyroptosis cells form a large number of vesicles, namely focal apoptotic bodies. The cell membrane then forms pores and breaks, meaning the content flows out. The biochemical characteristics of pyroptosis involve the formation of inflammatory bodies, the activation of caspase and gasdermin, and the release of a large number of inflammatory factors. Pyroptosis is also an important innate immune response and plays an important role in the fight against infection (Bergsbaken et al., 2009;Man et al., 2017;Robinson et al., 2019). Pyroptosis has been reported to be involved in different processes in different diseases, including cardiovascular diseases (Jia et al., 2019;Zhaolin et al., 2019), Parkinson's disease , Alzheimer's disease (Han et al., 2020), diabetic kidney disease Lin et al., 2020), inflammatory bowel disease , cancers (Pizato et al., 2018;Cui et al., 2019;Qiao et al., 2019;Yu et al., 2019;Zhang C. C. et al., 2019;Zhang T. et al., 2019;Teng et al., 2020), and so on.
Recently, many reports have outlined that pyroptosis-related genes and processes play an important role in cancers. For breast cancer, DRD2 which promotes macrophage M1 polarization, triggers gasdermin E (GSDME) and induces pyroptosis, which is then downregulated (Tan et al., 2021), and through increased levels of caspase-1, caspase-3, gasdermin D and E (GSDMD and GSDME), docosahexaenoic acid (DHA) and tetraarsenic hexoxide could promote pyroptosis and thus inhibit breast cancer progress (Pizato et al., 2018;An et al., 2021;Li, 2021). For colorectal cancer, lobaplatin (Yu et al., 2019), liver X receptors (LXPs) (Derangère et al., 2014), and arsenic trioxide (ATO) combined with ascorbic acid (AA) (Tian et al., 2020) could be used to better treat colorectal cancers through inducing pyroptosis. For hepatocellular carcinoma, Hage et al. (2019) reported that sorafenib could induce the pyroptosis of macrophages which were the key mediators of antitumor effect.  reported that miltirone, a phenanthraquinone derivative isolated from Radix Salviae Miltiorrhizae, could inhibit hepatocellular carcinoma cells proliferation and induce the proteolysis of GSDME and cleavage of caspase 3. Furthermore, gasdermin also attaches importance to gastric, lung, skin, and esophageal squamous cancer (Xia et al., 2019). Yang et al. (2020) have reported that NLRP3, caspase-1, and GSDMD were upregulated in EC, and activation in pyroptosis could lead to an anti-tumor effect. Since pyroptosis plays such an important role in cancers and to date there have been few studies exploring its role in EC, combined analysis of pyroptosis and EC is imminent.
This study used TCGA, a large-scale open database, for analyzing pyroptosis-related genes (PRGs) in EC and constructed a novel PRG panel to predict the overall survival and immune microenvironment in EC. This research confirmed that the panel can be used to predict the clinical outcomes of EC independently and could eventually instruct immunotherapy of EC in the future.

MATERIALS AND METHODS
The process of data analysis is shown in the flow chart in Figure 1.

Data Acquisition
The RNA expression and clinical data of UCEC patients were downloaded from the TCGA website on April 2, 2021 1 . After removing the cases with incomplete clinical information, 526 cases were remained, as shown in Table 1. The "limma" R package was used to normalize the RNA expression profiles. A total of 34 pyroptosis-related genes were achieved through literature mining and are provided in Table 2 (Man and Kanneganti, 2015;Wang and Yin, 2017;Karki and Kanneganti, 2019;Xia et al., 2019;Ye et al., 2021).

Identification of Prognostic Differentially Expressed Pyroptosis-Related Genes
Differential expression genes (DEGs) were searched for with the "limma" R package by comparing tumor tissues with

Construction and Confirmation of a Prognostic Panel by Differentially Expressed Pyroptosis-Related Genes
The "glmnet" R package was used for constructing the DEPRGs panel by lasso regression analysis, and 10-fold cross confirmation and a P value of less than 0.05 were used to simplify and rectify the model . The gene expression levels should be corrected by demographics in multivariate regression models.
At the same time, the risk score of patients was generated [risk score = (expression value of each gene × and its coefficient)]. We used the median risk score to divide patients into two groups and visualized the results with the "Rtsne" R package. We then used the "stats" R package to conduct a principal components analysis (PCA) to study the different gene expression patterns of samples.
To validate this panel, we conducted a Kaplan-Meier analysis to show the clinical outcome differences between two groups with or without other clinical characteristics and visualized the survival curve using the "survive" and "survminer" R packages. The 1-, 2-, 3-, 4-, and 5-year ROC curves and the ROC curves compared with panels were drawn to predict the survival status with the "survivalROC" R packages. Univariate and multivariate Cox regression analyses then showed whether our panel could predict the prognosis of endometrial cancer independently.
A Chi-square test was then used to analyze the relationship between the panel and clinicopathological characteristics. These were shown in a strip chart and labeled as follows: p < 0.001 = * * * , p < 0.01 = * * , Frontiers in Cell and Developmental Biology | www.frontiersin.org and p < 0.05 = * by "ComplexHeatmap" R package. Wilcoxon signed-rank test compared the risk scores among different groups of these clinicopathological characteristics.
The scatter diagram showed the analysis results.
We then conducted operating characteristic curve (ROC) and decision curve analysis (DCA) to assess the sensitivity and specificity of the prognostic panel for UCEC in comparison to other clinicopathological factors (Vickers and Elkin, 2006). Frontiers in Cell and Developmental Biology | www.frontiersin.org A nomogram was constructed, integrating the prognostic panel, for prediction of 1, 3, and 5-year OS of UCEC patients.

Gene Ontology, Kyoto Encyclopedia of Genes and Genomes Functional Enrichment, and Gene Set Enrichment Analysis Analysis Between High-Risk Versus Low-Risk Groups
Gene ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analyses predicted the potential functions of DEPRGs using the "clusterProfiler" R package. We then sorted and presented the top items of GO and KEGG pathways based on P-value < 0.05 through the statistical software R (version 4.0.2) and visualized the results in barplots using the "ggplot2" R package. Gene Set Enrichment Analysis (GSEA) software (Version 4.0.3) 2 was also used to explore the potential biological function difference between the two groups. GSEA (version 3.0) was run for the "c2.cp.kegg.v.7.2.symbols.gmt" gene sets. The number of permutations was set to 1,000 and the phenotype labels were high-risk and low risk. FDR < 0.25 and NOM P < 0.05 indicated statistical significance.

Immune Characteristics Analysis
Using TIMER, CIBERSORT, CIBERSORT-ABS, QUANTISEQ, MCPcounter, XCELL, and EPIC, we analyzed the correlation 2 http://www.broad.mit.edu/gsea between the riskScore and immune-cell characteristics in UCEC patients of TCGA database by Spearman correlation analysis, and the results were visualized in a diagram. A p value < 0.05 was considered as statistically significant.
The "gsva" R package compared normalized gene expression data with the gene sets that had common biological functions, chromosomal localization, and physiological regulation in ssGSEA (Subramanian et al., 2005). The infiltrating score of 16 types of immune cells and the activity of 13 immune-related pathways were concluded. To study the expression value of immune checkpoint-related genes in two groups, we performed a "ggstatsplot" R package and violin plot visualization.

Survival Analysis of Four Prognostic Differentially Expressed Pyroptosis-Related Genes
The Kaplan-Meier plotter 3 could assess the effect of 54,675 genes on survival using 18,674 cancer samples (Győrffy et al., 2013). These cover 543 UCEC cancer patients with relapse-free and overall survival information (Nagy et al., 2021). Our study analyzed the relapse-free and overall survival of four prognostic DEPRGs through the Kaplan-Meier plotter in UCEC. Patients were classified into two groups according to the median of each prognostic DEPRG expression in the Kaplan-Meier plotter for relapse-free and overall survival. This classification method could show the survival probability differences between the high-expression group and the lowexpression group.

Gene Expression Profiling Interactive Analysis Database Analysis of Four Prognostic Differentially Expressed Pyroptosis-Related Genes
GEPIA is an online database that facilitates the standardized analysis of a tremendous amount of RNA sequencing data in the TCGA and GTEx data sets (Tang et al., 2017)

Ethics Statement
This study was carried out in accordance with the standards of the Helsinki Declaration of the World Medical Association and approved by the Ethics Committee of China Medical University. All clinical samples were collected from the Shengjing Hospital of China Medical University, with informed consent from all patients.

Tissue Collection
There were 56 clinical samples used in this study, including 33 cases of primary UCEC tissue and 23 cases of normal endometrial tissue. All samples were collected from the patients undergoing surgical excision at the Department of Obstetrics and Gynecology, Shengjing Hospital of China Medical University. No patient received radiotherapy, chemotherapy, or hormone therapy before surgery. The histopathological diagnosis was obtained from the Pathology Department according to the criteria of the World Health Organization.

Statistical Analysis
Statistical analysis was conducted using the R package and Graphpad Prism 8. Variance homogeneous and normal distributed continuous variables were analyzed using the unpaired student's t-test, otherwise, the Mann-Whitney U-test or Kruskal-Wallis H-test was used. DEGs were Frontiers in Cell and Developmental Biology | www.frontiersin.org 9 October 2021 | Volume 9 | Article 705828 A P-value of less than 0.05 was considered to be statistically significant (p < 0.001 = * * * , p < 0.01 = * * , and p < 0.05 = * ).

Identification of Prognostic Differentially Expressed Pyroptosis-Related Genes of TCGA Database
In total, 5 of the 34 pyroptosis-related genes were related to prognosis and four of them (ELANE, GPX4, GSDMD, and TIRAP) had different expression values between tumor and normal tissues, which were, thus, chosen as prognostic DEPRGs (Figure 2A). The heatmap in endometrial cancer showed that ELANE and GSDMD were upregulated and GPX4 and TIRAP were downregulated ( Figure 2B). The univariate Cox analysis of OS of the prognostic pyroptosis-related DEGs indicated that the four genes were low-risk factors (HR < 1) ( Figure 2C). Therefore, these four genes were constituted into the optimal prognostic risk model of DEPRGs. The details of the genes are shown in Table 2.

Construction of a Prognostic Panel by Differentially Expressed Pyroptosis-Related Genes
A prognostic panel was constructed using the four DEPRGs by lasso regression analysis (Figures 2D,E) and multivariate Cox analysis ( Figure 2F). The risk score of each patient was calculated according to the formula: risk score = expression value of ELANE × (−2.41387833412429) + expression value of GPX4 × (−0.298331831704336) + expression value of GSDMD × (−0.337701574941294) + expression value of TIRAP × (−0.58979182292574).

Prognostic Confirmation of the Prognostic Panel by Differentially Expressed Pyroptosis-Related Genes
The median risk score was taken for dividing patients into two groups (high-and low-risk group) and visualizing the results, as seen in Figure 3A. The survival status of each group was also visualized, as seen in Figure 3B. Patients in the two groups were distributed in different directions according to the PCA and t-SNE analysis (Figures 3C,D). Then, KM analysis showed that the patients with lower risk scores had higher survival potential with a p value < 0.01 ( Figure 3E). The 1-, 2-, 3-, 4-, and 5-year ROC curves showed all areas under the curve (AUC) were over 0.6, and AUC under 5-year ROC curve was 0.692 ( Figure 3F).

Independence Confirmation and Clinical Evaluation of the Panel by Differentially Expressed Pyroptosis-Related Genes
We carried univariate Cox analyses among the common clinical features first, and multivariate Cox analysis next, which showed our panel can predict the prognosis of endometrial cancer independently (Figures 4A,B).
The chi-square test, analyzing the relationship between the panel and clinicopathological characteristics, suggested that age, clinical stage, and tumor grade were significantly associated with the risk Score (p < 0.001) (Figure 5A). Wilcoxon signed-rank test comparing the risk score differences among different groups of these clinicopathological characteristics indicated age clinical stage, and tumor grade showed a positive relationship to risk scores (Figures 5B-D). The AUC of the panel, traditional clinicopathological characteristics, and DCA are shown in Figures 6A,B. There was also a nomogram incorporating the panel and clinicopathological characteristics in Figure 6C.

Gene Ontology, Kyoto Encyclopedia of Genes and Genomes Functional Enrichment, and Gene Set Enrichment Analysis Analysis Between High-Risk Versus Low-Risk Groups
To understand the molecular mechanism related to the risk score, the DEGs between two groups were used for GO and KEGG pathway analysis (Figures 7A,B). GO functional enrichment indicated that processes related to immunoglobin, immune response, and phagocytosis were enriched ( Figure 7A). KEGG pathways analysis also indicated DEGs are enriched in some classic pathways including Wnt and cAMP signaling pathways and some necessary substance metabolism ( Figure 7B). The details are shown in Tables 3, 4

Immune Characteristics Analysis
Through TIMER, CIBERSORT, CIBERSORT-ABS, QUANTISEQ, MCPcounter, XCELL, and EPIC conjoint analysis, we indicated that patients with high riskScore accumulated tumor-infiltrating immune cells such as macrophage, NK cells, B cells, and T cells (Figure 8A). In addition, we quantified the scores of different immune cell subsets, related functions, or pathways with ssGSEA for analyzing the relationship between risk scores and immune status. Interestingly, the immune function of CCR, check point, cytolytic activity, T cell co inhibition, T cell co stimulation, type I IFN response, type II IFN response and immune cell proportion of CD8+ T cells, DCs, neutrophils, pDCs, T helper cells, Th1 cells, and TIL were significantly different among the two groups (Figures 8B,C). Furthermore, patients with high-risk scores had low expression of immune checkpoint related genes ( Figure 8D).

External Validation of the Four Prognostic Differentially Expressed Pyroptosis-Related Gene Panel
We first analyzed the relapse-free and overall survival of four prognostic DEPRGs through the Kaplan-Meier plotter in UCEC. The results showed that the survival probability differences between the high-expression group and the lowexpression group were consistent with the former analysis (Figures 9A-H). We then used the GEPIA database and RT-PCR experiment to analyze the expression levels of the four prognostic DEPRGs, which indicated that GPX4 and GSDMD were highly expressed in UCEC compared to normal endometrial tissue, and TIRAP and ELANE were downregulated (Figures 10A-H).

DISCUSSION
Recently, pyroptosis related genes, proteins, and related molecular processes in carcinogenesis, cancer novel therapy have attracted much attention in the research. There is a correlation with EC, as reported in another previous study . There is also another PRG panel that can help predict the overall survival and immune microenvironment in gastric and ovarian cancer (Shao et al., 2021;Ye et al., 2021). However, a PRG panel in EC has never been reported.
We used the TAGA database and 34 PRGs were analyzed in cases of UCEC and healthy people. We found four PRGs (ELANE, GPX4, GSDMD, and TIRAP) were significantly differentially expressed and related to overall survival in EC. As a member of the gasdermin family, GSDMD cleaved by caspase-1 and caspase-11/4/5 can trigger pyroptosis (Kayagaki et al., 2015;Shi et al., 2015Shi et al., , 2017. The role of GSDMD has been reported in EC . As a member of serine protease, ELANE can hydrolyze many proteins besides elastin. As Kambara et al. (2018) have reported, ELANE could mediate GSDMD cleavage and activation in a caspase-independent manner in neutrophils, and could in the future be used to promote macrophage pyroptosis. GPX4, an enzyme that can effectively reduce lipid peroxidation, regulates lipid peroxidation-dependent caspase-11 activation and GSDMD cleavage (Kang et al., 2018). TIRAP regulates the expression of caspase-11, a caspase-1-related protease, which is crucial for the activation of inflammatory bodies (Gurung et al., 2012). All in all, these 4 PRGs are included in pyroptosis and affect cancer cell progression, but whether it takes action in EC needs more exploration due to a few being reported in EC.
These four PRGs were used to construct the prediction panel. In univariate Cox analysis, there may be some false correlation or indirect correlation between the independent variable and the dependent variable. For example, factor A has no impact on the outcome event, while factor B is an influencing factor for the outcome event. However, since factor A is only simple and has a strong correlation with factor B, there is collinearity between them, so in univariate Cox analysis, there may be significant differences in factor A, which leads to factor a being mistaken as an influencing factor and included in the multi-factor analysis. In multivariate analysis, by adjusting the influence of factor B, the "false correlation" between factor A and the dependent variable disappears. At this time, it can be considered that factor A is not an influencing factor for the outcome event. Therefore, we conducted univariate, lasso, and multivariate Cox analyses to create the panel and avoid bias. KM plotter analysis, AUC of ROC curves, and univariate and multivariate cox analysis suggest that it can be used to predict the overall survival of EC patients independently. In addition, the risk score calculated by this panel in each patient was positively related to clinical characteristics, such as age, stage, and grade. The DEGs between two groups then conducted functional enrichment in GO, KEGG, GSEA, and immune characteristics analysis. In GO analysis, "phagocytosis, " "immunoglobin, " and "immune response" had the most frequent occurrences among the top 10 results, which suggested that the immune microenvironment may play a part in EC progression. In KEGG analysis, popular pathways, such as the Wnt signaling pathway, and cAMP signaling pathway, and some necessary substant metabolism, stressed the roles of PRGs in endometrial cancer. In GSEA analysis, pathways related to carcinogeneses such as ubiquitin-mediated proteolysis, basal transcription factors, cell cycle, regulation of autophagy, ERBB, TGFβ, Wnt pathways, and pathways in cancer were enriched. We next undertook immune function and immune cell distribution among the two groups. The results indicated that antigen-presenting cells, including DCs, pDCs, and T helper cells can help present pyroptosis cells to the T cells, and costimulate T cells and para secretion interferon to trigger a subsequent reaction. As reported, gasdermin-mediated cell pyroptosis was related to T cells, dendritic cells, neutrophils, and NK cells. Natural killer cells and cytolytic T cells release granzyme A or B mediating caspase activation and gasdermin family protein cleavage in macrophages, which increases pore-forming activity and formation of pore-forming protein, and thus cell pyroptosis (Wang et al., 2013;Jorgensen et al., 2016;Kambara et al., 2018;Xi et al., 2019;Liu et al., 2020;Zhang Z. et al., 2020;Zhou et al., 2020). However, how the pyroptosis-related immune system plays a role in EC needs more exploration and confirmation in vitro and in vivo.
Finally, we used an external database and low-throughput experiments to validate the four prognostic DEPRGs in EC. The results of the KM-plotter showed that low expression of these four prognostic DEPRGs had a longer survival time, which is consistent with the former analysis. Furthermore, the results of the GEPIA and RT-PCR experiment showed the expression level difference in EC and normal endometrial tissues.
In summary, we identified fpur prognosis-associated pyroptosis-related genes (ELANE, GPX4, GSDMD, and TIRAP). The model can also predict prognosis prediction and microenvironment in human endometrial cancer. However, to date, there have been few studies in the molecular mechanisms in pyroptosis related pathways and immunotherapy, and future large-scale clinical cohorts need to be undertaken to correct this panel.

DATA AVAILABILITY STATEMENT
The original contributions presented in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding author.

ETHICS STATEMENT
The studies involving human participants were reviewed and approved by This study was carried out in accordance with the standards of the Helsinki Declaration of the World Medical Association, and approved by the Ethics Committee of China Medical University. The patients/participants provided their written informed consent to participate in this study.

AUTHOR CONTRIBUTIONS
XZ downloaded the dataset, analyzed the data, and was a major contributor in writing the manuscript. QY reviewed the manuscript. Both authors read and approved the final manuscript.