Construction of a Ferroptosis-Related Long Non-coding RNA Prognostic Signature and Competing Endogenous RNA Network in Lung Adenocarcinoma

Ferroptosis-related genes play an important role in the progression of lung adenocarcinoma (LUAD). However, the potential function of ferroptosis-related lncRNAs in LUAD has not been fully elucidated. Thus, to explore the potential role of ferroptosis-related lncRNAs in LUAD, the transcriptome RNA-seq data and corresponding clinical data of LUAD were downloaded from the TCGA dataset. Pearson correlation was used to mine ferroptosis-related lncRNAs. Differential expression and univariate Cox analysis were performed to screen prognosis related lncRNAs. A ferroptosis-related lncRNA prognostic signature (FLPS), which included six ferroptosis-related lncRNAs, was constructed by the least absolute shrinkage and selection operator (LASSO) Cox regression. Patients were divided into a high risk-score group and low risk-score group by the median risk score. Receiver operating characteristic (ROC) curves, principal component analysis (PCA), and univariate and multivariate Cox regression were performed to confirm the validity of FLPS. Enrichment analysis showed that the biological processes, pathways and markers associated with malignant tumors were more common in high-risk subgroups. There were significant differences in immune microenvironment and immune cells between high- and low-risk groups. Then, a nomogram was constructed. We further investigated the relationship between six ferroptosis-related lncRNAs and tumor microenvironment and tumor stemness. A competing endogenous RNA (ceRNA) network was established based on the six ferroptosis-related lncRNAs. Finally, we detected the expression levels of ferroptosis-related lncRNAs in clinical samples through quantitative real-time polymerase chain reaction assay (qRT-PCR). In conclusion, we identified the prognostic ferroptosis-related lncRNAs in LUAD and constructed a prognostic signature which provided a new strategy for the evaluation and prediction of prognosis in LUAD.

Introduction

Lung cancer is the disease with the highest morbidity and mortality, among which lung adenocarcinoma (LUAD) is accounted for 40–50% of all lung cancer cases (Bray et al., 2018). At present, local resection is still the first choice for the treatment of LUAD. In recent years, molecular targeted therapy and immunotherapy have developed rapidly, which has brought success in the clinical treatment of a series of malignant tumors. Although the prognosis of some patients has improved significantly with the application of targeted therapy and immunotherapy, the survival rate is far from satisfaction (Miller et al., 2019). There is evidence showing that identification and application of novel biomarkers could bring benefits for the effective treatment of patients (Xu et al., 2021).

Ferroptosis is a form of regulatory cell death with iron dependence, caused by excessive lipid peroxidation and related to the occurrence of a variety of tumors and therapeutic response (Chen et al., 2021). Ferroptosis also plays an important role in LUAD. CAMP responsive element binding protein 1 (CREB) can directly bind to the promoter region of GPX4 to promote its expression, thereby inhibiting potential ferroptosis and promote the growth of LUAD (Wang Z. et al., 2021). P53 is a well-known tumor suppressor, which mediates apoptosis and cell-cycle arrest. In addition, it also can modulate ferroptosis in an ALOX12-dependent way and inhibit the proliferation of H1299 (Chu et al., 2019).

The abnormal expression of lncRNAs in LUAD is widely involved in the process of tumor proliferation and metastasis (Li and Chen, 2016). PD-L1-lnc was reported to directly bind to c-Myc and enhance its transcriptional activity, ultimately promoting proliferation and invasion of LUAD (Qu et al., 2021). Knockdown of lncRNA MSC-AS1 could inhibit LUAD cell growth and accelerate apoptosis through miR-33b-5p/GPAM axis (Li S. et al., 2021). What’s more, accumulating studies demonstrated that lncRNAs also regulated ferroptosis of lung cancer by interacting with ELAVL1 or miR-365a-3p (Wang et al., 2019; Gai et al., 2020). Previous studies have shown that ferroptosis-related lncRNAs could be prognostic risk factors and the signatures constructed with them can be used to distinguish high-risk patients from low-risk patients (Cai et al., 2021; Tang et al., 2021). However, the full role of ferroptosis-related lncRNAs in LUAD needs to be further explored.

Therefore, we retrieved ferroptosis-related genes from previous studies and mined potential dysregulated ferroptosis-related lncRNAs through bioinformatic analysis, based on the LUAD dataset of TCGA. Our study showed 33 lncRNAs were differentially expressed in LUAD and correlated with prognosis. Furthermore, a ferroptosis-related lncRNA prognostic signature (FLPS) was constructed based on LASSO regression. We confirmed the validity of FLPS through ROC, PCA, and Cox regression. Finally, a nomogram was constructed to predict overall survival (OS) of LUAD patients and a ceRNA network was built to further clarify the function of those lncRNAs in LUAD.

Materials and Methods

Data Collection

The RNA-seq transcriptome data, including 497 LUAD samples and 54 adjacent normal tissues, and corresponding clinical data were extracted from TCGA database¹ for differential expression analysis. All patients included in the prognostic analysis fit the following criteria: (1) histologically confirmed LUAD and (2) available information on gene expression and survival. Lastly, 468 patients with corresponding clinicopathological information were enrolled for further study. A total of 468 patients with LUAD were randomly classified into a training cohort (235 patients) and a test cohort (233 patients) at a 1:1 ratio by using the “caret” package, according to previous studies (Yang et al., 2019; Geng et al., 2021). Then 283 ferroptosis-related genes were retrieved from the previous study, including 132 drivers and 151 suppressors (Liu et al., 2020). TCGA pan-cancer data, including RNA-seq and clinical data and stemness scores based on mRNA (RNAss) and DNA-methylation (DNAss) of LUAD were extracted from xena browser². TCGA pan-cancer data include 33 cancer types. Among them, 18 cancer types had more than five associated normal tissue samples and were used to be further investigated.

Establishment and Validation of the Ferroptosis-Related Long Non-coding RNA Prognostic Signature

The annotation of lncRNA was performed as a previous study (Tu et al., 2020). Pearson correlation was used to mine ferroptosis-related lncRNAs (with the | Person R| > 0.4 and p < 0.001).

Differentially expressed genes (DEGs), miRNAs (DEMs) and lncRNAs between tumor and adjacent normal tissues were identified by using the “limma” or “edgeR” R package (| log2FC| > 1, FDR < 0.05). Univariate Cox analysis of overall survival (OS) was performed to screen ferroptosis-related genes and lncRNAs with prognostic value. Then, R package “glmnet” was used to conduct least absolute shrinkage and selection operator (LASSO) Cox regression. A ferroptosis-related lncRNA prognostic signature (FLPS) was constructed for the LUAD patients, which involved six ferroptosis-related lncRNAs. The risk score was calculated as the formula:

$$\mathrm{Risk}\mathrm{score}=\sum _{\mathrm{i}=1}^{\mathrm{n}}\mathrm{\beta }\mathrm{i}\mathrm{\chi }\mathrm{i}$$

where βimeans the coefficients, whereas χi is the FPKM value of each ferroptosis-related lncRNAs.

Risk scores were calculated for all patients in our study. Kaplan–Meier curve was used to compare the OS between the high-risk and low-risk group through R package “survival” and “survminer.” Time-dependent ROC analyses and the area under the curve (AUC) were performed to assess the model by using “timeROC.” Principal component analysis (PCA) and scatter diagrams were performed by R package “ggplot2.”

Functional Enrichment and Immune-Related Scores Analysis

Differentially expressed genes between the high-risk subgroup and low-risk group were identified based on the standards of | log2FC| > 1 and p < 0.05 using the R package “limma.” Then the DEGs were inputted into the “Metascape” website for functional and pathway enrichment analysis (Zhou et al., 2019). In addition, GSEA software was used to explore the hallmarks which were highly enriched in the high-risk group. The infiltrating score of 16 immune cells and the activity of 13 immune-related pathways were calculated with single-sample gene set enrichment analysis (ssGSEA) in the “gsva” R package. The immune, stromal and estimate score for each patient was calculated by the R “estimate” package.

Competing Endogenous RNA Network Construction

The online tool LncBase V2.0³ (Paraskevopoulou et al., 2016) was used to predict the downstream target miRNA with a 0.6 threshold (Paraskevopoulou et al., 2016). TargetScan⁴ (Agarwal et al., 2015), miRDB⁵ (Chen and Wang, 2020), and miRWalk⁶ (Dweep and Gretz, 2015) were used to predict the target genes. Then, the predicted target genes and miRNAs were intersected with the DEGs and DEMs, respectively. MiRNAs or mRNAs that were not consistent with the expression of lncRNAs or miRNAs were eliminated. Then the network of lncRNAs, miRNAs, and mRNAs was constructed by Cytoscape (version 3.8.2) (Shannon et al., 2003).

Screening of Hub Genes and Survival-Related Genes

The STRING database⁷ (Szklarczyk et al., 2021) was applied to construct PPI network and discover the relationship among the target genes with high confidence (0.7). Then, the gene relationship was imported into Cytoscape and the top 10 hub genes were calculated by cytoHubba. Kaplan–Meier curve and log-rank test were used to mine prognostic genes.

Samples and Quantitative Real-Time Polymerase Chain Reaction

A total of 20 paired LUAD tissues and corresponding adjacent non-tumorous tissues were obtained from patients who underwent radical resection of lung cancer in Changhai Hospital, from June 2020 to September 2020. This study was approved by Ethics Committee of Changhai Hospital (Shanghai, China), and informed consent was signed before the operation for each patient. Total RNA was extracted from samples with TRIzol^TM Reagent (Invitrogen). RNA was reverse-transcribed using PrimeScript^TM RT Master Mix (Takara). Real-time PCR was performed with TB Green^TM Premix EX Taq^TM II (Takara). The expression levels of lncRNAs were normalized by GAPDH. The sequences of primers were as follows: AC021016.1 forward 5′-GGGTCAAGCACACTGAGGGT-3′, and reverse 5′-ACCAGGTGTGAACCCTTGGG-3′; KTN1-AS1 forward 5′- AGGGAAATTTGGGCAGAAGT-3′, and reverse 5′-GTTACCC GTGTGAGCCTGAT-3′; GAPDH forward 5′GTCTCCTCT GACTTCAACAGCG-3′, and reverse 5′-ACCACCCTGTTGCT GTAGCCAA-3′.

Statistical Analysis

All statistical analyses were performed using R language 4.0.4 version and attached packages. Differentially expressed genes (DEGs), miRNAs (DEMs), and lncRNAs between tumor and adjacent normal tissues were identified by Wilcox test. Kaplan–Meier curve and log-rank test were used to compare OS between subgroups. Time-dependent ROC analyses and the area under the curve (AUC) were used to assess the performance of the model. Independent risk factor was screened by univariate and multivariate COX regression analysis. Spearman correlation was used to test the correlation between gene expression and stemness score, stromal score, immune score, and estimate score. Comparison of gene expression between the normal and tumors were performed in 18 cancer types which had more than five associated adjacent normal samples using linear mixed effects models. Mann–Whitney U test was used to compare the ssGSEA scores of immune cells or pathways between the high-risk and low-risk group.

Results

Identification of Ferroptosis-Related Long Non-coding RNA in Lung Adenocarcinoma Patients

The workflow for construction of risk model and subsequent analyses is shown in Supplementary Figure 1. Firstly, a total of 283 ferroptosis-related gene expression matrices were extracted from the TCGA dataset. Among 283 genes, 18 genes were differentially expressed between tumor and adjacent normal tissues (Figure 1A) and consistent with univariate Cox regression analysis results (Figure 1B). A lncRNA whose expression value was correlated with one or more of the 18 ferroptosis-related genes was defined as a ferroptosis-related lncRNA. Through Pearson correlation analysis (| Pearson R| > 0.4 and p < 0.001), 590 ferroptosis-related lncRNAs were uncovered. Through analysis of differential expression and univariate Cox regression, we finally obtained 33 lncRNAs which were highly expressed in tumor and predicted a worse prognosis or were lowly expressed and predicted a better prognosis (Figures 1C,D).

FIGURE 1

Figure 1. Identification of ferroptosis-related lncRNAs in LUAD patients and construction of a prognostic signature. (A) Heatmap for 18 differentially expressed genes, with red indicating high expression and blue indicating low expression. (B) Forest map for 18 prognostic genes. (C) Heatmap for 33 differentially expressed lncRNAs, with red indicating high expression and blue indicating low expression. (D) Forest map for 33 prognostic lncRNAs. (E) The optimal log λ value is indicated by the vertical black line in the plot. (F) The LASSO coefficient profile of ferroptosis-related lncRNAs, each line represents an independent lncRNA.

Construction and Validation of the Prognostic Signature for Ferroptosis-Related Long Non-coding RNA

To build the ferroptosis-LPS for predicting the OS of LUAD patients, a total of 468 patients with LUAD were randomly classified into a training cohort (235 patients) and a test cohort (233 patients) at a 1:1 ratio firstly. Then we performed the least absolute shrinkage and selection operator (LASSO) regression analysis based on the expression values of 33 ferroptosis-related lncRNA in the training cohort (Figures 1E,F). Six lncRNAs, namely, AC021016.1, AC068228.2, MIR223HG, AC009275.1, AL049555.1, and KTN1-AS1, were identified. The risk scores of each patient in TCGA training and validation cohorts were calculated based on the coefficient for each lncRNA, and the formula is as follows: risk score = -0.1889 × AC021016.1 expression level) + (0.2885 × AC068228.2 expression level) - (0.2782 × MIR223HG expression level) + (0.1150 × AC009275.1 expression level) + (0.1759 × AL049555.1 expression level) + (0.1607 × KTN1-AS1 expression level). Patients in the TCGA training and test cohorts were divided into high-risk-score and low-risk-score subgroups based on the median value of risk scores. Risk score and survival status distributions are plotted in Figures 2A,F. The heatmap results suggested that AC021016.1 and MIR223HG were downregulated in the high-risk group, whereas the expression of AC068228.2, AC009275.1, AL049555.1, and KTN1-AS1 were upregulated in the high-risk group. Kaplan–Meier survival curves indicated that LUAD patients with low-risk scores had better clinical outcomes in either training or validation cohort (Figures 2B,G), and the ROC curves showed that ferroptosis-LPS had a promising ability to predict OS in the training and validation cohorts (Figures 2C,H). PCA analysis indicated the patients in two risk groups were distributed in two directions (Figures 2D,E,I,J).

FIGURE 2

Figure 2. Validation of the prognostic signature. (A) Distribution of ferroptosis-related lncRNAs model based on risk score for the training set, patterns of the survival time, and survival status between the high- and low-risk groups for the training set and clustering analysis heatmap shows the display levels of the six ferroptosis-related lncRNAs for each patient in the training set. (B) Kaplan–Meier survival curves of the OS of patients in the high- and low-risk cohorts for the training set. (C) Time-dependent ROC analysis of accuracy of the model in the training set. (D) 2D and (E) 3D plots of the PCA of the training set. (F) Distribution of ferroptosis-related lncRNAs model based on risk score for the testing set, patterns of the survival time, and survival status between the high- and low-risk groups for the testing set and clustering analysis heatmap shows the display levels of the six ferroptosis-related lncRNAs for each patient in the testing set. (G) Kaplan–Meier survival curves of the OS of patients in the high- and low-risk cohorts for the testing set. (H) Time-dependent ROC analysis of accuracy of the model in testing set. (I) 2D and (J) 3D plots of the PCA of the testing set.

Independence of the Ferroptosis-LPS Considering Other Clinical Factors

Univariate and multivariate Cox regression were performed to assess whether ferroptosis-LPS was an independent prognostic factor for patients with LUAD. In the TCGA training dataset, univariate Cox analysis indicated that ferroptosis-LPS and stage were remarkably associated with OS (p < 0.001, Figure 3A) and multivariate Cox analysis further confirmed that ferroptosis-LPS was a prognostic risk factor (Figure 3B). The same result was found in the testing cohort, which confirmed that ferroptosis-LPS was an independent risk factor for LUAD patients (univariate: p < 0.001; multivariate: p < 0.05; Figures 3C,D). These results suggested that our ferroptosis-LPS might be useful for clinical prognosis evaluation.

FIGURE 3

Figure 3. Ferroptosis-LPS was an independent prognostic factor for LUAD patients. (A,B) Univariate and multivariate Cox regression analyses in the training set. (C,D) Univariate and multivariate Cox regression analyses in the testing set.

Functional Analyses in the TCGA Training and Testing Cohort

To investigate the biological functions and pathways which were associated with risk score in the training and test cohorts, the DEGs between the low- and high-risk groups were used to performed functional and pathway enrichment analysis. As expected, these DEGs were significantly enriched in cell cycle, retinoblastoma gene in cancer, and meiotic cell cycle (p < 0.05, Figures 4A,C). Gene set enrichment analysis showed that the genes in high-risk group of both train and test cohorts were significantly enriched in several hallmarks, such as MTORC1 signaling, MYC targets, G2M checkpoint, E2F targets, and so on (Figures 4B,D).

FIGURE 4

Figure 4. Enrichment analysis of differentially expressed genes (DEGs) between the high- and low-risk subgroups. (A,C) Bar graph of enriched terms across input gene lists, colored by p-values, and network of enriched terms, colored by cluster ID, where nodes that share the same cluster ID are typically close to each other. (B,D) Gene set enrichment analysis (GSEA) indicating that hallmarks were enriched in the high-risk subgroup of the training and testing sets. “Log10(P)” is the p-value in log base 10.

To further explore the relationship between the risk score and immune status, ssGSEA was performed to calculate the infiltrating score of immune cells and immune-related pathways. The scores of CCR, Check-point, HLA, T_cell_co-stimulation, and Type_II_IFN_Response were lower in the high-risk group, which were confirmed by the results of testing cohort. Moreover, the infiltrating scores of aDCs, B cells, DCs, iDCs, pDCs, neutrophils, T helper cells, and TIL were obviously higher in the low-risk groups of training and testing cohorts (Figures 5A,C). We further investigated the relationship between risk score with tumor microenvironment. Interestingly, the immune, stromal, and estimate score were higher in the low risk-score group (Figures 5B,D).

FIGURE 5

Figure 5. Comparison of immune status and tumor microenvironment between the high- and low-risk subgroups. (A,C) The scores of 16 immune cells and 13 immune-related functions of the training and testing set are displayed in boxplots. (B,D) The relationship between risk score with tumor microenvironment of the training and testing set. (cluster1: low-risk-score group; cluster2: high-risk-score group).

Stratification Analysis of the Ferroptosis-LPS and Construction of the Ferroptosis-LPS-Based Nomogram

The clinicopathological features and risk score of each patient in the TCGA dataset are shown in the heatmap (Figure 6A). For better assessment of the prognostic ability of the FLPS, we performed a stratification analysis to confirm its ability to predict OS in various subgroups. Compared to low-risk groups, high-risk LUAD patients had worse OS in the male and female subgroup. The same results were found in patients with age ≤ 65 or > 65, T1-2 or T3-4, and stage I-II, N0 and M0 (Figure 6B).

FIGURE 6

Figure 6. Stratification analysis of the ferroptosis-LPS and construction of the ferroptosis-LPS-based nomogram. (A) Heatmap of clinicopathological features and risk score of each patient in the entire TCGA dataset. (B) Survival analysis in various subgroups. (C) The nomogram predicts the probability of the 1-, 3-, and 5-year OS.

To create a clinically applicable quantitative tool to predict the OS of LUAD patients, we established a nomogram using the risk score, age, and stage in the TCGA dataset (Figure 6C).

Association of Six Ferroptosis-Related Long Non-coding RNAs With Tumor Microenvironment and Tumor Stemness

As the risk score was significantly correlated with tumor microenvironment, we investigated the association between the expression levels of the six ferroptosis-related lncRNAs and the tumor microenvironment. We found that AC021016.1 and MIR223HG had a positive correlation with stromal scores, immune score and estimate score, whereas AC068228.2 had a negative correlation (p < 0.001, Figure 7A).

FIGURE 7

Figure 7. Association of expression of six lncRNAs with tumor microenvironment and tumor stemness and the pan-cancer expression analysis. (A) Correlation matrices between six lncRNAs expression and RNAss, DNAss, stromal score, immune score, and estimate score. (B) The expression level of six lncRNAs in tumors compared with normal tissue in 18 cancer types which were composed of more than five normal samples.

In the process of cancer progression, tumor cells can gradually lose the phenotype of differentiation and acquire the characteristics of progenitor cells and stem cells, which might increase tumor resistance. Since those six lncRNAs were prognostic risk factors, we explored the correlation between them and RNA stemness score (RNAss) and DNA stemness score (DNAss), which were indicators of tumor stemness. The results showed that AC021016.1 and MIR223HG were negatively correlated with tumor stemness in both RNAss and DNAss while AC068228.2 and KTN1-AS1 had the positive correlation (p < 0.001, Figure 7A).

Analysis of Six Ferroptosis-Related Long Non-coding RNAs in Pan-Cancer

To further understand the roles of six ferroptosis-related lncRNAs in tumors, we examined their expression levels in 18 cancer types in TCGA pan-cancer data. A striking inter-tumor heterogeneity of expression of six ferroptosis-related lncRNAs was observed. For instance, the expression levels of AC021016.1 showed the largest inter-tumor heterogeneity with some tumors having high levels of AC021016.1 (CHOL, KIRP, LIHC), while others were characterized with low levels of AC021016.1 expression (BLCA, BRCA, COAD, ESCA, GBM, HNSC, KICH, KIRC, LUAD, LUSC, PRAD, READ, STAD, UCEC). In addition, except MIR223HG, other ferroptosis-related lncRNAs had high expression levels in most tumor groups (Figure 7B).

Then we further explored the prognostic value of lncRNAs in tumors. Median was used to distinguish between high and low expression groups, whereas KM curves were performed to compare outcomes between the two groups. The results showed that the altered expression of ferroptosis-related lncRNAs were generally associated with patients’ overall survival. In most tumors, high expression of AC009275.1 and AC068228.2 were associated with poor prognosis, while other lncRNAs’ prognostic effect depended on the type of tumors (Supplementary Figure 2).

Construction of the Competing Endogenous RNA Network and Functional Enrichment Analysis

Several lncRNAs were found to act as miRNA sponges to participate in the regulation of gene expression. Then we explored the DEMs and DEGs based on the TCGA-LUAD dataset (|log2FC| > 1, FDR < 0.05); 369 DEMs and 1,517 DEGs were identified by the “edgeR” package. Top 20 up- and downregulated DEMs and DEGs are shown in Supplementary Figure 3. Using the LncBase Predicted v.2, we predicted the targeted miRNA of the six ferroptosis-related lncRNAs with a 0.6 threshold. Then the intersection of target miRNAs for downregulated lncRNAs with upregulated miRNAs of TCGA dataset and target miRNAs for upregulated lncRNAs with downregulated miRNAs of LUAD were taken. To further explore the target genes of the miRNAs, we applied TargetScan, miRDB, and miRWalk to predict the target genes. We used the negative correlation criteria [(a) miRNAs should be targeted to the mRNA; (b) the expression level of mRNAs should be opposite to miRNAs; and (c) the mRNAs belong to DEGs] to screen the downstream mRNAs. Ultimately, five ferroptosis-related lncRNAs, 45 miRNAs, and 107 mRNAs were included in our ceRNA network (Figure 8A).

FIGURE 8

Figure 8. Construction of the ceRNA network and screening of hub genes. (A) The ceRNA Network which consists of five ferroptosis-related lncRNAs, 45 miRNAs and 107 mRNAs. (B–D) Bubble plots showing GO analysis of 107 genes for biological process, cellular component and molecular function. (E) Cnetplot for KEGG signal pathway of 107 genes. (F,G) The PPI network of 41 genes and hub genes among them. (H,I) Kaplan–Meier survival curves of BIRC5 and CKS1B.

Furthermore, Go and KEGG enrichment analyses were performed to clarify the biological function of 107 downstream genes. Through BP analysis, response to peptide hormone, cell growth, regulation of cell growth, and regulation of cell-cell adhesion were obviously enriched (Figure 8B). The results of CC analysis contained membrane raft, membrane microdomain, membrane region, and focal adhesion (Figure 8C). Moreover, MF analysis revealed that downstream genes were mostly enriched in protein serine/threonine kinase activity (Figure 8D). KEGG analysis demonstrated that axon guidance, hippo signaling pathway, TGF-beta signaling pathway, cellular senescence, and steroid biosynthesis were the top five enriched pathways (Figure 8E).

Screening of Hub Genes and Survival-Related Genes

Through STRING database, a total of 41 genes were screened out of 107 target genes to construct the PPI network, which contained 41 nodes and 42 edges (Figure 8F). The 10 hub genes (CAV1, BIRC5, CKS1B, SMAD7, BMPR2, PTPN1, BMP2, EDN1, TGFBR2, and CXCL12) were screened out by the rank of degree which was calculated by cytoHubba in cytoscape (Figure 8G). The KM curve showed that the expression of BIRC5 and CKS1B were negatively correlated with survival prognosis (p < 0.01, Figures 8H,I).

Validation of the Expression Level of Ferroptosis-Related Long Non-coding RNAs in Lung Adenocarcinoma Samples

To validate the expression levels of ferroptosis-related lncRNAs, we detected two ferroptosis-related lncRNAs expression levels in 20 LUAD samples and 20 adjacent normal tissues through qRT-PCR. The results showed that AC021016.1 was significantly downregulated in LUAD while KTN1-AS1 was highly expressed in LUAD (Supplementary Figure 4).

Discussion

Long non-coding RNAs (lncRNAs) are a kind of non-protein coding RNAs with over 200 nucleotides (Wang J. et al., 2021). They are widely expressed in human cells and play vital roles in regulation of pathological and physiological processes in various types of malignant tumors (Zhou et al., 2020). Previous studies have constructed a few of lncRNA prognosis signatures, like m6A-related lncRNA signature (Xu et al., 2021), immune-related lncRNA signature (Wang J. et al., 2021), and autophagy-related lncRNA signature (Zhou et al., 2020), which can be used to identify high-risk patients and judge the prognosis of them. Ferroptosis is a form of regulatory cell death and plays an important role in the occurrence and development of LUAD. Some studies have shown that lncRNAs can participate in the regulation of ferroptosis and growth of LUAD (Sánchez et al., 2014; Johnson et al., 2017; Zhang et al., 2017; Wang et al., 2019; Gai et al., 2020). Thus, ferroptosis-related lncRNAs are potential markers for prognosis and a prognostic model based on them for LUAD is still lacking.

To construct a prognostic model, we firstly uncovered 18 ferroptosis-related genes through gene expression differential analysis and Cox regression. Among them, cysteine dioxygenase 1 (CDO1), toll-like receptor 4 (TLR4), and dual oxidase 1 (DUOX1) were downregulated in tumors and played a protective role, while Fanconi anemia complementation group D2 (FANCD2) and others were with high expression in tumors and were risk factors for prognosis. CDO1 was found to transform cysteine to taurine while cysteine was an indispensable substrate of glutathione peroxidase 4 (GPX4), a lipid repair enzyme of ferroptosis. The suppression of CDO1 increased the expression of GSH and inhibited ferroptosis in gastric cancer (Hao et al., 2017). In another study (Chen et al., 2019), Chen et al. reported that the silence of TLR4 and NOX4 significantly retarded the autophagy and ferroptosis in rats with heart failure. FANCD2, a nuclear protein that responds to DNA damage repair, negatively regulated ferroptosis of bone marrow injury in cancer treatment and could be a potential target of anticancer therapies (Song et al., 2016). In addition, a 10 ferroptosis-related genes signature and a ferroptosis-related gene signature with five genes were constructed to predict the prognosis of patients with LUAD (Wang S. et al., 2021; Zhu et al., 2021). In general, ferroptosis-related genes play important roles in the progression of LUAD.

It has been reported that ferroptosis-related lncRNAs were prognostic risk factors and the signatures constructed with ferroptosis-related lncRNAs can be used to distinguish high-risk patients from low-risk patients (Cai et al., 2021; Tang et al., 2021). However, the full role of ferroptosis-related lncRNAs in LUAD remains to be further explored. Here, we mined ferroptosis-related lncRNAs through Pearson correlation and constructed a prognostic Signature with six ferroptosis-related lncRNAs. Through Kaplan–Meier survival curves, ROC curve, PCA, and univariate and multivariate COX regression, we further confirmed the validity of our signature. Among the six lncRNAs, AC021016.1 and MIR223HG were downregulated in LUAD and were protective factors for prognosis, while KTN1-AS1 and others were upregulated and were risk factors for prognosis. LncRNA KTN1-AS1 was upregulated in several tumors, such as non-small cell lung cancer, bladder cancer, head and neck squamous cell carcinoma, and hepatocellular carcinoma and promoted the progression of tumors (Zhang et al., 2019; Jiang et al., 2020; Hu et al., 2021; Li Y. Y. et al., 2021). In our study, we find that lncRNA KTN1-AS1 was significantly upregulated in LUAD and correlated with tumor stemness and tumor microenvironment. In another study (Geng et al., 2021), MIR223HG was found to be a genome instability-related lncRNA and used to construct a genome instability-related lncRNA signature, which further confirmed the potential function of MIR223HG in LUAD. AL049555.1 was demonstrated to be highly expressed and associated with high risk of PAAD (Wang et al., 2020). The function of AC021016.1, AC068228.2, and AC009275.1 has not been elucidated and our research uncovered their potential roles in LUAD.

To further elucidate the potential function of six ferroptosis-related lncRNAs, we constructed ceRNA network to find the downstream genes. Through STRING and Cytoscape, we found the 10 hub genes of PPI network. Among them, the high expression of BIRC5 and CKS1B were significantly correlated with poor prognosis. BIRC5 is a well-known therapeutic target of cancers (Li et al., 2019) and plays an important role in cell division and inhibits cell death (Wheatley and Altieri, 2019). Previous research revealed that BIRC5 was elevated and promoted the development of lung cancer (Han et al., 2020). CKS1B, a member of Cks/Suc1 family, is involved in the regulation of cell cycle and chemotherapeutic resistance in cancers (Wang et al., 2017). In lung cancer, it also promotes cell growth and induces drug resistance (Shi et al., 2020). In addition, although CAV1, SMAD7, BMP2, EDN1, and CXCL12 were not significantly different in the prognostic analysis in our study, they also play important roles in the occurrence and development of lung cancer (Katsura et al., 2018; Yan et al., 2019; Huang et al., 2020; Pulido et al., 2020; Tong et al., 2020). In general, the functions of target genes further confirmed the potential roles of those lncRNAs. Finally, we detected two ferroptosis-related lncRNAs expression levels through qRT-PCR with 20 pairs of samples, which increased the reliability of our results.

Ferroptosis is a new form of cell death that might potentially provide a new strategy in clinical cancer therapy. However, the crosstalk among ferroptosis and other cell death processes, auto-immune microenvironment, and maintenance and transformation of tumor stem cells remain unsolved. In this study, we constructed a ferroptosis-related lncRNA signature and confirmed its validity in predicting survival outcomes of LUAD. We hope that these findings will provide some useful insights for clinical practice. However, several limitations were in our study. Firstly, multicenter LUAD datasets are needed to verify the validity of our model. Secondly, although we confirmed the expression in our own samples through qRT-PCR, the application of this prognostic prediction model in LUAD need to be further explored.

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/s.

Ethics Statement

The studies involving human participants were reviewed and approved by the Medical Ethics Committee of Changhai Hospital, Navy Military Medical University. The patients/participants provided their written informed consent to participate in this study.

Author Contributions

XF, CH, and XW constructed this study. XF, CH, XW, CJL, and HC collected and analyzed the data. BS and CL wrote and revised the manuscript and were responsible for the critical reading of the manuscript. All authors contributed to the article and approved the submitted version.

Funding

This study was sponsored by Natural Science Foundation of Shanghai (19ZR1456400) and National Natural Science Foundation of China (81802288).

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s Note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Acknowledgments

We sincerely appreciate the contribution of the TCGA project.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fcell.2021.751490/full#supplementary-material

Footnotes

1. ^ https://portal.gdc.cancer.gov/
2. ^ https://xenabrowser.net/datapages/
3. ^ http://carolina.imis.athena-innovation.gr
4. ^ http://www.targetscan.org
5. ^ http://www.mirdb.org/miRDB/
6. ^ http://mirwalk.umm.uni-heidelberg.de/
7. ^ https://string-db.org/

References

Agarwal, V., Bell, G. W., Nam, J. W., and Bartel, D. P. (2015). Predicting effective microRNA target sites in mammalian mRNAs. Elife 4:e05005.

Google Scholar

Bray, F., Ferlay, J., Soerjomataram, I., Siegel, R. L., Torre, L. A., and Jemal, A. (2018). Global cancer statistics 2018: GLOBOCAN estimates of incidence and mortality worldwide for 36 cancers in 185 countries. CA Cancer J. Clin. 68, 394–424. doi: 10.3322/caac.21492

PubMed Abstract | CrossRef Full Text | Google Scholar

Cai, H. J., Zhuang, Z. C., Wu, Y., Zhang, Y. Y., Liu, X., Zhuang, J. F., et al. (2021). Development and validation of a ferroptosis-related lncRNAs prognosis signature in colon cancer. Bosn. J. Basic Med. Sci. 21, 569–576. doi: 10.17305/bjbms.2020.5617

PubMed Abstract | CrossRef Full Text | Google Scholar

Chen, X., Kang, R., Kroemer, G., and Tang, D. (2021). Broadening horizons: the role of ferroptosis in cancer. Nat. Rev. Clin. Oncol. 18, 280–296. doi: 10.1038/s41571-020-00462-0

PubMed Abstract | CrossRef Full Text | Google Scholar

Chen, X., Xu, S., Zhao, C., and Liu, B. (2019). Role of TLR4/NADPH oxidase 4 pathway in promoting cell death through autophagy and ferroptosis during heart failure. Biochem. Biophys. Res. Commun. 516, 37–43. doi: 10.1016/j.bbrc.2019.06.015

PubMed Abstract | CrossRef Full Text | Google Scholar

Chen, Y., and Wang, X. (2020). miRDB: an online database for prediction of functional microRNA targets. Nucleic Acids Res. 48, D127–D131.

Google Scholar

Chu, B., Kon, N., Chen, D., Li, T., Liu, T., Jiang, L., et al. (2019). ALOX12 is required for p53-mediated tumour suppression through a distinct ferroptosis pathway. Nat. Cell Biol. 21, 579–591. doi: 10.1038/s41556-019-0305-6

PubMed Abstract | CrossRef Full Text | Google Scholar

Dweep, H., and Gretz, N. (2015). miRWalk2.0: a comprehensive atlas of microRNA-target interactions. Nat. Methods 12:697. doi: 10.1038/nmeth.3485

PubMed Abstract | CrossRef Full Text | Google Scholar

Gai, C., Liu, C., Wu, X., Yu, M., Zheng, J., Zhang, W., et al. (2020). MT1DP loaded by folate-modified liposomes sensitizes erastin-induced ferroptosis via regulating miR-365a-3p/NRF2 axis in non-small cell lung cancer cells. Cell Death Dis. 11:751. doi: 10.1038/s41419-020-02939-3

PubMed Abstract | CrossRef Full Text | Google Scholar

Geng, W., Lv, Z., Fan, J., Xu, J., Mao, K., Yin, Z., et al. (2021). Identification of the prognostic significance of somatic mutation-derived LncRNA signatures of genomic instability in lung adenocarcinoma. Front. Cell Dev. Biol. 9:657667. doi: 10.3389/fcell.2021.657667

PubMed Abstract | CrossRef Full Text | Google Scholar

Han, F., Yang, S., Wang, W., Huang, X., Huang, D., and Chen, S. (2020). Silencing of lncRNA LINC00857 enhances BIRC5-dependent radio-sensitivity of lung adenocarcinoma cells by recruiting NF-κB1. Mol. Ther. Nucleic Acids 22, 981–993. doi: 10.1016/j.omtn.2020.09.020

PubMed Abstract | CrossRef Full Text | Google Scholar

Hao, S., Yu, J., He, W., Huang, Q., Zhao, Y., Liang, B., et al. (2017). Cysteine dioxygenase 1 mediates erastin-induced ferroptosis in human gastric cancer cells. Neoplasia 19, 1022–1032. doi: 10.1016/j.neo.2017.10.005

PubMed Abstract | CrossRef Full Text | Google Scholar

Hu, X., Xiang, L., He, D., Zhu, R., Fang, J., Wang, Z., et al. (2021). The long noncoding RNA KTN1-AS1 promotes bladder cancer tumorigenesis via KTN1 cis-activation and the consequent initiation of Rho GTPase-mediated signaling. Clin. Sci. (Lond.) 135, 555–574. doi: 10.1042/CS20200908

PubMed Abstract | CrossRef Full Text | Google Scholar

Huang, F., Cao, Y., Wu, G., Chen, J., Wang, C., Lin, W., et al. (2020). BMP2 signalling activation enhances bone metastases of non-small cell lung cancer. J. Cell Mol. Med. 24, 10768–10784.

Google Scholar

Jiang, Y., Wu, K., Cao, W., Xu, Q., Wang, X., Qin, X., et al. (2020). Long noncoding RNA KTN1-AS1 promotes head and neck squamous cell carcinoma cell epithelial-mesenchymal transition by targeting miR-153-3p. Epigenomics 12, 487–505. doi: 10.2217/epi-2019-0173

PubMed Abstract | CrossRef Full Text | Google Scholar

Johnson, G. S., Li, J., Beaver, L. M., Dashwood, W. M., Sun, D., Rajendran, P., et al. (2017). A functional pseudogene, NMRAL2P, is regulated by Nrf2 and serves as a coactivator of NQO1 in sulforaphane-treated colon cancer cells. Mol. Nutr. Food Res. 61:10. doi: 10.1002/mnfr.201600769

PubMed Abstract | CrossRef Full Text | Google Scholar

Katsura, M., Shoji, F., Okamoto, T., Shimamatsu, S., Hirai, F., Toyokawa, G., et al. (2018). Correlation between CXCR4/CXCR7/CXCL12 chemokine axis expression and prognosis in lymph-node-positive lung cancer patients. Cancer Sci. 109, 154–165. doi: 10.1111/cas.13422

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, C. H., and Chen, Y. (2016). Insight into the role of long noncoding RNA in cancer development and progression. Int. Rev. Cell Mol. Biol. 326, 33–65. doi: 10.1016/bs.ircmb.2016.04.001

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, F., Aljahdali, I., and Ling, X. (2019). Cancer therapeutics using survivin BIRC5 as a target: what can we do after over two decades of study? J. Exp. Clin. Cancer Res. 38:368. doi: 10.1186/s13046-019-1362-1

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, S., Yang, S., Qiu, C., and Sun, D. (2021). LncRNA MSC-AS1 facilitates lung adenocarcinoma through sponging miR-33b-5p to up-regulate GPAM. Biochem. Cell. Biol. 99, 241–248. doi: 10.1139/bcb-2020-0239

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, Y. Y., Li, W., Chang, G. Z., and Li, Y. M. (2021). Long noncoding RNA KTN1 antisense RNA 1exerts an oncogenic function in lung adenocarcinoma by regulating DEP domain containing 1 expression via activating epithelial-mesenchymal transition. Anticancer Drugs 32, 614–625. doi: 10.1097/CAD.0000000000001035

PubMed Abstract | CrossRef Full Text | Google Scholar

Liu, Y., Zhang, X., Zhang, J., Tan, J., Li, J., and Song, Z. (2020). Development and validation of a combined ferroptosis and immune prognostic classifier for hepatocellular carcinoma. Front. Cell Dev. Biol. 8:596679. doi: 10.3389/fcell.2020.596679

PubMed Abstract | CrossRef Full Text | Google Scholar

Miller, K. D., Nogueira, L., Mariotto, A. B., Rowland, J. H., Yabroff, K. R., Alfano, C. M., et al. (2019). Cancer treatment and survivorship statistics, 2019. CA Cancer J. Clin. 69, 363–385.

Google Scholar

Paraskevopoulou, M. D., Vlachos, I. S., Karagkouni, D., Georgakilas, G., Kanellos, I., Vergoulis, T., et al. (2016). DIANA-LncBase v2: indexing microRNA targets on non-coding transcripts. Nucleic Acids Res. 44, D231–D238. doi: 10.1093/nar/gkv1270

PubMed Abstract | CrossRef Full Text | Google Scholar

Pulido, I., Ollosi, S., Aparisi, S., Becker, J. H., Aliena-Valero, A., Benet, M., et al. (2020). Endothelin-1-mediated drug resistance in EGFR-mutant non-small cell lung carcinoma. Cancer Res. 80, 4224–4232. doi: 10.1158/0008-5472.CAN-20-0141

PubMed Abstract | CrossRef Full Text | Google Scholar

Qu, S., Jiao, Z., Lu, G., Yao, B., Wang, T., Rong, W., et al. (2021). PD-L1 lncRNA splice isoform promotes lung adenocarcinoma progression via enhancing c-Myc activity. Genome Biol. 22:104. doi: 10.1186/s13059-021-02331-0

PubMed Abstract | CrossRef Full Text | Google Scholar

Sánchez, Y., Segura, V., Marín-Béjar, O., Athie, A., Marchese, F. P., González, J., et al. (2014). Genome-wide analysis of the human p53 transcriptional network unveils a lncRNA tumour suppressor signature. Nat. Commun. 5:5812. doi: 10.1038/ncomms6812

PubMed Abstract | CrossRef Full Text | Google Scholar

Shannon, P., Markiel, A., Ozier, O., Baliga, N. S., Wang, J. T., Ramage, D., et al. (2003). Cytoscape: a software environment for integrated models of biomolecular interaction networks. Genome Res. 13, 2498–2504. doi: 10.1101/gr.1239303

PubMed Abstract | CrossRef Full Text | Google Scholar

Shi, W., Huang, Q., Xie, J., Wang, H., Yu, X., and Zhou, Y. (2020). CKS1B as drug resistance-inducing gene-a potential target to improve cancer therapy. Front. Oncol. 10:582451. doi: 10.3389/fonc.2020.582451

PubMed Abstract | CrossRef Full Text | Google Scholar

Song, X., Xie, Y., Kang, R., Hou, W., Sun, X., Epperly, M. W., et al. (2016). FANCD2 protects against bone marrow injury from ferroptosis. Biochem. Biophys. Res. Commun. 480, 443–449. doi: 10.1016/j.bbrc.2016.10.068

PubMed Abstract | CrossRef Full Text | Google Scholar

Szklarczyk, D., Gable, A. L., Nastou, K. C., Lyon, D., Kirsch, R., Pyysalo, S., et al. (2021). The STRING database in 2021: customizable protein-protein networks, and functional characterization of user-uploaded gene/measurement sets. Nucleic Acids Res. 49, D605–D612.

Google Scholar

Tang, Y., Li, C., Zhang, Y. J., and Wu, Z. H. (2021). Ferroptosis-related long non-coding RNA signature predicts the prognosis of Head and neck squamous cell carcinoma. Int. J. Biol. Sci. 17, 702–711. doi: 10.7150/ijbs.55552

PubMed Abstract | CrossRef Full Text | Google Scholar

Tong, L., Shen, S., Huang, Q., Fu, J., Wang, T., Pan, L., et al. (2020). Proteasome-dependent degradation of Smad7 is critical for lung cancer metastasis. Cell Death Differ. 27, 1795–1806. doi: 10.1038/s41418-019-0459-6

PubMed Abstract | CrossRef Full Text | Google Scholar

Tu, Z., Wu, L., Wang, P., Hu, Q., Tao, C., Li, K., et al. (2020). N6-methylandenosine-related lncRNAs are potential biomarkers for predicting the overall survival of lower-grade glioma patients. Front. Cell Dev. Biol. 8:642. doi: 10.3389/fcell.2020.00642

PubMed Abstract | CrossRef Full Text | Google Scholar

Wang, H., Sun, M., Guo, J., Ma, L., Jiang, H., Gu, L., et al. (2017). 3-O-(Z)-coumaroyloleanolic acid overcomes Cks1b-induced chemoresistance in lung cancer by inhibiting Hsp90 and MEK pathways. Biochem. Pharmacol. 135, 35–49. doi: 10.1016/j.bcp.2017.03.007

PubMed Abstract | CrossRef Full Text | Google Scholar

Wang, J., Xiang, J., and Li, X. (2020). Construction of a competitive endogenous RNA network for pancreatic adenocarcinoma based on weighted gene co-expression network analysis and a prognosis model. Front. Bioeng. Biotechnol. 8:515. doi: 10.3389/fbioe.2020.00515

PubMed Abstract | CrossRef Full Text | Google Scholar

Wang, J., Yin, X., Zhang, Y. Q., and Ji, X. (2021). Identification and validation of a novel immune-related four-lncRNA signature for lung adenocarcinoma. Front. Genet. 12:639254. doi: 10.3389/fgene.2021.639254

PubMed Abstract | CrossRef Full Text | Google Scholar

Wang, M., Mao, C., Ouyang, L., Liu, Y., Lai, W., Liu, N., et al. (2019). Long noncoding RNA LINC00336 inhibits ferroptosis in lung cancer by functioning as a competing endogenous RNA. Cell Death Differ. 26, 2329–2343.

Google Scholar

Wang, S., Wu, C., Ma, D., and Hu, Q. (2021). Identification of a ferroptosis-related gene signature (FRGS) for predicting clinical outcome in lung adenocarcinoma. PeerJ 9:e11233. doi: 10.7717/peerj.11233

PubMed Abstract | CrossRef Full Text | Google Scholar

Wang, Z., Zhang, X., Tian, X., Yang, Y., Ma, L., Wang, J., et al. (2021). CREB stimulates GPX4 transcription to inhibit ferroptosis in lung adenocarcinoma. Oncol. Rep. 45:88. doi: 10.3892/or.2021.8039

PubMed Abstract | CrossRef Full Text | Google Scholar

Wheatley, S. P., and Altieri, D. C. (2019). Survivin at a glance. J. Cell Sci. 132:jcs223826. doi: 10.1242/jcs.223826

PubMed Abstract | CrossRef Full Text | Google Scholar

Xu, F., Huang, X., Li, Y., Chen, Y., and Lin, L. (2021). m⁶A-related lncRNAs are potential biomarkers for predicting prognoses and immune responses in patients with LUAD. Mol. Ther. Nucleic Acids 24, 780–791. doi: 10.1016/j.omtn.2021.04.003

PubMed Abstract | CrossRef Full Text | Google Scholar

Yan, Y., Xu, Z., Qian, L., Zeng, S., Zhou, Y., Chen, X., et al. (2019). Identification of CAV1 and DCN as potential predictive biomarkers for lung adenocarcinoma. Am. J. Physiol. Lung Cell Mol. Physiol. 316, L630–L643. doi: 10.1152/ajplung.00364.2018

PubMed Abstract | CrossRef Full Text | Google Scholar

Yang, G., Zhang, Y., and Yang, J. (2019). A five-microRNA signature as prognostic biomarker in colorectal cancer by bioinformatics analysis. Front. Oncol. 9:1207. doi: 10.3389/fonc.2019.01207

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhang, D., Zhang, G., Hu, X., Wu, L., Feng, Y., He, S., et al. (2017). Oncogenic RAS regulates long noncoding RNA Orilnc1 in human cancer. Cancer Res. 77, 3745–3757. doi: 10.1158/0008-5472.CAN-16-1768

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhang, L., Wang, L., Wang, Y., Chen, T., Liu, R., Yang, W., et al. (2019). LncRNA KTN1-AS1 promotes tumor growth of hepatocellular carcinoma by targeting miR-23c/ERBB2IP axis. Biomed. Pharmacother. 109, 1140–1147. doi: 10.1016/j.biopha.2018.10.105

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhou, M., Shao, W., Dai, H., and Zhu, X. (2020). A robust signature based on autophagy-associated LncRNAs for predicting prognosis in lung adenocarcinoma. Biomed. Res. Int. 2020:3858373. doi: 10.1155/2020/3858373

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhou, Y., Zhou, B., Pache, L., Chang, M., Khodabakhshi, A. H., Tanaseichuk, O., et al. (2019). Metascape provides a biologist-oriented resource for the analysis of systems-level datasets. Nat. Commun. 10:1523. doi: 10.1038/s41467-019-09234-6

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhu, G., Huang, H., Xu, S., Shi, R., Gao, Z., Lei, X., et al. (2021). Prognostic value of ferroptosis-related genes in patients with lung adenocarcinoma. Thorac. Cancer 12, 1890–1899.

Google Scholar