Abstract

Background. Lung cancer is considered to be the second most aggressive and rapidly fatal cancer after breast cancer. Necroptosis, a novel discovered pattern of cell death, is mediated by Receptor-interacting serine/threonine-protein kinase 1 (RIPK1), Receptor-interacting serine/threonine-protein kinase 3 (RIPK3), and Mixed Lineage Kinase Domain Like Pseudokinase (MLKL). Methods. For the purpose of developing a prognostic model, Least absolute shrinkage and selection operator (LASSO) Cox regression analysis was conducted. Using Pearson’s correlation analysis, we evaluated the correlation between necroptosis-related markers and tumor immune infiltration. A bioinformatics analysis was conducted to construct a necroptosis-related regulatory axis. Results. There was a downregulation of most of necroptosis-related genes in lung adenocarcinoma (LUAD) versus lung tissues but an increase in PGAM5, HMGB1, TRAF2, EZH2 levels. We also summarized the Single Nucleotide Variant (SNV) and copy number variation (CNV) of necroptosis-related genes in LUAD. Consensus clustering identified two clusters in LUAD with distinct immune cell infiltration and ESTIMATEScore. Genes related to necroptosis were associated with necroptosis, Tumor necrosis factor (TNF) signaling pathway, and apoptosis according to Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways. Four prognostic genes (ALDH2, HMGB1, NDRG2, TLR2) were combined to develop a prognostic gene signature for LUAD patients, which was highly accurate in predicting prognosis. Univariate and multivariate analysis identified HMGB1, pT stage, and pN stage as independent factors impacting on LUAD patients’ prognosis. A significant correlation was found between the level of TLR2 and NDRG2 and clinical stage, immunity infiltration, and drug resistance. Additionally, the progression of LUAD might be regulated by lncRNA C5orf64/miR-582-5p/NDRG2/TLR2. Conclusion. The current bioinformatics analysis identified a necroptosis-related prognostic signature for LUAD and their relation to immunity infiltration. This result requires further investigation.

1. Introduction

Lung cancer will account for d 228,820 new cases and 135,720 deaths in 2020 all globally [1]. Adenocarcinoma of the lung (LUAD) accounts for more than 40% of all cases of lung cancer [2]. Due to the absence of typical clinical symptoms, patients usually have disseminated metastatic tumors when initially diagnosed with LUAD. Worse still, LUAD is characterized by high aggression and rapidly fatality with overall survival (OS) less than 3 years [3]. Even though smoking has been identified as a risk factor for LUAD, the molecular mechanism is still not understood. And exploration and identification of the potential tumorigenesis molecular mechanism and prognostic markers for LUAD are urgent.

Necroptosis, a novel discovered pattern of cell death, is mediated by Receptor-interacting serine/threonine-protein kinase 1 (RIPK1), Receptor-interacting serine/threonine-protein kinase 3 (RIPK3), and Mixed Lineage Kinase Domain Like Pseudokinase (MLKL) [4, 5]. Studies have implicated necroptosis in the pathogenesis of Parkinson, Alzheimer, vascular atherosclerosis, and infectious disease [57]. The death of T cells and cancer metastasis can also be accelerated by necroptosis, according to recent evidence [8]. On the other hand, when immunotherapy is used to treat malignancies, necroptosis may trigger and amplify antitumor immunity [4]. Necroptosis regulators may also provide prognostic information for cancer and other diseases [9, 10]. Moreover, some necroptosis-related prognostic signature had been identified for types of cancer, including stomach adenocarcinoma, breast cancer and cervical squamous cell carcinoma, and endocervical adenocarcinoma [1113]. As such, necroptosis-related genes may also contribute to LUAD prognosis.

Big data mining has been proposed as a promising tool for examining tumorigenesis mechanisms, associated prognosis markers, and therapy targets, following the development of the Cancer Genome Atlas (TCGA). Herein, we systematically investigated the expression of genes associated with necroptosis, their prognostic significance, and their association with immune infiltration. The data may offer another evidence about the vital functions of necroptosis in LUAD.

2. Materials and Methods

2.1. Datasets and Preprocessing

The RNA sequencing data (Fragments Per Kilobase of exon model per Million mapped fragments (FPKM) value) and copy number variation (CNV) data of LUAD patients were downloaded from TCGA database. After requiring the clinic characters of LUAD cohort from TCGA, we rearranged it and Supplementary Table 1 showed the detail. Using R (version 4.0.5) and R Bioconductor packages, dataset processing and further analysis were carried out. All the patients losing some clinical information were rejected and a total of 486 cases were obtained. In order to analyze the expression profile, we normalized it to transcripts per kilobase million values.

2.2. Expression, Genetic Mutation, and Functional Enrichment Analysis

By reviewing the previous literatures, we identified 17 necroptosis-related genes, including RIPK1, RIPK3, MLKL, TLR2, TLR3, TLR4, TNFRSF1A, PGAM5, ZBP1, NR2C2, HMGB1, CXCL1, USP22, TRAF2, ALDH2, EZH2, and NDRG2 [1422]. Using “limma” and “reshape2” package [23], the expression patterns of genes associated with necroptosis were generated. A mutation frequency of the gene was calculated using “maftools” package [24]. The “RCircos” package also allowed us to visualize the chromosome locations of CNVs associated with necroptosis [25]. We then performed functional enrichment analysis (Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG)) with “clusterProfiler “package [26].

2.3. Consensus Cluster Analysis

We then conducted consensus cluster analysis using the “ConsensusClusterPlus” package, which was calculated 1,000 times [27]. In the following step, survival and gene expression were analyzed using the packages “survival” and “pheatmap”. ESTIMATE algorithm was applied for evaluating the difference of Immunoscore, StromaScore, ESTIMATEScore, and immune cell infiltration in each cluster of LUAD. The immune cell landscape in each cluster of LUAD was generated with “vioplot” package.

2.4. Development of Necroptosis-Related Prognostic Gene Signature

The log-rank test was used to calculate the -values, and hazard ratio for the necroptosis-related prognostic gene found by Kaplan–Meier analysis. In the following step, the development of the prognostic model was counted on the Least absolute shrinkage and selection operator (LASSO) Cox regression analysis based on necroptosis-related prognostic genes. The risk score of each LUAD patient was calculated by a computational equation (sum of coefficients × necroptosis-related gene expression). TCGA-LUAD patients were separated into two subgroups (low- and high-risk), and the cut-off was the median value of the gene expression. We used Kaplan-Meier analysis to determine the OS curve, and time ROC analysis to determine how well this prognostic signature predicted outcomes. The prognostic risk factor was further analyzed using univariate and multivariate Cox analysis using clinical characteristics and gene signatures. A predicted nomogram was developed to evaluate the predictive performance in OS rate (1-, 3-, and 5-year).

2.5. Prognostic Gene Individual Analysis

Wilcox test was utilized to evaluated expression difference of prognostic signature in different pathological stage of LUAD patients. A correlation was then determined between prognostic signature and immune cell infiltration using Tumor Immune Estimation Resource (TIMER) (https://cistrome.shinyapps.io/timer/). This was followed by tumor mutation burden (TMB) and microsatellite instability (MSI) analysis using Spearman’s correlation method. We then collected the IC50 of 265 small molecules in 860 cell lines, and its corresponding necroptosis-related prognostic gene expression from Genomics of Drug Sensitivity in Cancer (GDSC). Using Pearson correlation analysis, the relation between gene expression and drug IC50 was determined.

2.6. lncRNA–miRNA–MRNA Network Construction

The miRNA targets of necroptosis-related prognostic gene were predicted with miRDB (http://mirdb.org/) [28]. Moreover, we evaluated the expression and prognosis significance of these miRNAs in LUAD to further screen the most promising targets. This was followed by lncRNA targets predicting using LncBase (https://carolina.imis.athena-innovation.gr/) [29] and RNAInter (http://www.rna-society.org/) [30] with a coefficient>0.7. A similar approach to identifying promising lncRNA targets was used to investigate the expression and prognosis significance of these lncRNAs.

3. Results

3.1. Expression and Mutation Landscape of Necroptosis-Related Genes in LUAD

Supplementary Figure S1 revealed the work flow of the current study. Figure 1(a) displays the expression landscape, revealing that downregulation of most of necroptosis-related genes were shown in LUAD versus lung tissues while the level of PGAM5, HMGB1, TRAF2, EZH2 was raised (all ). Figures 1(b) and 1(c) showed the genetic landscape of necroptosis-related gene in LUAD. To be more specific, 64.59% (106/138) of LUAD cases presented with genetic mutations (Figure 1(b)). In terms of mutation rate, TLR4 ranked first, followed by ZBP1 (Figure 1(b)). Missense mutations ranked highest, and C > A was the most common Single Nucleotide Variant (SNV) type (Figure 1(c)). For CNV analysis, the result indicated widespread copy number amplification of CXCL1, TNFRSF1A, NDRG2, RIPK1, USP22, TRAF2, TLR4, and MLKL, as well as widespread CNV deletion of RIPK3, PGAM5, EZH2, HMGB1, TLR3, ZBP1, NR2C2, ALDH2, and TLR2 (Figure 1(d)). The location of CNV alteration of necroptosis-related genes on chromosomes was showed in Figure 1(e).

3.2. GO and KEGG Analysis

We then further confirmed whether these genes were associated with necroptosis in LUAD by performing GO and KEGG pathways analysis. Accordingly, these genes were widely associated with programmed necrotic cell death, necrotic cell death, necroptotic process and NF-kappa B signaling, membrane raft, CD40 receptor complex, transcription coregulator activity, cytokine binding, and Tumor necrosis factor (TNF) receptor superfamily binding in GO analysis (Supplementary Figure 2(a)). Furthermore, these necroptosis-related genes showed widespread association with necroptosis, TNF signaling pathway, NF-kappa B signaling pathway, Nucleotide oligomerization domain (NOD)-like receptor signaling pathway, and apoptosis in KEGG pathway analysis (Supplementary Figure 2(b)).

3.3. Identification of Two Clusters Based on Necroptosis-Related Genes in LUAD

To clarify whether these LUAD patients can be divided into multiple subgroups to achieve precise treatment for patients, we then performed consensus clustering analysis. LUAD patients were differentiated using consensus clustering based on gene patterns. And two clusters (cluster 1/2) were suggested as the optimal clustering stability based on the similarity (Figure 2(a)). No significant difference was generated in OS rate between these two clusters in LUAD (Figure 2(b), ). Interestingly, cluster 1 showed distinctly different age and pM stage with cluster 2 (Figure 2(c), ). In Figure 3(a), the immune cell infiltration landscape was shown for two clusters of LUAD, with cluster 1 correlated with high abundance of plasma cells (), CD4 memory resting T cells (), Tresgs (), and neutrophils () while cluster 1 was correlated with high abundance of CD8 T cells () and follicular helper T cells () in LUAD. Moreover, the data suggested a higher Immunoscore (), StromaScore (), and ESTIMATEScore () in cluster 1 versus cluster 2 (Figures 3(b), 3(c), and 3(d)).

3.4. A Prognostic Signature Based on Necroptosis-Related Genes

Prognosis analysis demonstrated that high level of ALDH2, NDRG2, TLR2, TLR4, and low HMGB1 level had a better OS rate in LUAD (Supplementary Figures 3(a), 3(b), 3(c), 3(d), and 3(e)) and Supplementary Table 2). Based on these five necroptosis-related prognostic genes, we then developed a prognostic gene signature with LASSO Cox regression analysis. The coefficients of each LUAD case was calculated with the followed computational equation: risk score = sum of coefficients × the expression of necroptosis-related genes. Finally, based on the result of LASSO Cox regression analysis, the best module was obtained. TLR4 was ejected, and the risk score of patients was calculated by including four other genes in this prognostic signature (Riskscore = (−0.1017) × ALDH2 expression + (0.1559) × HMGB1 expression + (−0.0698) × NDRG2 expression + (−0.0845) × TLR2 expression). The coefficient and partial likelihood deviance of prognostic signature were shown in Figures 4(a) and 4(b). LUAD cohort was divided into high- and low-risk group, and the riskscore, survival status of patients, and gene expression were shown in Figure 4(c). Compared with low-risk group, high-risk group had a worse OS rate () with the area under the curves in 1-year, 3-year, and 5-year periods being 0.684, 0.592, and 0.584, respectively (Figures 4(d) and 4(e)).

3.5. Predictive Nomogram Based on Prognostic Signature and Clinical Characters

Considering clinical characters and above four necroptosis-related prognostic genes, we identified HMGB1, pT stage, and pN stage as independent factors impacting on LUAD patients’ prognosis in further analysis (univariate and multivariate analysis) (Figures 5(a) and 5(b)). These factors were used to construct a predictive nomogram to predict survival probability, which showed a good prediction ability (Figures 5(c) and 5(d)).

3.6. Individual Analysis of Necroptosis-Related Prognostic Signature

As shown in Figure 6(a), a noteworthy correlation was obtained between pathological stage and TLR2 expression () and NDRG2 expression (), suggesting that TLR2 and NDRG2 may be correlated with the progression of LUAD. And we select TLR2 and NDRG2 for further study. Previous study had suggested that immune infiltration was involved in tumor development and progression in LUAD [31]. The current result demonstrated a positive correlation between the expression of TLR2 and NDRG2 and the immune infiltration level of B cells, CD8+ T cells, CD4+ T cells, macrophage, neutrophils, and dendritic cells (Figure 6(b), all ). Moreover, some somatic copy number alterations of TLR2 and NDRG2 could inhibit immune cell infiltration level (Figure 6(c)). TMB and MSI were suggested as predictive marker for cancer immunotherapy, including lung cancer [3234]. Interestingly, the expression of TLR2 and NDRG2 showed significant correlation with TMB score () and MSI score () (Figure 6(c)). In MSI analysis, MSI score decreased as NDRG2 expression increased (, Figure 6(d)). To identify cancer immunotherapy target, a vital way is to evaluate the correlation between gene expression and existed drug targets. In order to clarify whether TLR2 and NDRG2 could serve as potential drug screening targets, we explored the correlation between TLR2 and NDRG2 and existed drug targets in GDSC database. Interestingly, the expression of TLR2 and NDRG2 showed positive or negative correlation with GDSC drug sensitivity, including methotrexate and vinblastine (Figure 6(e)).

3.7. Construction of Necroptosis-Related Regulatory Axis

We then constructed a necroptosis-related regulatory axis to clarify the potential molecular mechanism of TLR2 and NDRG2 in LUAD. Using miRDB, we identified miR-582-5p and miR-5699-5p as the miRNA targets of NDRG2 and TLR2 (Figure 7(a)). Moreover, the expression of miR-582-5p () and miR-5699-5p () were upregulated in LUAD versus lung tissues (Figures 7(b) and 7(c)). OS analysis suggested a poor survival in LUAD patients with high miR-582-5p expression (Figure 7(d), ). Accordingly, miR-582-5p may be the most potential miRNA target of NDRG2 and TLR2. To explore its upstream lncRNA targets, we submitted miR-582-5p to RNAInter and lncBase, and the result suggested three lncRNA targets (MALAT1, C5orf64, and SNHG16) interacting to miR-582-5p (Figure 7(e)). Further analyses indicated downregulation of MALAT1, C5orf64, and upregulation of SNHG16 in LUAD versus lung tissues (Figure 7(f), all ). However, only lncRNA C5orf64 was correlated with OS rate in LUAD (Figure 7(g), ), suggested C5orf64 as the most promising lncRNA target. Therefore, we identified lncRNA C5orf64/miR-582-5p/NDRG2/TLR2 regulatory axis in the progression in LUAD. This result requires further investigation.

4. Discussion

Previous study had suggested the involvement of necroptosis migration and invasion regulation of tumor [35]. The necroptosis mechanism was proposed as an effective way for eradicating cancer cells [36]. The identification of prognostic value and potential regulatory axis of necroptosis-related genes will allow necroptosis to be leveraged for therapeutic benefits and prognosis improvement of LUAD.

To confirm whether these necroptosis-related genes were associated with necroptosis in LUAD, we then performed GO and KEGG pathways analysis. As expected, these necroptosis-related genes showed widespread association with necroptosis, programmed necrotic cell death, necroptotic process, and TNF signaling pathway. These functions or pathways could mediate necroptosis and cancer progression. NF-κB signaling may influenced inflammation and the progression of tumor [37]. Signaling mediated by TNF is also crucial for homeostasis and immunity in mammals [38]. Interestingly, TNF was referred as a key mediator in balancing cell survival and necroptosis [39].

LUAD patients were differentiated using consensus clustering based on gene patterns and we identified two clusters, which showed conspicuous difference in pM stage and immune cell characterization. Cluster 1 of LUAD was linked to high Immunoscore, StromaScore, and ESTIMATEScore and abundant immune cell infiltration, referring to hot tumor [40]. Further moreover, high Immunoscore was significantly correlated with better prognosis in LUAD [40].

Our study also developed a prognostic model based on four prognostic necroptosis-related genes (ALDH2, HMGB1, NDRG2, TLR2), which was highly accurate in predicting LUAD prognosis. Univariate and multivariate analysis identified HMGB1, pT stage, and pN stage as independent factors impacting on LUAD patients’ prognosis. Though certain prognostic signatures had been identified for LUAD, our study firstly developed a prognostic model with necroptosis-related markers in LUAD, providing another biomarker in LUAD. A machine learning strategy had constructed and validated a prognostic signature using 12 immune-related genes for LUAD [41]. Another prognostic signature showed good performance in predicting prognosis and reflecting tumor immune microenvironment in LUAD [42]. Jin et al. also constructed a 7-lncRNA prognosis signature and a predictive nomogram in LUAD [43].

The result of our study identified a lncRNA C5orf64/miR-582-5p/NDRG2/TLR2 regulatory axis in the progression in LUAD. Interestingly, previous study demonstrated lncRNA C5orf64 as a novel biomarker associated with tumor microenvironment and mutation pattern remodeling in LUAD [44]. Moreover, miR-582-5p served as a prognostic biomarker in LUAD and inhibited tumor cell proliferation and invasion [45]. High mRNA level of TLR2 could accelerate tumor progression in LUAD [46]. NDRG2 acted a prognostic marker in LUAD and associated with depth of invasion, vascular invasion, and better OS [47]. Thus, lncRNA C5orf64/miR-582-5p/NDRG2/TLR2 regulatory axis may be involved in the progression in LUAD. This result requires further investigation.

Our study also had some limitations. Not all of 17 necroptosis-related genes were specific to necroptosis. Moreover, the result of consensus clustering analysis is barely satisfactory. This result requires further investigation.

5. Conclusion

The current bioinformatics analysis identified a necroptosis-related prognostic signature for LUAD and their relation to immunity infiltration. This result requires further investigation.

Data Availability

The analyzed data sets generated during the study are available from the corresponding author on reasonable request.

Conflicts of Interest

The authors declare that they have no competing interests.

Authors’ Contributions

Zhiwu Han designed the study and wrote the manuscript. Libo Sun performed data analysis work. Wenwen Li, Zhenhuan Zhao, and Yanhua Zuo edited the manuscript. All authors read and approved the final manuscript.

Acknowledgments

This study was funded by Libo Sun and Zhiwu Han.

Supplementary Materials

Supplementary Materials. Supplementary Table 1. The clinical characters of lung adenocarcinoma patients in TCGA cohort. Supplementary Table 2. The result of necroptosis-related genes in prognosis analysis.

Supplementary Materials. Supplementary Figure 1. The work flow of the current study.

Supplementary Materials. Supplementary Figure 2. The enriched item in gene ontology (a) and Kyoto Encyclopedia of Genes and Genomes analysis (b). BP, biological process; CC, cellular component; MF molecular function.

Supplementary Materials. Supplementary Figure 3. The prognostic value of necroptosis-related genes in LUAD. Overall survival curve in LUAD patients with high/low expression of ALDH2 (a), NDRG2 (b), TLR2 (c), TLR4 (d), and HMGB1 (e).