N6-Methyladenosine-Related Long Non-Coding RNAs Are Identified as a Potential Prognostic Biomarker for Lung Squamous Cell Carcinoma and Validated by Real-Time PCR

Currently, the precise mechanism by which N⁶-methyladenosine (m⁶A) modification of long non-coding RNAs (lncRNAs) promotes the occurrence and development of lung squamous cell carcinoma (LUSC) and influences tumor microenvironment (TME) remains unclear. Therefore, we studied the prognostic value of m⁶A-related lncRNAs and their relationship with TME in 495 LUSC samples from The Cancer Genome Atlas (TCGA) database. Pearson’s correlation and univariate Cox regression analysis identified 6 m⁶A-related lncRNAs with prognostic values for LUSC patients. LUSC patients were divided into two subgroups (clusters 1 and 2) using principal component analysis. The expression of PD-L1 was lower in tumor tissues and cluster 2 of LUSC patients. Cluster 2 of LUSC patients had a high immune score, stromal score, and unique immune cell infiltration. The focal adhesion kinase (FAK) pathway and cytokine receptor pathways are enriched in cluster 1. The m⁶A-related lncRNA prognostic markers (m⁶A-LPMs) were established using the least absolute shrinkage and selection operator (LASSO) Cox regression analysis. The risk score was calculated by 4 m⁶A-LPMs and associated with OS, TME, clinicopathological characteristics of LUSC patients. After adjusting for age, gender, and stage, the risk score was also an independent prognostic factor for LUSC patients. Real-time PCR results showed that the expression of 4 m⁶A-LPMs was consistent with our prediction results. Our study found that 4 m⁶A-LPMs (AC138035.1, AC243919.2, HORMAD2-AS1, and AL122125.1) are closely associated with LUSC prognosis, in future, they may as novel diagnostic biomarkers for LUSC and provide new immunotherapy targets for LUSC patients.

1 Introduction

Lung cancer is currently the second leading cause of global cancer incidence and the primary cause of cancer death in the worldwide. In 2020, about 2.2 million new lung cancer cases were diagnosed, and another 1.8 million died from lung cancer worldwide (Sung et al., 2021). Non-small cell lung cancer (NSCLC) is the most common type of lung cancer, accounting for about 85%, which mainly includes lung adenocarcinoma (LUAD) and lung squamous cell carcinoma (LUSC) (Woodman et al., 2021). Although significant progress has been made in lung cancer treatment modalities, including surgery, radiotherapy, chemotherapy, and targeted therapy, the 5-years survival rate of NSCLC patients is still low, at about 16.6% (Torre et al., 2016; Liu et al., 2020). Various factors determine the development and progression of NSCLC. The genetic makeup of tumor cells and the tumor microenvironment (TME) have critical regulatory roles (Forde et al., 2014). The immune system, in particular, plays a dual role in the evolution of tumors. It can identify and control new tumor cells through immune surveillance and promote tumor progression through various mechanisms of immunosuppression (Kunimasa and Goto, 2020). As a result, immunotherapy has gradually shifted from second-line therapy to first-line therapy since 2013, and patients without targeted oncogene mutations have benefited from immunotherapy (Zhang et al., 2019). Paclitaxel and carboplatin combined with immune checkpoint inhibitors such as pembrolizumab is the first-line treatment for LUSC at present (Fan et al., 2019; Lin et al., 2020). This method can significantly prolong the overall survival and progression-free survival of LUSC patients (Paz-Ares et al., 2018). However, the long-term follow-up shows that patients treated with immunotherapy may eventually develop resistance to immunotherapy, and in some patients tumor progression may even be accelerated after immunotherapy (Huang et al., 2015; Ribas, 2015), which may be related to the imbalance of the immune system in LUSC patients. Therefore, it is necessary to investigate the regulatory network of the immune TME to identify new biomarkers and targets for the prognosis and immunotherapy of patients with LUSC.

N⁶-methyladenosine (m⁶A) methylation is one of the most common post-transcriptional modifications. It is reported to play an essential role in normal physiological processes and pathological processes leading to various diseases (Ma et al., 2019). M⁶A methylation modification is the methylation of the sixth N atom of adenine (A) in RNA, referred to as m⁶A, which mainly occurs in mRNA and accounts for about 80% of all RNA methylation (Ping et al., 2014). Recent studies have found that m⁶A modifications can also regulate the generation and function of non-coding RNA, such as micro-RNAs (miRNAs), long non-coding RNAs (lncRNAs), and circular RNAs (circRNAs) (Yang et al., 2017; Wu et al., 2018; Yang et al., 2018). The m⁶A modification of RNA is a dynamic and reversible mode of regulation. It is co-regulated by m⁶A methyltransferase (writers) and demethylase (erasers), and can be selectively recognized and bound by m⁶A recognition proteins (readers) to regulate gene expression post-transcription (Zaccara et al., 2019). Some studies have found that m⁶A promotes tumor progression by upregulating oncogenes in various tumors (He et al., 2019). Liu et al. (2020) performed a reliable cluster analysis on the gene expression of 13 m⁶A-related proteins and divided LUAD patients into two groups with significant differences in age, race, tumor size, stage, and distant metastasis, and their expression profile analysis showed that the expression of m⁶A-related genes could be used as the prognostic criterion for patients with advanced LUAD. Some studies have also found that the m⁶A demethylase FTO promotes the tumor progression of LUSC by regulating the expression of MZF1 and provides a potential target for LUSC treatment (Liu et al., 2018). These findings indicate that m⁶A methylation plays a vital role in the progression and prognosis of LUSC.

LncRNAs are a class of non-coding RNAs with a length of more than 200 nucleotides, which do not encode proteins but have some characteristics of mRNA (Cabili et al., 2011). Their abnormal expression in lung cancer can promote tumor progression. For example, MALAT1, which has carcinogenic effects, is the commonest lncRNA in lung cancer, which can promote NSCLC tumor growth and metastasis by regulating the miR124/STAT3 and miR204/SLUG axis (Li et al., 2016; Li et al., 2018). The interaction between m⁶A and lncRNAs has two modes: lncRNAs regulating m⁶A-modified RNAs and m⁶A-modified lncRNAs (Dai et al., 2020). M⁶A has extensive modifications to lncRNAs that functionally regulate the eukaryotic transcriptome, affecting RNA splicing, export, localization, translation, and stability (Linder et al., 2015). Some studies have found that m⁶A modification can regulate the stability or expression of lncRNAs, which regulate tumor progression through the lncRNA-mediated ceRNA network (Fazi and Fatica, 2019; Chen J. et al., 2020; Huang et al., 2020). For instance, there is a triple helix structure at the 3 ‘terminal of lncRNA MALAT-1, and the deletion or mutation of this structure affects the stability and expression of MALAT1 (Brown et al., 2014). Methyltransferase METTL16 can bind to the triple helix structure at the 3 ‘terminal of lncRNA MALAT1 and influence its structural stability and functional expression, thus participating in the occurrence and development of tumor (Brown et al., 2016). Studies have reported that METTL3 can improve the methylation level of m⁶A and enhance the transcriptional stability of lncRNA ABHD11-AS1 to increase its expression. At the same time, the overexpression of METTL3 is closely related to the poor prognosis of patients with NSCLC (Xue et al., 2021). In recent years, some studies have also found that m⁶A-related lncRNAs are potential biomarkers for predicting the prognosis and immune response in LUAD patients (Xu et al., 2021; Zheng et al., 2021). However, the mechanism how the lncRNAs modificated by m⁶A promote the occurrence and development of LUSC remains unclear. Therefore, to understand the regulating role of m⁶A modificated lncRNAs in the progression of LUSC may help us to get new biomarkers and valuable therapeutic targets.

This study uses Pearson correlation analysis to identify potential m⁶A-related lncRNAs. Then, univariate Cox regression analysis identified m⁶A-related lncRNAs with prognostic values. We analyzed their expression in lung tumor tissues and normal tissues datasets through bioinformatics and statistical methods. We also established cluster subtypes of m⁶A-related lncRNAs with prognostic values and systematically assessed the correlation between these lncRNAs in different clusters with clinicopathological characteristics, prognosis, PD-L1, and TME. Gene enrichment analysis was performed to evaluate the possible reasons for the differences in TME of the different clusters. Additionally, we identified 4 m⁶A-related lncRNA prognostic markers (m⁶A-LPMs) based on the prognostic value of m⁶A-related lncRNAs using the LASSO Cox analysis. We randomly divided 495 LUSC patients into the training and test groups. In addition, we constructed a risk model and calculated the risk score of each patient based on the expression and regression coefficients of the four m⁶A-LPMs in the different groups. According to the median risk score, all patients were divided into high-risk and low-risk groups. We evaluate the relationship between the clinicopathological characteristics, prognosis, immune score, and immune cell infiltration in the high- and low-risk groups. Real-time PCR was used to further verify the expression of m⁶A-LPMs in LUSC cell line HCC1588 and human bronchial epithelial cell line 16HBE. The Study flowchart of this research is shown in Figure 1.

FIGURE 1

FIGURE 1. Study flowchart of this research.

2 Materials and Methods

2.1 Datasets and m⁶A-Related Genes

The RNA-seq transcriptome data and related clinical data of LUSC patients were downloaded from the TCGA data portal (https://portal.gdc.cancer.gov/). Data from 504 LUSC samples and 49 adjacent normal tissues were downloaded from the TCGA database. After excluding samples with unknown survival time and survival status from clinical features, we ultimately obtained 498 LUSC patient samples and their corresponding clinicopathological information, including survival time, survival status, age, sex, and TNM staging for further analysis. We extracted 23 m⁶A-related gene expression matrices from the TCGA dataset based on the data in previous publications, including the expression data of writers (METTL3, METTL14, METTL16, WTAP, VIRMA, ZC3H13, RBM15, and RBM15B), readers (YTHDC1, YTHDC2, YTHDF1, YTHDF2, YTHDF3, HNRNPC, FMR1, LRPPRC, HNRNPA2B1, IGFBP1, IGFBP2, IGFBP3, RBMX) and erasers (FTO and ALKBH5) (Wojtas et al., 2017; Tang et al., 2018; Yu et al., 2021; Wang et al., 2021; and Xu et al., 2021.).

2.2 Annotation of lncRNAs

Obtain The annotation file of the lncRNAs of the Genome Reference Consortium Human Build 38 (GRCh38) from the GENCODE website (https://www.gencodegenes.org/human) to annotate the lncRNAs in the TCGA database. 14,086 lncRNAs were identified in the LUSC datasets from the TCGA database by identifying the Ensemble IDs of the gene.

2.3 Identification of m⁶A-Related lncRNAs and Construction of a Co-expression Network

We acquired detailed clinical information for 498 LUSC patients from the extracted clinical data in the TCGA database. We assessed survival outcomes of the LUSC patients based on the overall survival (OS). Through the “limma” package (version 3.42.2) of the R language, we use Pearson correlation analysis to identify potential m⁶A-related lncRNAs based on a correlation coefficient |R| = 0.4 (p = 0.001). Then, the “igraph” R package was used to draw a co-expression network of the m⁶A-related lncRNAs and m⁶A-related genes (Figure 1).

2.4 Screening of Prognostic-Related lncRNAs and Their Expression Heatmap

First, the “limma” R package (version 3.42.2) merged the expression of m⁶A-related lncRNAs with the survival time and status of 498 LUSC patients. We obtained the survival outcomes and the expression of m⁶A-related lncRNAs of 495 LUSC patients. Then, univariate Cox regression analysis identified m⁶A-related lncRNAs with prognostic values (p < 0.05) (Figure 1). Finally, the “pheatmap,” “reshape2,” and “ggpubr” R packages (versions 1.0.12, 1.4.4, and 0.4.0) were used to draw the expression heatmap of prognostic-related lncRNAs in lung cancer tumor tissues and normal tissues. We also classified the LUSC patients into different groups using the “limma” and the “ConsensusClusterPlus” R packages (versions 3.42.2 and 1.58.0) based on the expression levels of m⁶A-related lncRNAs with prognostic values.

2.5 Survival Analysis and Clinical Correlation Analysis Between Clusters 1 and 2

We used the “survival” and “survminer” R packages (versions 3.2–11 and 0.4.9) to compare survival outcomes between the two clusters of LUSC patients. Then, the “pheatmap” R package (version 1.0.12) drew expression heatmaps of different LUSC clusters and clinical features.

2.6 Comparison of Immune Score, Stromal Score, and TME Between the Two Clusters

The immune score of each patient with LUSC was calculated using the estimation algorithm in the “estimate” R package (version 1.0.13). Then, we estimated the cell type identification of the relative subset of RNA transcripts and obtained the ratio of 22 immune cell types for each sample using the “e1071”, “preprocessCore” packages (versions 1.7-6 and 1.48.0), and " CIBERSORT” R script (version 1.03). In this process, we used 1,000 sorting algorithms, and only performed follow-up analysis on samples with a CIBERSORT value < 0.05, and compared the differential immune penetration levels between two clusters.

2.7 Gene Set Enrichment Analysis

The differences in TME between two clusters of LUSC were analyzed using the gene set enrichment analysis (GSEA). p < 0.05 and a false discovery rate of <0.05 were considered statistically significant.

2.8 Construction of m⁶A-Related lncRNA Prognostic Markers

The 495 LUSC patients were randomly divided into the training group (248 patients) and the test group (247 patients) with a 1:1 ratio by using the “caret” R package (version 6.0–86). We established 4 m⁶A-related lncRNA prognostic markers (m⁶A-LPMs) by Least Absolute Shrinkage and Selection Operator (LASSO) Cox regression analysis (using 10-fold cross-validation to estimate the penalty parameters) using the “glmnet” R package (version 4.1–1). Using the coefficients generated by the LASSO regression algorithm, we obtained the risk score equation: risk score = sum of coefficients×m⁶A-related-LPMs expression level, which was used to calculate the risk score of each patient in the training and test group. Then, the median of the risk score was set as the standard, and the patients were divided into high-risk and low-risk groups (Figure 1).

2.9 Survival Analysis Between High-Risk and Low-Risk Groups

According to the risk score, using the “survival,” “survminer,” “timeROC,” and “pheatmap” R packages (versions 3.2–11, 0.4.9, 0.4 and 1.0.12) to compare survival outcomes between high-risk and low-risk groups.

2.10 Correlation Analysis of the Risk Score and Clinicopathological Characteristics

The risk scores and clinicopathological characteristics of all patients were integrated using the “survival,” “survminer,” and “pheatmap” R packages (versions 3.2–11, 0.4.9, and 1.0.12) to draw heatmaps and survival curves of the correlation between different risk scores and clinicopathological characteristics.

2.11 Independent Prognostic Analysis

After integrating the risk scores of the training group and the test group with clinical data, we used univariate and multivariate Cox regression analysis to assess whether the risk score can be an independent prognostic factor for LUSA patients. The “survival” R package (version 3.2–11) was used to draw a forest map of Cox regression analysis.

2.12 Correlation Analysis Between Risk Score and Immunity

We used the “limma,” “ggplot2”, “ggpubr,” and “ggExtra” R packages (versions 3.42.2, 3.3.5, 0.4.0, and 0.9) to integrate the risk score and the content of immune cells in order to analyze the correlation between the risk score and the expression level of immune cells.

2.13 Cell Culture

The human bronchial epithelial cell line 16HBE was cultured in Dulbecco’s modified Eagle medium (DMEM, Gibco), and LUSC cell line HCC1588 was cultured in RPMI1-1,640 (Solarbio), which contained 10% heat-inactivated fetal bovine serum (FBS) and 1% penicillin/streptomycin. All cells are incubated in an incubator containing 5% CO₂ at 37°C.

2.14 16HBE Cells and HCC1588 Cells RNA Extraction and Real-Time PCR

The expression of lncRNAs was detected by real-time PCR in different cell lines. Using Trizol method to extract total RNA from different cells. Primers for AC138035.1 (forward, 5′-TCA​GAT​GAG​CAG​CAG​CGT​TAG​ATT​C-3′, reverse, 5′- GAG​CGT​GAT​GGG​AGA​AAG​TGA​CAG-3′), Primers for HORMAD2-AS1 (forward, 5′-GCA​GGA​GAA​CCA​CAG​TGA​CAA​CC-3′, reverse, 5′- TGC​TGA​CAG​TAG​TGC​TTG​CCT​ATT-3′), Primers for AC243919.2 (forward, 5′-CTC​GCC​AAC​TGG​TCC​TTC​ATC​TTC-3′, reverse, 5′- CTC​CTT​CAT​GCT​AAG​CCT​CCT​CTT​G-3′), Primers for AL122125.1 (forward, 5′-CAA​GCA​TGT​GGC​ACT​AGA​GGA​GAC-3′, reverse, 5′- CAG​AGA​TGG​AGG​CAG​AGG​TTG​AAT​G-3′) and Primers for GAPDH (forward, 5′-TGA​CTT​CAA​CAG​CGA​CAC​CCA-3′, reverse, 5′- CAC​CCT​GTT​GCT​GTA​GCC​AAA-3′) by real-time PCR and it was performed with three replicates of samples from three independent experiments.

2.15 Statistical Analysis

In this study, the R programming language (version 3.6.1, http://www.R-project.org) was used for analysis. The Kaplan-Meier method was used to generate the survival curves, and the log-rank test compared the differences between groups. The receiver operating characteristic (ROC) curves were used to estimate the predictive efficiency of m⁶A-LPMs for the 1-year overall survival (OS). Real-time PCR experimental data are expressed as the mean ± standard deviation (SD) of at least three independent experiments. All statistical analysis was performed using GraphPad Prism 8 (GraphPad Software, La Jolla, CA, USA). The student’s t-test (two-tailed) compares the two groups in real-time PCR. A p < 0.05 indicates statistical significance.

3 Results

3.1 Clinical Characteristics of LUSC Patients

We obtained prognostic information and clinicopathological features in a total of 498 LUSC patients after excluding patients with unknown survival outcomes. The clinicopathological characteristics of the patients are shown in Table 1.

TABLE 1

TABLE 1. Clinicopathological characteristics of LUSC patients in the TCGA cohort.

3.2 Prognosis of m⁶A-Related lncRNAs in LUSC Patients

We defined lncRNAs with expression levels that correlated with one or more of the 23 m⁶A-related genes as m⁶A-related lncRNAs by Pearson’s correlation analysis, according to the selection criteria R > |0.4|, p = 0.001. We identified that the expression of 351 lncRNAs was significantly correlated with m⁶A-related genes and constructed the co-expression network of m⁶A regulatory genes and m⁶A-related lncRNA (Figure 2A and Supplementary Tables S1, S2). Subsequently, we combined the expression levels of 351 m⁶A-related lncRNAs with the survival period and status of 498 LUSC patients to identify lncRNAs with prognostic value. Univariate Cox regression analysis was used to analyze prognostic information and expression of m⁶A-related lncRNAs in the TCGA database (p < 0.05). We finally identified the expression of 6 m⁶A-related lncRNAs: AC138035.1, AP001469.3, AC243919.2, PRC1−AS1, AL122125.1, and HORMAD2−AS1 was significantly impact on the OS of LUSC patients. (Figure 2B and Supplementary Table S3). Meanwhile, we analyzed the expression of prognostic m⁶A-related lncRNAs in LUSC tissues and normal tissues. As shown in Figure 2C, an increased expression of AC138035.1, AP001469.3, AC243919.2, PRC1−AS1, and AL122125.15 were observed in LUSC tissues compared to normal tissues.

FIGURE 2

FIGURE 2. Identification of m⁶A-related lncRNAs with prognostic value in LUSC patients. (A) Co-expression network of m⁶A-related genes and m⁶A-related lncRNAs. (B) m⁶A-associated lncRNAs with prognostic values were validated by univariate Cox analysis. (C) Expression of 6 prognostic m6A-related lncRNAs in tumor and normal tissues. LUSC, lung squamous cell carcinoma; m⁶A, N⁶-methyladenosine; lncRNAs, long non-coding RNAs; Statistical methods: Wilcoxon test **p < 0.01 and ***p < 0.001.

3.3 Consensus Clustering of m⁶A-Related lncRNAs in LUSC Patients

To further estimate the clinical relevance of m⁶A-related lncRNAs with prognostic value in LUSC patients, we used the “ConsensusClusterPlus” R package to divide LUSC tissue samples based on the expression of prognostic lncRNAs into two subgroups. In consensus clustering analysis, a k represented the cluster count, and k = 2-9 in the cumulative distribution function (CDF) (Figure 3A). The stability of consensus clustering analysis results should meet three conditions: the curve area does not increase significantly under the CDF; there was a close correlation among different clusters; the number of samples in a particular cluster cannot be too small. Figure 3B showed that the CDF curves of the consensus matrix did not significantly increase when k = 2. In addition, we selected k = 2 as the consensus matrix k value to divide the LUSC cohorts into two independent subgroups, cluster 1 (n = 418) and cluster 2 (n = 77); there was a significant correlation within the same cluster and a very low correlation between different clusters (Figure 3C). Heat maps showed the expression of 6 prognostic m⁶A-related lncRNAs in different clusters and their relationship with clinicopathological features (Figure 3D). We found no statistical significance among the two clusters in OS of LUSC patients (Figure 3E).

FIGURE 3

FIGURE 3. Clinicopathological characteristics and survival analysis of cluster 1/2 subgroups of LUSC. (A,B) Consensus clustering cumulative distribution function (CDF) for k = 2 to 9. (C) Consensus clustering matrix for k = 2. (D) Heatmap of the correlation between the expression levels of 6 m⁶A-related lncRNAs with prognostic value and clinicopathological features in the TCGA-LUSC dataset (red indicates high expression, blue indicates low expression). (E) Kaplan-Meier curves of overall survival for two clusters in LUSC. LUSC, lung squamous cell carcinoma; m⁶A, N⁶-methyladenosine; lncRNAs, long non-coding RNAs; TCGA, the cancer genome atlas; Statistical methods: Chi-square test, log-rank test.

Because previous studies revealed that multiple lncRNAs indirectly regulated PD-L1 expression to impact the survival of cancer patients, immune checkpoint inhibitors (ICIs) targeting programmed death-1/programmed death ligand 1 (PD-1/PD-L1) had been integrated into standard-of-care regimens for patients with advanced LUSC (Tian P. et al., 2021; Tian Y. et al., 2021; Mu et al., 2021; and; Yuan et al., 2021). Some studies also reported that the effect of first-line treatment with the PD-1 inhibitor pembrolizumab was superior to platinum-doublet chemotherapy in patients with non-small-cell lung cancer and high PD-L1 expression (Aguilar et al., 2019). Therefore, we wanted to estimated the expression of PD-L1 among different clusters and which cluster would more benefit from pembrolizumab, a PD-1 inhibitor. As shown in Figure 4, we detected that the expression of PD-L1 in LUSC tissue samples was lower than normal tissue samples (p < 0.001, Figure 4A), as well as the expression of PD-L1 in cluster 2 was also lower than cluster 1 (p < 0.05, Figure 4B). However, the expression of PD-L1 was not significantly correlated with the expression of lncRNA in LUSC (Figure 4C).

FIGURE 4

FIGURE 4. The expression level of PDL-1 in LUSC and the relationship between PD-L1 with prognostic-related lncRNAs. (A) The expression of PDL-1 in LUSC tissue samples. (B) The expression of PDL-1 in two clusters. (C) The relationship between PDL-1 and prognostic-related lncRNAs. LUSC, lung squamous cell carcinoma; lncRNAs, long non-coding RNAs; PD-L1, programmed cell death 1 ligand 1. Statistical methods: Wilcoxon test. *p < 0.05, ***p < 0.001.

3.4 Immune Score, Stromal Score, and TME in Different Clusters

Malignant solid tumor tissue is composed of tumor-associated epithelial cells, immune cells, stromal cells, vascular cells, and tumor cells (Yoshihara et al., 2013). Stromal cells play an essential role in tumor progression and drug resistance, and infiltrating immune cells can regulate the tumor microenvironment (TME) to influence tumor growth, invasion, and metastasis. (Hanahan and Weinberg, 2011; Mlecnik et al., 2011; Straussman et al., 2012,). The stromal score was planned to capture the presence of stroma in tumor tissues, and the immune score was designed to represent the infiltration of immune cells in tumor tissues. The estimated score was inferred tumor purity (Yoshihara et al., 2013). We use the ESTIMATE algorithm to calculate the immune, stromal, and ESTIMATE scores of each patient with LUSC patients. As is shown in Figures 5A–C, the ESTIMATE score, immune score, and stromal score of cluster 2 were lower than cluster 1. Additionally, we also assessed the expression of 22 immune cells in the two clusters of LUSC using the CIBERSORT algorithm. In cluster 1, the numbers of CD4 memory resting T cells, activated M2 macrophages, neutrophils, regulatory T cells (Tregs), and resting NK cells were higher than cluster 2; however, the numbers of naïve B cells, follicular helper T cells, and activated NK cells were higher in cluster 2 compared with cluster 1 (Figure 5D).

FIGURE 5

FIGURE 5. Immune function and immune cell infiltration in the two clusters of LUSC. (A) ESTIMATE score in the two clusters. (B) Immune score in the two clusters. (C) Stroma score in the two clusters. (D) The infiltration level of 22 immune cells in the two clusters. LUSC, lung squamous cell carcinoma. Statistical methods: Wilcoxon test.

3.5 GSEA Analysis Between the Two Clusters of LUSC

We performed GSEA to elucidate the potential regulatory mechanisms leading to the difference in TME between the two clusters. The enrichment analysis revealed that the malignant tumors correlated pathways, such as focal adhesion kinase (FAK), extracellular matrix (ECM) receptor interaction, and cytokine receptors, were mainly enriched in cluster 1 (Figure 6). Therefore, the activation of the FAK pathway and cytokine receptors may be the main reason for the difference in TME in LUSC cluster 1/2.

FIGURE 6

FIGURE 6. Gene enrichment analysis (GSEA) among different two clusters. (A–C) GSEA shows that focal adhesion kinase (A), extracellular matrix receptor interaction (B), and cytokine receptor (C) are enriched in cluster 1 of LUSC. LUSC, lung squamous cell carcinoma.

3.6 Identification and Evaluation Prognosis of 4 m⁶A-LPMs

First, univariate Cox regression analysis showed that the differential expression of the six m⁶A-related lncRNAs were closely associated with OS in LUSC patients. Next, the six lncRNAs were performed LASSO regression analysis, and we finally obtained four lncRNAs (AC138035.1, HORMAD2-AS1, AC243919.2 and AL122125.1) and coefficients of each of the lncRNAs to identify the most potent prognostic m⁶A-LPMs (Figures 7A,B). We randomly divided 495 LUSC patients into the training group (248) and the test group (247) in a 1:1 ratio. The clinical characteristics of LUSC patients in the training and test groups were shown in Table 2. Based on the expression and regression coefficients of the 4 m⁶A-LPMs in the training and test groups, we constructed a risk model and calculated the risk score of each patient. Subsequently, all patients were divided into the high-risk group and the low-risk group according to the median risk score. The risk score curve and survival status and the expression heatmap of 4 m⁶A-LPMs in different groups were shown in Figures 7C,D. These results revealed that the survival period of patients gradually decreased, and the number of deaths gradually increased with the increase of the risk score in the training group and the test group (green dots represented live patients, red dots represented dead patients). The heatmap results showed that the HORMAD2-AS1 was highly expressed in the high-risk group compared with the low-risk group, whereas the expression of AC138035.1, AC243919.2, and AL122125.1 was downregulated in the high-risk group compared with the low-risk group (Figures 7C,D). In the LUSC training and test groups, the OS of the low-risk group was significantly longer than the high-risk group (Figures 7E,F). ROC curves was used to assess the prognostic accuracy of this risk model. The area under the curve (AUC) values for 1-year of this risk model in the training and the test groups were 0.617 and 0.572, respectively (Figures 7G,H). These results indicate that this risk model can better predict the prognosis of LUSC patients.

FIGURE 7

FIGURE 7. Verifying the risk model. (A,B) LASSO Cox regression analysis of 6 m⁶A-related lncRNAs with prognostic value. (C) Survival status and expressive heatmap of 4 m⁶A-LPMs in the training group (red indicates high expression, green indicates low expression). (D) Survival status and expressive heatmap of 4 m⁶A-LPMs in the test group (red indicates high expression, green indicates low expression). (E,F) Kaplan-Meier curve of 4 m⁶A-LPMs in the training group (E) and test group (F). (G,H) time-ROC curves of 4 m⁶A-LPMs in the training group (G) and test group (H). LASSO, the least absolute shrinkage and selection operator; ROC, the receiver operating characteristic curve; m⁶A-LPMs, m⁶A-related lncRNA prognostic markers. Statistical methods: Chi-square test, log-rank test.

TABLE 2

TABLE 2. Clinicopathological characteristics of LUSC patients in training and test group.

3.7 Correlation Analysis of m⁶A-LPMs With Clinicopathological Features and Immune Score

We further evaluated the relationship between the risk score, clinicopathological features, and immune score. As shown in Figure 8A, the risk score in the females is higher than in males; we also found that compared with the high immune score group and cluster 1, the low immune score group and cluster 2 had a lower risk score (Figures 8B,C). As well as, patients with LUSC stage I-II had lower risk scores than those with stage III-IV (Figure 8D). These results showed that the risk score was correlated with TME and tumor malignant stage.

FIGURE 8

FIGURE 8. Correlation analysis of risk score with clinicopathological characteristics and immune scores. (A) Boxplots of risk score with different genders. (B) Boxplots of risk scores and different immune scores. (C) Boxplots of risk scores and different clusters. (D) Boxplots of risk scores and different Stages. Statistical methods: Chi-square test. ***p < 0.001.

To better evaluate the prognostic value of m⁶A-LPMs, we further assessed whether they retained the ability to predict OS of patients in different subgroups by stratified analysis. Our analysis found that high-risk patients had poor OS than low-risk patients in all age groups (Figures 9A,B). Gender influenced the ability of m⁶A-LPMs to predict OS in LUSC patients. As is shown in Figures 9C,D, males in the low-risk group had better OS than those in the high-risk group, whereas risk scores did not affect OS in females. In addition, m⁶A-LPMs also predicts OS in LUSC patients with stage I-II but does not predict OS in LUSC patients with stage III-IV patients (Figures 9E,F). Overall, these results indicate that m⁶A-LPMs may be potential predictors for patients with LUSC.

FIGURE 9

FIGURE 9. Stratified analysis of survival based on clinicopathological features. (A–F) The predictive ability of m⁶A-LPMs in multiple subgroups of LUSC patients, including the different ages: >65 (A) or≤65 (B), different sexes: female (C) or male (D), different stages: I-II (E) or III-IV (F). LUSC, lung squamous cell carcinoma; m⁶A-LPMs, m⁶A-related lncRNA prognostic markers. Statistical methods: Log-rank test.

3.8 Independent Prognostic Factors for Patients With LUSC

We also performed univariate and multivariate cox analyses to determine whether the risk score are independent prognostic factors for LUSC patients. As is shown in Figure 10A, univariate Cox analysis showed that the risk score, age, and stage as optimal prognostic factors for LUSC patients in the entire TCGA dataset (the hazard ratio (HR) of risk score was 1.690 and 95% confidence interval (CI) was 1.168–2.445, p = 0.005; the HR of age was 1.017 and 95% CI was 1.000–1.034, p = 0.044; and the HR of the stage was 1.256 and 95% CI was 1.064–1.482, p = 0.007, Supplementary Table S4). Moreover, after adjusting for age, sex, and stage, we found that risk scores was independent prognostic factors for LUSC patients in the TCGA entire set (Figure 10B, Supplementary Table S5).

FIGURE 10

FIGURE 10. Univariate and multivariate cox analyses of the risk score. (A–B) Univariate (A) and multivariate (B) Cox analyses of the risk score in the TCGA entire dataset. TCGA, the cancer genome atlas.

3.9 Correlation Analysis of m⁶A-LPMS in LUSC With Immune Cells

We also analyzed the relationship between the risk score and immune cells to evaluate the impact of m⁶A-LPMs on the TME of LUSC. We found that the infiltration level of some immune cells (resting dendritic cells, eosinophils, activated M2 macrophages, neutrophils, activated CD4 memory-activated cells, and gamma delta T cells) in LUSC is positively correlated with the risk score (Figures 11A–F). The infiltration level of other immune cells (naïve B cells, non-activated M0 macrophages, activated NK cells, and follicular helper T cells) is negatively correlated with the risk score (Figures 11G–J). These results suggest that m⁶A-LPMS might influence cancer progression by regulating the level of immune cells in LUSC.

FIGURE 11

FIGURE 11. The correlation analysis between the risk score and immune cells. (A) resting dendritic cells, (B) eosinophils, (C) M2 macrophages, (D) neutrophils, (E) activated CD4 memory activated cells, (F) gamma delta T cells, (G) naïve B cells, (H) M0 macrophages, (I) activated NK cells, (J) follicular helper T cells.

3.10 Verifying the Expression of m⁶A-LPMs in 16HBE Cell and HCC1588 Cell

Real-time PCR was used further to verify the expression of m⁶A-LPMs in different cells. We found that the expression of HORMAD2-AS1 was not statistically significant between the two groups; however, the trend of expression was consistent with our predicted results (Figure 12A). The expression of AC138035.1, AC243919.2, and AL122125.1 in the HCC1588 cell was statistically higher than in the 16HBE cell (Figures 12B–D). We will collect more clinical specimens and cell lines to verify our results.

FIGURE 12

FIGURE 12. The expression of 4 m⁶A-LPMs in 16HBE cell and HCC1588 cell. (A–D) The expression of HORMAD2-AS1 (A), AC138035.1 (B), AC243919.2 (C), and AL122125.1 (D) in 16HBE cells and HCC1588 cells were verified by real-time PCR. m⁶A-LPMs, m⁶A-related lncRNA prognostic markers. Statistical methods: t-test, *p < 0.05, **p < 0.01.

4 Discussion

LUSC has a high rate of metastasis and recurrence. The therapy for LUSC is limited, so that the prognosis of LUSC patients is still poor (Qian et al., 2016). Different subtypes of lung cancer patients have different clinical characteristics and clinical outcomes; therefore, many researchers tend to study non-coding RNA signatures to predict the prognosis and the response to immunotherapy of NSCLC patients (Song et al., 2019; Tian et al., 2020). Some studies have found that m⁶A modification plays an essential role in tumor progression. In several tumors, m⁶A modulators can promote tumor growth and metastasis by regulating the expression of the corresponding lncRNAs. Zheng et al. (2019) found that m⁶A methylation in nasopharyngeal carcinoma reduces the rate of RNA degradation and increases the stability of methylated FAM225A transcripts; it combines with miR-590–3p and miR-1275 to regulate the expression of ITGB3 and activate the FAK/PI3K/AKT pathway to promote the proliferation and invasion of nasopharyngeal carcinoma. Wu et al. (2019) also found that m⁶A-induced lncRNA RP11-triggered epithelial-mesenchymal transition (EMT) by the upregulation of ZEB1 promotes invasion of colorectal cancer cells. Recently, a few studies found that m⁶A regulation of lncRNAs can promote the progression and resistance of NSCLC. For instance, Jin et al. (2019) found that m⁶A mRNA methylation initiated by the m⁶A transferase METL3 promotes YAP mRNA translation by recruiting YTHDF1/3 and iodine trifluoride to the translation initiation complex and increases expression of YAP mRNA by regulating the MALAT1-miR-1914-3p-YAP axis stability. The increase in YAP mRNA expression and activity induces resistance and metastasis of NSCLC. These studies have confirmed that the modification of lncRNAs by m⁶A may affect the occurrence and progression of various cancers, including NSCLC, and that lncRNAs can also target m⁶A regulators to promote tumor progression. However, there have only been a few studies on the mechanism by which m⁶A regulates lncRNAs in LUSC to promote tumor progression. Therefore, more research focused on the interaction between m⁶A modification and lncRNAs is needed to identify potential new targets for the prognosis and treatment of LUSC.

Our study confirmed the expression network of m⁶A-related lncRNAs in LUSC, its prognostic value, its impact on the TME, and its relationship with clinicopathological characteristics. We identified 6 m⁶A-related prognostic lncRNAs from 481 LUSC patients. The expression of 5 lncRNAs (AC138035.1, AP001469.3, AC243919.2, PRC1−AS1, and AL122125.1) in tumor tissues were higher than those in normal adjacent tissues, but the expression levels of lncRNAs HORMAD2−AS1 were lower than those in normal adjacent tissues. Among them, the high expression of AC138035.1 is closely related to the poor prognosis of ovarian cancer patients (Lin H. et al., 2021). Zhao et al. (2021). found that a natural plant peptide aldehyde inhibitor MG-132 (carbobenzoxy-Leu-Leu-leucinal) to inhibit human pterygium fibroblasts (HPFs) proliferate and increase apoptosis by regulating the expression of lncRNAs such as AL122125.1. In addition, we also discovered other lncRNAs for the first time. According to the consensus cluster analysis of the prognostic-related lncRNAs, we identified two subtypes of LUSC, namely cluster 1 and cluster 2. It was found that the expression of PD-L1 in LUSC tumor tissues was lower than that in normal tissues, and the expression of PD-L1 in cluster 2 in LUSC was lower than that in cluster 1. PD-L1 protein expression level has been used as a biomarker to predict which patients are more likely to respond to immunotherapy; and studies have found that the expression level of PD-L1 in NSCLC patients is 24%–60%, which may result from differences in patients’ treatments methods, detection antibodies or demographics (Yu et al., 2016). We found that the immune score and matrix score of cluster 2 in LUSC were lower than cluster 1. However, the immune and matrix scores were not significantly correlated with clinicopathological characteristics and prognosis. Our findings are similar to some previous studies. For example, Qu et al. (2020) found that the matrix score had no significant correlation with clinical characteristics and prognosis for LUSC cases. In cluster 1 of LUSC, there were higher infiltration levels of resting CD4 memory T cells, activated M2 macrophages, neutrophils, regulatory T cells (Tregs), and resting NK cells than in cluster 2. However, the infiltration levels of naïve B cells, follicular helper T cells, and activated NK cells were higher in cluster 2 than in cluster 1. The results of the GSEA showed that focal adhesion, ECM receptor interaction, and cytokine receptors related to the malignant characteristics of tumors were mainly found in cluster 1. Focal adhesion is a non-receptor tyrosine kinase involved in cancer cell growth, proliferation, survival, migration, angiogenesis, invasion, and EMT, essential for maintaining cancer stem cells and macrophages (Soria et al., 2016). m⁶A modification directly controls RNA metabolism, including mRNA processing, output, translation initiation and maintenance stability, and the biogenesis of lncRNA, which may be involved in the regulation of cytokines, immune response, and the self-renewal ability of cancer stem cells and play a key role in cell proliferation and apoptosis (Chang et al., 2019). Therefore, we suggest that the differences in TME between different subgroups of LUSC are regulated by the FAK pathway and cytokine receptors.

Subsequently, we evaluated and verified the prognostic value of 4 m⁶A-related-LPMs in LUSC patients. Compared with bronchial epithelial cells, the expression of AC138035.1, AC243919.2, and AL122125.1 was increased in HCC1588 cells by real-time PCR, while the expression of HORMAD2−AS1 was not statistically significant. The limitation of cells number, lacking verified tests by tissues, and tumor specificity may be the main reasons. More clinical specimens and cell lines should be used to enhance the research credibility by further experimental proof as the study of these lncRNAs remains unknown. We divide LUSC patients into high- and low-risk groups by risk score. In the LUSC training group and test group, the OS of low-risk patients was longer than that of high-risk patients. In addition, the AUC values for 1 year OS were 0.617 and 0.572 in the training and test groups, respectively. These results indicate that this risk model can better predict the prognosis of LUSC patients. Previous studies have shown that the root mean square difference between the estimated and true metrics is significant in the case of a small sample size, which may affect the real performance of the ROC curve and lead to a low AUC value (Hanczar et al., 2010; Berrar and Flach, 2012) Since our data come from the TCGA database, the possible reason for the low AUC value is related to the small sample size in our different groups. The risk score of the male, low immune score group, and cluster 2 was lower than those of the female, high immune score group, and cluster 1. The stratified analysis revealed that the prognostic characteristics of m⁶A-LPMs retained the ability to predict the OS rate of different ages, male patients, and early tumor patients. We also confirmed that the risk score could be an independent prognostic factor for LUSC patients by univariate and multivariate Cox regression analysis. Previous studies have confirmed that the increased expression level of AC138035.1 is closely related to the poor prognosis of patients with ovarian cancer (Lin N. et al., 2021). These studies show that m⁶A-LPMs may be effective in predicting the prognosis of LUSC patients.

The TME plays a vital role in the development of tumors and affects the treatment and prognosis of patients (Fridman et al., 2012; Galon et al., 2014; Hui and Chen, 2015). However, the mechanism of the interaction between immune infiltration and tumor in LUSC remains unclear, and the understanding of m⁶A-modified lncRNAs in the regulation of TME is still limited. Previous studies have confirmed a large number of neutrophil infiltration in the high-grade invasive histological subtypes of LUAD, which is closely related to the poor prognosis of LUAD patients (Kurebayashi et al., 2016). Neutrophils exist in the tumor immune microenvironment of the tumor and can boost tumor progression by promoting cell growth, angiogenesis, metastasis, and immune evasion (Mollaoglu et al., 2018). Macrophages are classified into M1 and M2 macrophages according to their functions. Studies have found that the M1 macrophages exert pro-inflammatory effects and anti-tumor immunity, while M2 macrophages participate in anti-inflammatory response and promote tumor progression effects (Chanmee et al., 2014). Some studies have found that increased eosinophils in tumor tissues can improve the prognosis of patients, but eosinophils are related to poor prognosis in Hodgkin lymphoma (Davis and Rothenberg, 2014). Our study discovered that the infiltration level of some immune cells (resting dendritic cells, eosinophils, activated M2 macrophages, neutrophils, activated CD4 memory T cells, and gamma delta T cells) is positively correlated with the risk score and poor prognosis of LUSC patients. At the same time, we also found that other immune cells (naïve B cells, follicular helper T cells, activated NK cells, non-activated M0 macrophages) infiltration level is negatively correlated with a risk score. Chen et al., 2020b have found that naïve-like B cells not only inhibit the proliferation of NSCLC cells but are also associated with a good prognosis of NSCLC; the expression of NK cells in tumor tissues is positively correlated with the prognosis of cancer patients. NK cells can enhance the ability of antibodies and T cells to recognize tumors, and activated NK cells can eliminate the cancer cells in the circulating blood and overcome tumor resistance (Shimasaki and Jain, 2020). Infiltrating follicular helper T cells play a protective role in breast and colon cancer, and its increased expression level in tumors indicates that it has a better prognosis for patients (Bindea et al., 2013; Gu-Trantien et al., 2013). The above results indicated that m⁶A modified lncRNAs might affect patients with LUSC by regulating TME. Therefore, it is urgent to study further the regulatory role of m⁶A-related lncRNAs in the TME of LUSC.

However, there are still some defects in our research. Firstly, our data mainly comes from the lung cancer cohort in the TCGA database, and more external validation is needed in the future to assess whether it can be used in clinical patients. Secondly, because of the limited number of cells and tumor tissues, there might be a slight deviation in the validation of m⁶A-LPMs; moreover, the molecular mechanism and the role of these lncRNAs in tumor growth, metastasis, and immune infiltration need to be further studied. Therefore, in subsequent studies, we will further verify our results by experiments in vivo and in vitro; and investigate the impact of m⁶A and its prognostic-related lncRNAs on the TME to improve the efficacy of immunotherapy LUSC patients.

5 Conclusion

In conclusion, we found that 4 m⁶A-LPMs (AC138035.1, AC243919.2, HORMAD2-AS1 and AL122125.1) are closely associated with LUSC patients’ prognoses. We also divided LUSC patients into high-risk and low-risk groups according to the risk score of the m⁶A-related LPMs. This risk score was found to be closely related to the immune score, clinicopathological characteristics, and the level of immune cell infiltration of LUSC patients, and it was also found to be an independent prognostic indicator of LUSC patients. These biomarkers are beneficial to the treatment and prognosis of LUSC patients. Further study of their potential regulatory mechanisms may identify novel targets for the prognosis and immunotherapy of LUSC patients.

Data Availability Statement

The datasets generated and analyzed during this study are available in the TCGA database (https://portal.gdc.cancer.gov). The original contributions presented in the study are included in the Article/Supplementary Material; further inquiries can be directed to the corresponding authors.

Author Contributions

SL: designed the study, WZ: analyzed, interpreted the data, and drafted the manuscript. QZ, ZX collected the data. LC, TX, and XC: assisted in manuscript revision and figure organization. All authors have read and approved the manuscript.

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.

Supplementary Material

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

References

Aguilar, E. J., Ricciuti, B., Gainor, J. F., Kehl, K. L., Kravets, S., Dahlberg, S., et al. (2019). Outcomes to First-Line Pembrolizumab in Patients with Non-small-cell Lung Cancer and Very High PD-L1 Expression. Ann. Oncol. 30, 1653–1659. doi:10.1093/annonc/mdz288

PubMed Abstract | CrossRef Full Text | Google Scholar

Berrar, D., and Flach, P. (2012). Caveats and Pitfalls of ROC Analysis in Clinical Microarray Research (And How to Avoid Them). Briefings Bioinforma. 13, 83–97. doi:10.1093/bib/bbr008

PubMed Abstract | CrossRef Full Text | Google Scholar

Bindea, G., Mlecnik, B., Tosolini, M., Kirilovsky, A., Waldner, M., Obenauf, A. C., et al. (2013). Spatiotemporal Dynamics of Intratumoral Immune Cells Reveal the Immune Landscape in Human Cancer. Immunity 39, 782–795. doi:10.1016/j.immuni.2013.10.003

PubMed Abstract | CrossRef Full Text | Google Scholar

Brown, J. A., Bulkley, D., Wang, J., Valenstein, M. L., Yario, T. A., Steitz, T. A., et al. (2014). Structural Insights into the Stabilization of MALAT1 Noncoding RNA by a Bipartite Triple Helix. Nat. Struct. Mol. Biol. 21, 633–640. doi:10.1038/nsmb.2844

PubMed Abstract | CrossRef Full Text | Google Scholar

Brown, J. A., Kinzig, C. G., DeGregorio, S. J., and Steitz, J. A. (2016). Methyltransferase-like Protein 16 Binds the 3′-terminal Triple Helix of MALAT1 Long Noncoding RNA. Proc. Natl. Acad. Sci. U.S.A. 113, 14013–14018. doi:10.1073/pnas.1614759113

PubMed Abstract | CrossRef Full Text | Google Scholar

Cabili, M. N., Trapnell, C., Goff, L., Koziol, M., Tazon-Vega, B., Regev, A., et al. (2011). Integrative Annotation of Human Large Intergenic Noncoding RNAs Reveals Global Properties and Specific Subclasses. Genes Dev. 25, 1915–1927. doi:10.1101/gad.17446611

PubMed Abstract | CrossRef Full Text | Google Scholar

Chang, G., Leu, J.-S., Ma, L., Xie, K., and Huang, S. (2019). Methylation of RNA N6-Methyladenosine in Modulation of Cytokine Responses and Tumorigenesis. Cytokine 118, 35–41. doi:10.1016/j.cyto.2018.06.018

PubMed Abstract | CrossRef Full Text | Google Scholar

Chanmee, T., Ontong, P., Konno, K., and Itano, N. (2014). Tumor-associated Macrophages as Major Players in the Tumor Microenvironment. Cancers 6 (3), 1670–1690. doi:10.3390/cancers6031670

PubMed Abstract | CrossRef Full Text | Google Scholar

Chen, J., Tan, Y., Sun, F., Hou, L., Zhang, C., Ge, T., et al. (2020a). Single-cell Transcriptome and Antigen-Immunoglobin Analysis Reveals the Diversity of B Cells in Non-small Cell Lung Cancer. Genome Biol. 21, 152. doi:10.1186/s13059-020-02064-6

PubMed Abstract | CrossRef Full Text | Google Scholar

Chen, Y., Lin, Y., Shu, Y., He, J., and Gao, W. (2020b). Interaction between N6-Methyladenosine (m6A) Modification and Noncoding RNAs in Cancer. Mol. Cancer 19, 94. doi:10.1186/s12943-020-01207-4

PubMed Abstract | CrossRef Full Text | Google Scholar

Dai, F., Wu, Y., Lu, Y., An, C., Zheng, X., Dai, L., et al. (2020). Crosstalk between RNA m6A Modification and Non-coding RNA Contributes to Cancer Growth and Progression. Mol. Ther. - Nucleic Acids 22, 62–71. doi:10.1016/j.omtn.2020.08.004

PubMed Abstract | CrossRef Full Text | Google Scholar

Davis, B. P., and Rothenberg, M. E. (2014). Eosinophils and Cancer. Cancer Immunol. Res. 2, 1–8. doi:10.1158/2326-6066.CIR-13-0196

PubMed Abstract | CrossRef Full Text | Google Scholar

Fan, F. S., Yang, C.-F., and Chang, C.-L. (2019). Nivolumab Plus Carboplatin and Paclitaxel as the First-Line Therapy for Advanced Squamous Cell Carcinoma of the Lung with Strong Programmed Death-Ligand 1 Expression: A Case Report. Cureus 11, e5881. doi:10.7759/cureus.5881

PubMed Abstract | CrossRef Full Text | Google Scholar

Fazi, F., and Fatica, A. (2019). Interplay between N6-Methyladenosine (m6A) and Non-coding RNAs in Cell Development and Cancer. Front. Cell Dev. Biol. 7, 116. doi:10.3389/fcell.2019.00116

PubMed Abstract | CrossRef Full Text | Google Scholar

Forde, P. M., Kelly, R. J., and Brahmer, J. R. (2014). New Strategies in Lung Cancer: Translating Immunotherapy into Clinical Practice. Clin. Cancer Res. 20, 1067–1073. doi:10.1158/1078-0432.ccr-13-0731

PubMed Abstract | CrossRef Full Text | Google Scholar

Fridman, W. H., Pagès, F., Sautès-Fridman, C., and Galon, J. (2012). The Immune Contexture in Human Tumours: Impact on Clinical Outcome. Nat. Rev. Cancer 12, 298–306. doi:10.1038/nrc3245

PubMed Abstract | CrossRef Full Text | Google Scholar

Galon, J., Mlecnik, B., Bindea, G., Angell, H. K., Berger, A., Lagorce, C., et al. (2014). Towards the Introduction of the 'Immunoscore' in the Classification of Malignant Tumours. J. Pathol. 232, 199–209. doi:10.1002/path.4287

PubMed Abstract | CrossRef Full Text | Google Scholar

Gu-Trantien, C., Loi, S., Garaud, S., Equeter, C., Libin, M., de Wind, A., et al. (2013). CD4+ Follicular Helper T Cell Infiltration Predicts Breast Cancer Survival. J. Clin. Invest. 123, 2873–2892. doi:10.1172/JCI67428

PubMed Abstract | CrossRef Full Text | Google Scholar

Hanahan, D., and Weinberg, R. A. (2011). Hallmarks of Cancer: the Next Generation. Cell 144, 646–674. doi:10.1016/j.cell.2011.02.013

PubMed Abstract | CrossRef Full Text | Google Scholar

Hanczar, B., Hua, J., Sima, C., Weinstein, J., Bittner, M., and Dougherty, E. R. (2010). Small-sample Precision of ROC-Related Estimates. Bioinformatics 26, 822–830. doi:10.1093/bioinformatics/btq037

PubMed Abstract | CrossRef Full Text | Google Scholar

He, L., Li, H., Wu, A., Peng, Y., Shu, G., and Yin, G. (2019). Functions of N6-Methyladenosine and its Role in Cancer. Mol. Cancer 18, 176. doi:10.1186/s12943-019-1109-9

PubMed Abstract | CrossRef Full Text | Google Scholar

Huang, H., Weng, H., and Chen, J. (2020). m6A Modification in Coding and Non-coding RNAs: Roles and Therapeutic Implications in Cancer. Cancer Cell 37, 270–288. doi:10.1016/j.ccell.2020.02.004

PubMed Abstract | CrossRef Full Text | Google Scholar

Huang, Y.-H., Zhu, C., Kondo, Y., Anderson, A. C., Gandhi, A., Russell, A., et al. (2015). CEACAM1 Regulates TIM-3-Mediated Tolerance and Exhaustion. Nature 517, 386–390. doi:10.1038/nature13848

PubMed Abstract | CrossRef Full Text | Google Scholar

Hui, L., and Chen, Y. (2015). Tumor Microenvironment: Sanctuary of the Devil. Cancer Lett. 368, 7–13. doi:10.1016/j.canlet.2015.07.039

PubMed Abstract | CrossRef Full Text | Google Scholar

Jin, D., Guo, J., Wu, Y., Du, J., Yang, L., Wang, X., et al. (2019). m6A mRNA Methylation Initiated by METTL3 Directly Promotes YAP Translation and Increases YAP Activity by Regulating the MALAT1-miR-1914-3p-YAP axis to Induce NSCLC Drug Resistance and Metastasis. J. Hematol. Oncol. 12, 135. doi:10.1186/s13045-019-0830-6

PubMed Abstract | CrossRef Full Text | Google Scholar

Kunimasa, K., and Goto, T. (2020). Immunosurveillance and Immunoediting of Lung Cancer: Current Perspectives and Challenges. Ijms 21, 597. doi:10.3390/ijms21020597

PubMed Abstract | CrossRef Full Text | Google Scholar

Kurebayashi, Y., Emoto, K., Hayashi, Y., Kamiyama, I., Ohtsuka, T., Asamura, H., et al. (2016). Comprehensive Immune Profiling of Lung Adenocarcinomas Reveals Four Immunosubtypes with Plasma Cell Subtype a Negative Indicator. Cancer Immunol. Res. 4, 234–247. doi:10.1158/2326-6066.CIR-15-0214

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, J., Wang, J., Chen, Y., Li, S., Jin, M., Wang, H., et al. (2016). LncRNA MALAT1 Exerts Oncogenic Functions in Lung Adenocarcinoma by Targeting miR-204. Am. J. Cancer Res. 6, 1099–1107. PMID: 27294002; PMCID: PMC4889723.

PubMed Abstract | Google Scholar

Li, S., Mei, Z., Hu, H. B., and Zhang, X. (2018). The lncRNA MALAT1 Contributes to Non‐small Cell Lung Cancer Development via Modulating miR‐124/STAT3 axis. J. Cell Physiol. 233, 6679–6688. doi:10.1002/jcp.26325

PubMed Abstract | CrossRef Full Text | Google Scholar

Lin, C.-P., Sung, Y.-C., Lo, C.-Y., and Yen, M.-H. (2020). Pathological Complete Response of Initially Inoperable Lung Squamous Cell Carcinoma Treated by Immunochemotherapy: A Case Report. Asian J. Surg. 43, 393–395. doi:10.1016/j.asjsur.2019.10.008

PubMed Abstract | CrossRef Full Text | Google Scholar

Lin, H., Weng, J., Mei, H., Zhuang, M., Xiao, X., Du, F., et al. (2021a). 5‐Lipoxygenase Promotes Epithelial-Mesenchymal Transition through the ERK Signaling Pathway in Gastric Cancer. J. Gastroenterology Hepatology 36, 455–466. doi:10.1111/jgh.15184

CrossRef Full Text | Google Scholar

Lin, N., Lin, J.-z., Tanaka, Y., Sun, P., and Zhou, X. (2021b). Identification and Validation of a Five-lncRNA Signature for Predicting Survival with Targeted Drug Candidates in Ovarian Cancer. Bioengineered 12, 3263–3274. doi:10.1080/21655979.2021.1946632

PubMed Abstract | CrossRef Full Text | Google Scholar

Linder, B., Grozhik, A. V., Olarerin-George, A. O., Meydan, C., Mason, C. E., and Jaffrey, S. R. (2015). Single-nucleotide-resolution Mapping of m6A and m6Am throughout the Transcriptome. Nat. Methods 12, 767–772. doi:10.1038/nmeth.3453

PubMed Abstract | CrossRef Full Text | Google Scholar

Liu, J., Ren, D., Du, Z., Wang, H., Zhang, H., and Jin, Y. (2018). m 6 A Demethylase FTO Facilitates Tumor Progression in Lung Squamous Cell Carcinoma by Regulating MZF1 Expression. Biochem. Biophysical Res. Commun. 502, 456–464. doi:10.1016/j.bbrc.2018.05.175

CrossRef Full Text | Google Scholar

Liu, Y., Guo, X., Zhao, M., Ao, H., Leng, X., Liu, M., et al. (2020). Contributions and Prognostic Values of M 6 A RNA Methylation Regulators in Non‐small‐cell Lung Cancer. J. Cell Physiol. 235, 6043–6057. doi:10.1002/jcp.29531

PubMed Abstract | CrossRef Full Text | Google Scholar

Ma, S., Chen, C., Ji, X., Liu, J., Zhou, Q., Wang, G., et al. (2019). The Interplay between m6A RNA Methylation and Noncoding RNA in Cancer. J. Hematol. Oncol. 12, 121. doi:10.1186/s13045-019-0805-7

PubMed Abstract | CrossRef Full Text | Google Scholar

Mlecnik, B., Tosolini, M., Kirilovsky, A., Berger, A., Bindea, G., Meatchi, T., et al. (2011). Histopathologic-based Prognostic Factors of Colorectal Cancers Are Associated with the State of the Local Immune Reaction. Jco 29 (6), 610–618. doi:10.1200/JCO.2010.30.5425

PubMed Abstract | CrossRef Full Text | Google Scholar

Mollaoglu, G., Jones, A., Wait, S. J., Mukhopadhyay, A., Jeong, S., Arya, R., et al. (2018). The Lineage-Defining Transcription Factors SOX2 and NKX2-1 Determine Lung Cancer Cell Fate and Shape the Tumor Immune Microenvironment. Immunity 49, 764–779. doi:10.1016/j.immuni.2018.09.020

PubMed Abstract | CrossRef Full Text | Google Scholar

Mu, L., Wang, Y., Su, H., Lin, Y., Sui, W., Yu, X., et al. (2021). HIF1A-AS2 Promotes the Proliferation and Metastasis of Gastric Cancer Cells through miR-429/pd-L1 Axis. Dig. Dis. Sci. 66, 4314–4325. doi:10.1007/s10620-020-06819-w

PubMed Abstract | CrossRef Full Text | Google Scholar

Paz-Ares, L., Luft, A., Vicente, D., Tafreshi, A., Gümüş, M., Mazières, J., et al. (2018). Pembrolizumab Plus Chemotherapy for Squamous Non-small-cell Lung Cancer. N. Engl. J. Med. 379, 2040–2051. doi:10.1056/NEJMoa1810865

PubMed Abstract | CrossRef Full Text | Google Scholar

Ping, X.-L., Sun, B.-F., Wang, L., Xiao, W., Yang, X., Wang, W.-J., et al. (2014). Mammalian WTAP Is a Regulatory Subunit of the RNA N6-Methyladenosine Methyltransferase. Cell Res. 24, 177–189. doi:10.1038/cr.2014.3

PubMed Abstract | CrossRef Full Text | Google Scholar

Qian, L., Lin, L., Du, Y., Hao, X., Zhao, Y., and Liu, X. (2016). MicroRNA-588 Suppresses Tumor Cell Migration and Invasion by Targeting GRN in Lung Squamous Cell Carcinoma. Mol. Med. Rep. 14, 3021–3028. doi:10.3892/mmr.2016.5643

PubMed Abstract | CrossRef Full Text | Google Scholar

Qu, Y., Cheng, B., Shao, N., Jia, Y., Song, Q., Tan, B., et al. (2020). Prognostic Value of Immune-Related Genes in the Tumor Microenvironment of Lung Adenocarcinoma and Lung Squamous Cell Carcinoma. Aging 12, 4757–4777. doi:10.18632/aging.102871

PubMed Abstract | CrossRef Full Text | Google Scholar

Ribas, A. (2015). Adaptive Immune Resistance: How Cancer Protects from Immune Attack. Cancer Discov. 5, 915–919. doi:10.1158/2159-8290.CD-15-0563

PubMed Abstract | CrossRef Full Text | Google Scholar

Shimasaki, N., Jain, A., and Campana, D. (2020). NK Cells for Cancer Immunotherapy. Nat. Rev. Drug Discov. 19, 200–218. doi:10.1038/s41573-019-0052-1

PubMed Abstract | CrossRef Full Text | Google Scholar

Song, Q., Shang, J., Yang, Z., Zhang, L., Zhang, C., Chen, J., et al. (2019). Identification of an Immune Signature Predicting Prognosis Risk of Patients in Lung Adenocarcinoma. J. Transl. Med. 17, 70. doi:10.1186/s12967-019-1824-4

PubMed Abstract | CrossRef Full Text | Google Scholar

Soria, J. C., Gan, H. K., Blagden, S. P., Plummer, R., Arkenau, H. T., Ranson, M., et al. (2016). A Phase I, Pharmacokinetic and Pharmacodynamic Study of GSK2256098, a Focal Adhesion Kinase Inhibitor, in Patients with Advanced Solid Tumors. Ann. Oncol. 27, 2268–2274. doi:10.1093/annonc/mdw427

PubMed Abstract | CrossRef Full Text | Google Scholar

Straussman, R., Morikawa, T., Shee, K., Barzily-Rokni, M., Qian, Z. R., Du, J., et al. (2012). Tumour Micro-environment Elicits Innate Resistance to RAF Inhibitors through HGF Secretion. Nature 487 (7408), 500–504. doi:10.1038/nature11183

PubMed Abstract | CrossRef Full Text | Google Scholar

Sung, H., Ferlay, J., Siegel, R. L., Laversanne, M., Soerjomataram, I., Jemal, A., et al. (2021). Global Cancer Statistics 2020: GLOBOCAN Estimates of Incidence and Mortality Worldwide for 36 Cancers in 185 Countries. CA A Cancer J. Clin. 71, 209–249. doi:10.3322/caac.21660

CrossRef Full Text | Google Scholar

Tang, C., Klukovich, R., Peng, H., Wang, Z., Yu, T., Zhang, Y., et al. (2018). ALKBH5-dependent m6A Demethylation Controls Splicing and Stability of Long 3′-UTR mRNAs in Male Germ Cells. Proc. Natl. Acad. Sci. U.S.A. 115, E325–E333. doi:10.1073/pnas.1717794115

PubMed Abstract | CrossRef Full Text | Google Scholar

Tian, P., Wei, J. X., Li, J., Ren, J. K., and Yang, J. J. (2021a). LncRNA SNHG1 Regulates Immune Escape of Renal Cell Carcinoma by Targeting miR‐129‐3p to Activate STAT3 and PD‐L1. Cell Biol. Int. 45, 1546–1560. doi:10.1002/cbin.11595

PubMed Abstract | CrossRef Full Text | Google Scholar

Tian, Y., Li, L., Lin, G., Wang, Y., Wang, L., Zhao, Q., et al. (2021b). lncRNA SNHG14 Promotes Oncogenesis and Immune Evasion in Diffuse Large-B-Cell Lymphoma by Sequestering miR-152-3p. Leukemia Lymphoma 62, 1574–1584. doi:10.1080/10428194.2021.1876866

PubMed Abstract | CrossRef Full Text | Google Scholar

Tian, Y., Yu, M., Sun, L., Liu, L., Wang, J., Hui, K., et al. (2020). Distinct Patterns of mRNA and lncRNA Expression Differences between Lung Squamous Cell Carcinoma and Adenocarcinoma. J. Comput. Biol. 27, 1067–1078. doi:10.1089/cmb.2019.0164

PubMed Abstract | CrossRef Full Text | Google Scholar

Torre, L. A., Siegel, R. L., and Jemal, A. (2016). Lung Cancer Statistics. Adv. Exp. Med. Biol. 893, 1–19. doi:10.1007/978-3-319-24223-1_1

PubMed Abstract | CrossRef Full Text | Google Scholar

Wang, H., Meng, Q., and Ma, B. (2021). Characterization of the Prognostic m6A-Related lncRNA Signature in Gastric Cancer. Front. Oncol. 11, 630260. doi:10.3389/fonc.2021.630260

PubMed Abstract | CrossRef Full Text | Google Scholar

Wojtas, M. N., Pandey, R. R., Mendel, M., Homolka, D., Sachidanandam, R., and Pillai, R. S. (2017). Regulation of m6A Transcripts by the 3ʹ→5ʹ RNA Helicase YTHDC2 Is Essential for a Successful Meiotic Program in the Mammalian Germline. Mol. Cell 68, 374–387. doi:10.1016/j.molcel.2017.09.021

PubMed Abstract | CrossRef Full Text | Google Scholar

Woodman, C., Vundu, G., George, A., and Wilson, C. M. (2021). Applications and Strategies in Nanodiagnosis and Nanotherapy in Lung Cancer. Seminars Cancer Biol. 69, 349–364. doi:10.1016/j.semcancer.2020.02.009

CrossRef Full Text | Google Scholar

Wu, B., Su, S., Patil, D. P., Liu, H., Gan, J., Jaffrey, S. R., et al. (2018). Molecular Basis for the Specific and Multivariant Recognitions of RNA Substrates by Human hnRNP A2/B1. Nat. Commun. 9, 420. doi:10.1038/s41467-017-02770-z

PubMed Abstract | CrossRef Full Text | Google Scholar

Wu, Y., Yang, X., Chen, Z., Tian, L., Jiang, G., Chen, F., et al. (2019). m6A-induced lncRNA RP11 Triggers the Dissemination of Colorectal Cancer Cells via Upregulation of Zeb1. Mol. Cancer 18, 87. doi:10.1186/s12943-019-1014-2

PubMed Abstract | CrossRef Full Text | Google Scholar

Xu, F., Huang, X., Li, Y., Chen, Y., and Lin, L. (2021). m6A-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

Xue, L., Li, J., Lin, Y., Liu, D., Yang, Q., Jian, J., et al. (2021). m 6 A Transferase METTL3‐induced lncRNA ABHD11‐AS1 Promotes the Warburg Effect of Non‐small‐cell Lung Cancer. J. Cell Physiol. 236, 2649–2658. doi:10.1002/jcp.30023

PubMed Abstract | CrossRef Full Text | Google Scholar

Yang, D., Qiao, J., Wang, G., Lan, Y., Li, G., Guo, X., et al. (2018). N 6-Methyladenosine Modification of lincRNA 1281 Is Critically Required for mESC Differentiation Potential. Nucleic Acids Res. 46, 3906–3920. doi:10.1093/nar/gky130

PubMed Abstract | CrossRef Full Text | Google Scholar

Yang, Z., Li, J., Feng, G., Gao, S., Wang, Y., Zhang, S., et al. (2017). MicroRNA-145 Modulates N6-Methyladenosine Levels by Targeting the 3′-Untranslated mRNA Region of the N6-Methyladenosine Binding YTH Domain Family 2 Protein. J. Biol. Chem. 292, 3614–3623. doi:10.1074/jbc.M116.749689

PubMed Abstract | CrossRef Full Text | Google Scholar

Yoshihara, K., Shahmoradgoli, M., Martínez, E., Vegesna, R., Kim, H., Torres-Garcia, W., et al. (2013). Inferring Tumour Purity and Stromal and Immune Cell Admixture from Expression Data. Nat. Commun. 4, 2612. doi:10.1038/ncomms3612

PubMed Abstract | CrossRef Full Text | Google Scholar

Yu, H., Boyle, T. A., Zhou, C., Rimm, D. L., and Hirsch, F. R. (2016). PD-L1 Expression in Lung Cancer. J. Thorac. Oncol. 11, 964–975. doi:10.1016/j.jtho.2016.04.014

PubMed Abstract | CrossRef Full Text | Google Scholar

Yu, J., Mao, W., Sun, S., Hu, Q., Wang, C., Xu, Z., et al. (2021). Identification of an m6A-Related lncRNA Signature for Predicting the Prognosis in Patients with Kidney Renal Clear Cell Carcinoma. Front. Oncol. 11, 663263. doi:10.3389/fonc.2021.663263

PubMed Abstract | CrossRef Full Text | Google Scholar

Yuan, H., Liu, J., and Zhang, J. (2021). The Current Landscape of Immune Checkpoint Blockade in Metastatic Lung Squamous Cell Carcinoma. Molecules 26, 1392. doi:10.3390/molecules26051392

PubMed Abstract | CrossRef Full Text | Google Scholar

Zaccara, S., Ries, R. J., and Jaffrey, S. R. (2019). Reading, Writing and Erasing mRNA Methylation. Nat. Rev. Mol. Cell Biol. 20, 608–624. doi:10.1038/s41580-019-0168-5

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhang, C., Leighl, N. B., Wu, Y.-L., and Zhong, W.-Z. (2019). Emerging Therapies for Non-small Cell Lung Cancer. J. Hematol. Oncol. 12, 45. doi:10.1186/s13045-019-0731-8

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhao, D., Zhao, H., He, Y., Yang, Y., Du, Y., and Zhang, M. (2021). The Inhibitive Effects of Proteasome Inhibitor MG-132 on Pterygium Fibroblasts In Vitro and the Potential Key Regulators Involved. Life Sci. 270, 119088. Epub 2021 Jan 20. doi:10.1016/j.lfs.2021.119088

PubMed Abstract | CrossRef Full Text | Google Scholar

Zheng, J., Zhao, Z., Wan, J., Guo, M., Wang, Y., Yang, Z., et al. (2021). N‐6 Methylation‐related lncRNA Is Potential Signature in Lung Adenocarcinoma and Influences Tumor Microenvironment. J. Clin. Lab. Anal. 35, e23951. doi:10.1002/jcla.23951

PubMed Abstract | CrossRef Full Text | Google Scholar

Zheng, Z.-Q., Li, Z.-X., Zhou, G.-Q., Lin, L., Zhang, L.-L., Lv, J.-W., et al. (2019). Long Noncoding RNA FAM225A Promotes Nasopharyngeal Carcinoma Tumorigenesis and Metastasis by Acting as ceRNA to Sponge miR-590-3p/miR-1275 and Upregulate ITGB3. Cancer Res. 79, 4612–4626. doi:10.1158/0008-5472.CAN-19-0799

PubMed Abstract | CrossRef Full Text | Google Scholar