IER5L is a Prognostic Biomarker in Pan-Cancer Analysis and Correlates with Immune Infiltration and Immune Molecules in Non-Small Cell Lung Cancer

Introduction

Lung cancer is the leading cause of cancer-related deaths worldwide.¹ Based on pathological point of view, lung cancer is divided into two types: non-small cell lung cancer (NSCLC) and small cell lung cancer (SCLC),² and NSCLC accounts for more than 80% of all cases.³ Unfortunately, many NSCLCs are determined at advanced stages;⁴ thus, many of them have lost surgical opportunity. Moreover, although many therapies have been developed for lung cancer, including chemotherapy, targeted therapy, and immunotherapy, its average 5-year overall survival (OS) probability remains less than 20%.^5–7 Thus, identification of effective biomarkers to predict prognosis and even guide therapy is important for NSCLC.

Various biomarkers for diagnosis and prognosis have been identified for NSCLC using the Gene Expression Omnibus (GEO) and The Cancer Genome Atlas (TCGA) databases, such as ADAM-like decysin-1 (ADAMDEC1),⁸ long non-coding RNASEH1-AS1,⁹ sialic acid binding Ig-like lectin 15 (Siglec15),¹⁰ lysine acetyltransferase 2B (KAT2B),¹¹ and Erb-b2 receptor tyrosine kinase 4 (ERBB4).¹² Moreover, many gene signatures have been constructed to predict prognosis and therapeutic response for NSCLC, such as RNA-binding proteins related signature,¹³ telomere length-related gene signature,¹⁴ ferroptosis-related long non-coding RNAs signature,¹⁵ immune microenvironment signature,¹⁶ neutrophil extracellular traps-related lncRNA signature,¹⁷ and glycolysis-related gene signature.¹⁸ However, owing to complex molecular mechanisms and disease heterogenicity in NSCLC, more biomarkers having different targets and prognostic value need to be explored.

Immediate-early response family has three different members, which is activated under different conditions.¹⁹ IER2 is activated through activating P38 and JNK, whereas IER5 expression is increased via proteotoxic stress.¹⁹ In recent years, researchers have reported that activated immediate early response 5 like (IER5L) is correlated with tumorigenesis and its progression. In glioblastoma, there is increased IER5L expression in cancer tissues compared with control tissues.²⁰ Moreover, a growth on IER5L expression has been reported in metastatic colorectal cancer compared with primary cancer, and is negatively correlated with survival.²¹ However, there are few pan-cancer analyses of IER5L and reports on the role of IER5L in NSCLC, so this study aimed to perform a pan-cancer analysis of IER5L and evaluate its value as a prognostic biomarker and a potential immunotherapeutic target for NSCLC.

In the present study, data were downloaded from several public databases. Pan-cancer analysis of IER5L were first performed in terms of gene expression, prognosis, and disease progression. In NSCLC, relative expression, prognostic value, and correlation with clinical parameters were analyzed using R language. Quantitative polymerase chain reaction (qPCR) and single-cell RNA sequencing (scRNA seq) analysis were used to validate dysregulated IER5L expression in NSCLC. Upstream and downstream mechanisms were explored using transcriptional network that was constructed using microRNA (miRNA) and transcription factors (TFs) and functional enrichment analysis, a protein–protein interaction (PPI) network, and gene set enrichment analysis (GSEA), respectively. Methylation levels of IER5L and tumor mutation burden (TMB) between the two groups were analyzed using R language. Finally, correlations of IER5L with immune infiltration and immune molecules (inhibitors and stimulators) were analyzed using R language to evaluate its value in immunotherapy.

Methods and Materials

Data Collection and Processing

The mRNA expression data and clinical information for multiple cancers were obtained from the Gene Expression Profiling Interactive Analysis 2 (GEPIA2) database (http://gepia2.cancer-pku.cn/#index). Moreover, mRNA expression profiles (1037 cancer and 108 controls), methylation data (807 cancer and 71 controls), gene mutation data (1059 cancer), and clinical information (1027 cancer) were downloaded from the TCGA-GDC database (https://portal.gdc.cancer.gov/repository). GSE31210 (226 cancer and 20 controls) and GSE30219 (293 cancer and 14 controls) datasets were downloaded from the GEO database (https://www.ncbi.nlm.nih.gov/gds/). GSE31210 and GSE30219 were generated from the GPL570 [HG-U133_Plus_2] Affymetrix Human Genome U133 Plus 2.0 Array. Inclusion criteria were patients with 1) NSCLC; 2) IER5L mRNA expression; and 3) complete follow-up data. Therefore, 987 patients from TCGA, 226 patients from GSE31210, and 268 patients from GSE30219 were included in the current study. The clinical characteristics of patients are summarized in Table 1.

Table 1 Clinical Characteristics of Patients with NSCLC from Public Databases

Levels of IER5L mRNA Expression and Its Prognostic Value in a Pan-Cancer Analysis

Levels of IER5L mRNA expression were estimated between cancer tissues and control tissues in 23 types of cancer and among different International Cancer Control (UICC) stages in the GEPIA2 database. Moreover, patients with each kind of cancer were divided into high or low IER5L expression groups based on the median IER5L mRNA expression, and survival analysis was performed between two IER5L expression groups.

Levels of IER5L mRNA Expression in NSCLC and Controls

Levels of IER5L mRNA expression in NSCLC were compared with controls in TCGA, GSE31210, and GSE30219 datasets. The diagnostic value of IER5L was evaluated using a receiver operating characteristic (ROC) curve.

RNA Extraction and qPCR Validation

Thirty-two lung tissues (16 NSCLC and 16 controls) were obtained from patients with NSCLC in Zigong First People’s Hospital. The baseline characteristics of patients are listed in Table 2. Total RNA was first extracted from lung tissues using the E.Z.N.A. HP Total RNA Kit (Omega, GA, USA) and was then applied to synthesize cDNA using The PrimeScript^TM RT reagent kit (Takara, Japan) based on the manufacturer’s protocol. Finally, qPCR was performed using cDNA via iQ™ SYBR Green Supermix (Bio-Rad, USA) according to the manufacturer’s instructions. Levels of IER5L mRNA expression were normalized using β-actin Ct value in a 2^−ΔΔCt quantification method. The qPCR primers were as follows.

Table 2 The Clinical Characteristics of Patients with NSCLC from Our Cohort

IER5L-forward: 5′-AGACAAAGACTCGGCAAAGGG-3′.

IER5L-reverse: 5′-TGGTAGGCTCCTCGCTAAACA-3′.

β-actin′-forward: 5′-CCACGAAACTACCTTCAACTCC′-3′.

β-actin′-reverse: 5′-GTGATCTCCTTCTGCATCCTGT′-3′.

Single-Cell RNA Sequencing Analysis

Six single-cell RNA sequencing datasets for NSCLC were downloaded from the Tumor Immune Single-cell Hub database (TISCH, http://tisch1.comp-genomics.org/home/) including NSCLC_EMTAB6149, NSCLC_GSE117570, NSCLC_GSE127465, NSCLC_GSE127471, NSCLC_GSE139555, and NSCLC_GSE99254. A single-cell atlas was plotted, and IER5L mRNA expression was compared among different cells.

The Prognostic Value of IER5L and Its Correlation with Clinical Characteristics

The high or low IER5L expression groups were constructed using NSCLC patients based on the median IER5L mRNA expression. The OS probability was estimated between the two groups in TCGA, GSE31210, and GSE30219, and the DFS probability was calculated in GSE30219. Moreover, the levels of IER5L mRNA expression were compared between the different UICC stages in TCGA, GSE31210, and GSE30219, and were also compared between patients with relapse and patients without relapse in GSE30219. Independent prognostic factors for the NSCLC were screened using univariate and multivariate Cox regression analyses with the “survival” package in TCGA, GSE31210, and GSE30219.

Nomogram and Calibration Plots of Nomogram

A nomogram was constructed for NSCLC using independent prognostic factors screened by multivariate Cox regression analysis with the “rms” package in R language (version 4.0.2). Calibration plots of the nomogram were used to assess conformity of the nomogram-predicted with actual OS probability.

IER5L Methylation Levels and Gene Mutation Analysis

Methylation data of IER5L in NSCLC were obtained using Strawberry Perl (version 5.34.1). Correlation of IER5L methylation levels with IER5L mRNA expression levels was evaluated using Spearman rank correlation analysis. The “maftools” package in R language was applied to analyze gene mutation frequencies between the two groups in both LUAD and LUSC. Moreover, TMB was calculated using Strawberry Perl in both LUAD and LUSC, and correlation of IER5L expression levels with TMB was evaluated using Spearman rank correlation analysis.

Construction of the Upstream Regulatory Network of IER5L

TFs were obtained from the Harmonizome (https://maayanlab.cloud/Harmonizome/) and GRNdb (http://www.grndb.com/) databases. miRNAs of IER5L were acquired from the miRWalk (http://mirwalk.umm.uni-heidelberg.de/), TargetScan (https://www.targetscan.org/vert_80/), and mirDIP (http://ophid.utoronto.ca/mirDIP/) databases.

Identification of Differentially Expressed Genes (DEGs) and Construction of a PPI Network

DEGs were identified between the two IER5L expression groups in NSCLC using the “limma” package in R language, and the criteria were |log₂ fold change (FC)|≥1 and adjusted P< 0.05. A heatmap was generated between the two IER5L expression groups using the “pheatmap” package in R language. The STRING online database (https://string-db.org, version 11.0) was used to construct a PPI network.²² The top 200 DEGs were uploaded into the STRING online database to construct a PPI network, tailoring to “human”, “minimum required interaction score: 0.04”, and “network display options: hide disconnected nodes in the network”, and a “string_interactions_short.tsv” file was acquired. This file was then imported into Cytoscape (version 3.7.1) to identify hub genes using the MCC algorithm of cytoHubba app.

Functional Enrichment Analysis and GSEA

Functional enrichment analysis was performed using “clusterProfiler”, “org.Hs.eg.db”, “enrichplot” and “ggplot2” packages in R language, including biological process (BP), molecular function (MF), and cellular component (CC). Terms with adjusted P< 0.05 were considered of statistical significance. GSEA can determine whether a defined set of genes shows statistically significant and concordant differences between the two biological states.²³ GSEA software (version 4.0.3) and “c2.cp.kegg.v7.5.1.symbols.gmt” were obtained from the GSEA database (http://www.gsea-msigdb.org/gsea/index.jsp) and was used to perform pathway analysis between the two expression groups. Statistical significance was defined as false discovery rate (FDR)<0.25 and adjusted P<0.05.

Correlation of IER5L mRNA Expression with Immune Infiltration

Tumor micro-environment score was acquired using the “limma” and “e1071” packages in R language. The “vioplot” package in R language was used to calculate immune infiltration profiles between the two expression groups. In addition, correlations of IER5L expression with immune infiltration were assessed using Spearman rank correlation analysis in R language, and statistical significance was set at a P-value of <0.05 and |r|>0.1.

Correlation of IER5L mRNA Expression with Immune Molecules

Correlations of IER5L mRNA expression with multiple immune molecules (inhibitors and stimulators) in different cancers were obtained from the TISIDB database (http://cis.hku.hk/TISIDB/index.php). Correlations were assessed using Spearman rank correlation analysis.

Statistical Analysis

R language (version 4.0.2) and GraphPad Prism (version 7.00) were applied to perform statistical analysis and to draw all figures. Shapiro–Wilk normality test was used to evaluate data distribution type. Levels of IER5L relative expression were displayed as median (interquartile range). Comparisons between the two groups were performed using the Mann–Whitney test for nonparametric data. Kruskal–Wallis test was performed for multiple comparisons for nonparametric data. Correlation analysis was performed using Spearman rank correlation analysis. Survival curve was compared using Log rank test. Statistical significance was set at a P-value of <0.05.

Results

IER5L Was Associated with Prognosis and Progression in a Pan-Cancer Analysis

We firstly performed a pan-cancer analysis of IER5L in GEPIA2 database. Levels of IER5L mRNA expression were increased in bladder urothelial carcinoma, breast invasive carcinoma (BRCA), cholangiocarcinoma, colon adenocarcinoma, esophageal carcinoma (ESCA), head and neck squamous cell carcinoma, LUAD, LUSC, rectum adenocarcinoma (READ), stomach adenocarcinoma, thyroid carcinoma (THCA), and uterine corpus endometrial carcinoma when compared with controls (P<0.05, Figure 1A–C), while levels of IER5L mRNA expression were lower in kidney chromophobe (KICH) than corresponding controls (P< 0.05, Figure 1A). Kaplan–Meier survival curves showed that the high IER5L expression was correlated with a lower OS probability than the low IER5L expression in cervical and endocervical cancer (CESC) (Log-rank P= 0.072), kidney papillary cell carcinoma (KIRP) (Log-rank P= 0.064), brain lower grade glioma (Log-rank P< 0.001), liver hepatocellular carcinoma (LIHC) (Log-rank P= 0.017), NSCLC (Log-rank P< 0.001), ovarian serous cystadenocarcinoma (Log-rank P= 0.068), READ (Log-rank P= 0.025), and uveal melanoma (UVM) (Log-rank P= 0.008) (Figure 2A–H). Moreover, patients with the high IER5L expression displayed a worse DFS in comparison with patients with the low IER5L expression in BRCA (Log-rank P= 0.037), KIRP (Log-rank P= 0.02), NSCLC (Log-rank P= 0.007), pancreatic adenocarcinoma (Log-rank P= 0.063), prostate adenocarcinoma (Log-rank P< 0.001), READ (Log-rank P= 0.047), and UVM (Log-rank P= 0.026) (Figure 2I–O). Besides, levels of IER5L mRNA expression were associated with tumor progression in adrenocortical carcinoma (P= 0.070), CESC (P= 0.021), diffuse large B-cell lymphoma (P= 0.018), ESCA (P= 0.042), KICH (P= 0.025), LIHC (P= 0.021), testicular germ cell tumor (P= 0.051), THCA (P= 0.002), uterine carcinosarcoma (P= 0.019) (Figure 3A–I). Therefore, these findings show that IER5L expression is upregulated in multiple cancers and is associated with disease progression and prognosis.

Figure 1 Levels of IER5L expression in 23 kinds of cancers and corresponding controls from GEPIA2 database. (A) BLCA, BRCA, CESC, CHOL, COAD, ESCA, NHSC. (B) KICH, KIRC, KIRP, LIHC, LUAD, LUSC, PAAD, PCPG. (C) PRAD, READ, SARC, SKCM, STAD, THCA, THYM, UCEC. *P < 0.05. Abbreviations: IER5L, immediate early response 5 like; GEPIA2, Gene Expression Profiling Interactive Analysis 2.

Figure 2 Kaplan-Meier survival analysis of IER5L in different tumors. (A-H) Overall survival between the two IER5L groups in CESC, KIRP, LGG, LIHC, NSCLC, OV, READ, and UVM. (I-O) Disease free survival between the two IER5L groups in BRCA, KIRP, NSCLC, PAAD, PRAD, READ, and UVM. Abbreviation: IER5L, immediate early response 5 like.

Figure 3 Correlation of IER5L expression with tumor stages. (A) ACC. (B) CESC. (C) DLBC. (D) ESCA. (E) KICH. (F) LIHC. (G) TGCT. (H) THCA. (I) UCS. IER5L, immediate early response 5 like.

IER5L as a Diagnostic Biomarker in NSCLC

Because of its increasing prevalence and mortality, NSCLC was selected as the focus of this study. Levels of IER5L mRNA expression were significantly higher in NSCLC tissues than in control tissues in the TCGA (P< 0.001, Figure 4A). The area under the curve (AUC) value was 0.882 (0.857–0.906), and cutoff value was 2.595 with a sensitivity of 74.25% and a specificity of 91.67% (Figure 4B). In addition, there were increased levels of IER5L expression in NSCLC in comparison with controls in both GSE31210 (P= 0.005; Figure 4C) and GSE30219 (P< 0.001; Figure 4D) datasets. qPCR was performed to validate the findings of bioinformatic analysis, and the results showed that levels of IER5L mRNA expression in NSCLC tissues were higher than those in control tissues (P= 0.015, Figure 4E).

Figure 4 Levels of IER5L expression and its diagnostic value in NSCLC. (A) Levels of IER5L expression in NSCLC and corresponding controls from the TCGA. (B) Receiver operating characteristic curve of IER5L in NSCLC. (C-E) Levels of IER5L expression in NSCLC and corresponding controls from the GSE31210, GSE30219, and qPCR validation. IER5L, immediate early response 5 like; NSCLC, non-small cell lung cancer; TCGA, The Cancer Genome Atlas. *P < 0.05, **P < 0.01, ***P < 0.001. Abbreviation: qPCR, quantitative polymerase chain reaction.

Single-Cell RNA Sequencing Analysis

Single-cell RNA sequencing analysis was performed to explore source of overexpressed IER5L in NSCLC. A heatmap of six datasets is displayed in Figure 5A. IER5L was mainly expressed in malignant cells in NSCLC_GSE117570 (Figure 5A and B), in endothelial cells and fibroblasts in NSCLC_GSE127465 (Figure 5A and C), and in monocytes and macrophages in both NSCLC_GSE139555 (Figure 5A and D) and NSCLC_GSE99254 (Figure 5A and E).

Figure 5 Single cell-RNA sequencing analysis in NSCLC. (A) A heatmap of six datasets. (B-E) Single-cell atlas and violin plot of IER5L in NSCLC_GSE117570, NSCLC_GSE127465, NSCLC_GSE139555, and NSCLC_GSE99254. Abbreviations: NSCLC, non-small cell lung cancer; IER5L, immediate early response 5 like.

IER5L Was a Prognostic Biomarker and Was Correlated with Clinical Characteristics in NSCLC

Kaplan–Meier survival curves showed that patients with the high IER5L mRNA expression had a lower OS probability than patients with the low IER5L mRNA expression in the TCGA (P< 0.001, Figure 6A), GSE31210 (P< 0.001, Figure 6B), and GSE30219 (P< 0.001, Figure 6C). In addition, the DFS probability was significantly lower in the high IER5L expression group than in the low IER5L expression group in GSE30219 (P< 0.001, Figure 6D). Levels of IER5L mRNA expression were higher in patients with stage III–IV cancer than in patients with stage I–II cancer in TCGA (P= 0.001, Figure 6E) and were increased in patients with stage II cancer than in patients with stage I cancer in GSE31210 (P< 0.001, Figure 6F). However, levels of IER5L mRNA expression were not a significant difference between stage I–II and stage III–IV in GSE30219 (P> 0.05, Figure 6G). The patients with relapse had higher levels of IER5L mRNA expression compared with the patients without relapse in GSE30219 (P= 0.039, Figure 6H).

Figure 6 The clinical value of IER5L in NSCLC. (A-C) OS probability in the TCGA, GSE31210, and GSE30219. (D) DFS probability in the GSE30219. (E-G) Levels of IER5L expression in different UICC stages in TCGA, GSE31210, and GSE30219. (H) Levels of IER5L expression in different relapse groups in GSE30219. *P < 0.05, **P < 0.01. Abbreviations: NSCLC, non-small cell lung cancer; IER5L, immediate early response 5 like; OS, overall survival; DFS, disease free survival; ns: no significance.

IER5L as an Independent Prognostic Factor in NSCLC

Univariate Cox regression analysis reported that age (HR: 1.014; P=0.025), UICC stage (HR: 1.475; P<0.001), and IER5L expression (HR: 1.034; P<0.001) were correlated with NSCLC prognosis (Figure 7A), which was confirmed using multivariate Cox regression analysis (for age, HR:1.020, P=0.002; for UICC stage, HR: 1.499; P<0.001; for IER5L expression, HR: 1.037; P<0.001, Figure 7B) in TCGA. In GSE30219, age (HR: 1.034; P<0.001), T stage (HR: 1.568; P<0.001), N stage (HR: 1.820; P<0.001), and IER5L expression (HR:1.536; P=0.003) were correlated with prognosis using univariate Cox regression analysis (Figure 7C), and age (HR: 1.036; P<0.001), N stage (HR: 1.646; P<0.001), and IER5L expression (HR:1.540; P=0.003) were correlated with prognosis using multivariate Cox regression analysis (Figure 7D). In GSE31210, UICC stage (HR: 4.232; P<0.001) and IER5L expression (HR: 1.001; P=0.035) were correlated with NSCLC prognosis (Figure 7E), however, only UICC stage was correlated with prognosis using multivariate Cox regression analysis (HR: 3.928; P< 0.001, Figure 7F).

Figure 7 Univariate and multivariate Cox regression analysis, and a nomogram and its calibration plots. Univariate and multivariate Cox regression analyses in NSCLC in (A and B) TCGA, (C and D) GSE30219, and (E and F) GSE31210. (G) Nomogram. (H) Calibration plots of the nomogram. **P < 0.01, ***P < 0.001. Abbreviations: IER5L, immediate early response 5 like; OS, overall survival.

Nomogram and Calibration Plots of Nomogram

A nomogram was generated to predict OS probability, including age, UICC stage, and IER5L expression, in patients with NSCLC. The nomogram-predicted OS probability at 1, 3, and 5 years was calculated by plotting a vertical line between the total point axis and each prognostic axis (Figure 7G). Calibration plots of nomogram showed the high conformity of the nomogram-predicted and actual OS probability at 1, 3, and 5 years (Figure 7H).

IER5L Expression Was Correlated with TMB in NSCLC

Levels of IER5L methylation were negatively correlated with levels of IER5L mRNA expression without a statistical difference (rho=−0.05, P=0.13, Figure 8A). Genes with the highest mutation frequencies were tumor protein p53 (TP53) and titin (TTN) in both IER5L expression groups of LUAD and LUSC (Figure 8B-E). In addition, levels of IER5L mRNA expression were positively correlated with TMB in NSCLC (rho=0.065, P=0.043, Figure 8F).

Figure 8 Methylation and mutation analysis. (A) Correlation of IER5L expression with IER5L methylation. (B, C) Gene mutation frequencies in the high and low IER5L expression groups in lung adenocarcinoma. (D and E) Gene mutation frequencies in the high and low IER5L expression groups in lung squamous cell carcinoma. (F) Correlation of IER5L expression with TMB. Abbreviations: IER5L, immediate early response 5 like; TMB, tumor mutation burden.

An Upstream Regulatory Network of IER5L in NSCLC

There were 35 and 6 TFs from the Harmonizome and GRNdb databases, respectively, and early growth response 1 (EGR1) was the only overlapping TF between the two databases (Figure 9A). Levels of EGR1 mRNA expression were decreased in NSCLC compared with controls (P< 0.001, Figure 9B). There were 252 miRNAs, 394 miRNAs, and 10 miRNAs from miRWalk, TargetScan, and mirDIP databases, respectively, and hsa-miR-6789-3p and hsa-miR-8075 were overlapping miRNAs between the three databases (Figure 9C). hsa-miR-8075 has been associated with tumor cell proliferation and invasion.^24,25 Therefore, an EGR1-hsa-miR-8075-IER5L network was constructed in NSCLC (Figure 9D).

Figure 9 A miRNA-transcription factor-mRNA network. (A) Overlapping transcription factors between Harmonizome and GRNdb databases. (B) Levels of EGR1 relative expression in NSCLC. (C) Overlapping microRNAs among miRWalk, TargetScan, and mirDIP databases. (D) An EGR1-hsa-miR-8075-IER5L network. Abbreviations: NSCLC, non-small cell lung cancer; IER5L, immediate early response 5 like. EGR1, early growth response 1.

Identification of DEGs, Construction of a PPI Network, Functional Enrichment Analysis, and GSEA

A total of 191 DEGs (38 upregulated and 153 downregulated genes) were identified when comparing the high IER5L expression group with the low IER5L expression group (|log₂ fold FC|≥1 and adjusted P< 0.05). A heatmap is illustrated in Figure 10A. A PPI network was generated in the STRING online database and the top 30 hub genes were screened using the MCC method in Cytoscape, among which the top five hub genes were neurotrophic receptor tyrosine kinase 2 (NTRK2), contactin 2 (CNTN2), synaptoporin (SYNPR), glutamate ionotropic receptor AMPA type subunit 1 (GRIA1), and keratin 20 (KRT20) (Figure 10B). These DEGs were significantly enriched in BP such as action potential, sodium ion transport, and regulation of membrane potential; CC such as intermediate filament, intermediate filament cytoskeleton, and cation channel complex; and MF such as channel activity, passive transmembrane transporter activity, and monooxygenase activity (adjusted P< 0.05, Figure 10C). GSEA showed that pathways correlated with the high IER5L expression were primarily enriched in cell cycle, DNA replication, Notch signaling pathway, oxidative phosphorylation, base excision repair, and nucleotide excision repair (FDR<0.25 and nominal P<0.05; Figure 10D-I).

Figure 10 A PPI network, functional enrichment analysis, and gene set enrichment analysis. (A) A heatmap of two IER5L expression groups. (B) A PPI network. (C) Functional enrichment analysis, including BP, CC, and MF. (D) Cell cycle. (E) DNA replication. (F) Notch signaling pathway, (G) Oxidative phosphorylation. (H) Base excision repair. (I) Nucleotide excision repair. Abbreviations: PPI, protein–protein interaction; BP, biological process; CC, cellular component; MF, molecular function; IER5L, immediate early response 5 like.

IER5L Expression Was Correlated with Immune Infiltration

There was an increased number of regulatory T (Treg) cells (P<0.001), neutrophils (P<0.05), and activated natural killer (NK) cells (P<0.001) in the high IER5L expression group in comparison with the low IER5L expression group (Figure 11A); whereas a reduced number of naive B cells (P<0.001), resting memory CD4+ T cells (P<0.01), resting dendritic cells (P<0.001), resting mast cells (P<0.001), and eosinophils (P<0.01) in the high IER5L expression group compared with the low IER5L expression group (Figure 11A). Levels of IER5L expression were positively correlated with the number of Treg cells (rho=0.19, P<0.001) and activated NK cells (rho=0.14, P<0.001) (Figure 11B and C), and were negatively correlated with the number of naive B cells (rho=−0.13, P<0.001), resting memory CD4+ T cells (rho=−0.12, P<0.001), resting dendritic cells (rho=−0.17, P<0.001), and resting mast cells (rho=−0.17, P<0.001) (Figure 11D-G).

Figure 11 Correlations of IER5L expression with immune infiltration. (A) Immune infiltration profiles between the two IER5L expression groups. Correlations of IER5L expression with (B) Treg cells, (C) activated NK cells, (D) naive B cells, (E) resting memory CD4+ T cells, (F) resting dendritic cells, (G) resting mast cells. *P < 0.05, **P < 0.01, ***P < 0.001. Abbreviations: IER5L, immediate early response 5 like; Treg, regulatory T; NK, natural killer.

IER5L Expression Was Correlated with Immune Molecules

We next explored the correlation of IER5L with immune molecules including inhibitors and stimulators in the TISIDB database. IER5L mRNA expression was positively correlated with many immune inhibitors in NSCLC (Figure 12A). The top three immune inhibitors were nectin cell adhesion molecule 2 (PVRL2) (rho=0.274, P<0.001), lymphocyte activating 3 (LAG3) (rho=0.161, P<0.001), and programmed cell death 1 (PDCD1) (rho=0.154, P<0.001) in LUAD (Figure 12B), and the top three immune inhibitors were transforming growth factor beta 1 (TGFB1) (rho=0.337, P<0.001), PVRL2 (rho=0.336, P<0.001), and adenosine A2a receptor (ADORA2A) (rho=0.082, P=0.066) in LUSC (Figure 12C). Besides, IER5L expression was positively correlated with many immune stimulators in NSCLC (Figure 12D). The top three immune stimulators were TNF receptor superfamily member 4 (TNFRSF4) (rho=0.316, P<0.001), TNF receptor superfamily member 18 (TNFRSF18) (rho=0.26, P<0.001), and CD276 molecule (CD276) (rho=0.252, P<0.001) in LUAD (Figure 12E), and the top three immune stimulators were PVR cell adhesion molecule (PVR) (rho=0.276, P<0.001), TNF receptor superfamily member 25 (TNFRSF25) (rho=0.267, P<0.001), and retinoic acid early transcript 1E (RAET1E) (rho = 0.246, P< 0.001) in LUSC (Figure 12F).

Figure 12 Correlations of IER5L expression with immune molecules. (A) Correlations of IER5L expression with immune inhibitors in multiple cancers. (B and C) The top 3 immune inhibitors in lung adenocarcinoma and lung squamous cell carcinoma. (D) Correlations of IER5L expression with immune stimulators in multiple cancers. (E and F) The top 3 immune stimulators in lung adenocarcinoma and lung squamous cell carcinoma. Abbreviation: IER5L, immediate early response 5 like.

Discussions

In the current study, our findings showed that IER5L expression was dysregulated in multiple cancers and was also correlated with progression and prognosis in multiple cancers. In NSCLC. IER5L was a diagnostic biomarker and was an independent prognostic factor. An EGR1-hsa-miR-8075-IER5L network was constructed for NSCLC. DEGs are significantly enriched in the biological process of action potential, sodium ion transport, and regulation of membrane potential. Increased IER5L expression was primarily enriched in pathways of cell cycle, NOTCH signaling, and oxidative phosphorylation and was correlated with increased regulatory T cells and neutrophils, elevated levels of immune molecules, and high TMB.

Previous studies have shown that levels of IER5L expression were increased in both primary tumor tissues and metastatic tumor tissues and were negatively correlated with survival.^20,21 Our study showed that IER5L expression was upregulated in many kinds of tumors in various systems and was associated with UICC stages and prognosis. Thus, in combination with the previous findings, we could conclude that IER5L may be a diagnostic and prognostic biomarker in cancer. Similar results have also been confirmed in NSCLC using several datasets, and increased IER5L expression has been validated using qPCR. In addition, univariate and multivariate Cox regression analyses reported that IER5L expression was correlated with NSCLC prognosis after adjustment for age and UICC stage. Thus, IER5L is an independent prognostic factor for NSCLC. A nomogram was generated using independent prognostic factors, and calibration plots of the nomogram showed highly predictive accuracy.

We then explore the underlying downstream mechanisms between the two IER5L expression groups using GSEA. The oxidative phosphorylation, Notch signaling pathway, cell cycle, and DNA replication were the primary pathways in our study. Increasing evidences reveal that oxidative phosphorylation is more active in cancer cells than in normal cells,²⁶ Attenuating cell oxidative phosphorylation could suppress cancer cell migration and invasion.^27,28 The Notch signaling pathway has been implicated in the pathogenesis of many cancers,^29–31 and it influences tumor microenvironment,³² and predicts prognosis, recurrence, and response to therapy.^33–36 Deregulation in cell cycle progression result in the unconstrained tumor cell proliferation.³⁷ DNA replication must be precisely controlled to maintain genome stability,³⁸ and genome instability increases the risk of cancer.³⁹ Thus, IER5L may be involved in the tumorigenesis, progression, and prognosis of NSCLC through the above-mentioned pathways. However, these pathways may not function as an isolated state; instead, they may interact to influence roles of tumor cells. Thus, more researches were needed to explore potential mechanisms between these pathways and role of IER5L in these pathways. TMB has been associated with prognosis of disease and immunotherapeutic response in multiple cancers.⁴⁰ Our study has reported that increased IER5L expression was correlated with higher TMB. Thus, the worse prognosis in the high IER5L expression may be correlated with increased TMB. These results show that interfering IER5L may prevent progression and improve prognosis in NSCLC, and it deserves further investigation.

The tumor microenvironment is mainly composed of malignant cells, fibroblasts, immune cells, and endothelial cells,⁴¹ and they play critical roles in initiation, progression, and prognosis in cancer.⁴² Treg cells are abundant in NSCLC and have an increased risk of recurrence in early-stage cancer,⁴³ and could predict the efficacy of PD-1/PD-L1 inhibitors.⁴⁴ Additionally, increased tumor-infiltrating Treg cells are associated with lower OS and DFS probability in lung cancer.^45–47 Lung cancer tissues have more neutrophil infiltration and neutrophil extracellular traps (NETs) formation.⁴⁸ Neutrophils are higher in advanced NSCLC⁴⁹ and are independently associated with survival.⁵⁰ Moreover, NETs can enable A549 cells to invade and migrate in a macrophage-maintained inflammatory microenvironment.⁴⁸ Brain metastasis of lung cancer is correlated with nicotine-induced N2-neutrophils.⁵¹ Resting mast cells and eosinophils are significantly lower in lung cancer tissues.⁵² The increased Treg cells and neutrophils as well as reduced resting mast cells and eosinophils were reported in the high IER5L expression group in the current study. Thus, IER5L was correlated with a dysregulated immune microenvironment, and may be correlated with response on immunotherapy. Our results also showed that IER5L expression was positively correlated with many immune molecules. Increased PD1 levels were significantly associated with better response to immunotherapy and worse OS probability,^53,54 and PD1/PDL1 complex had a higher predictive value compared with single PD1.⁵⁵ QL1706, a mixture of anti-PD1 and anti-CTLA4, showed promising anti-tumor activities for immunotherapy-naive patients.⁵⁶ Thus, IER5L may function as an immunotherapeutic biomarker for NSCLC.

There are some highlights in our study. First, we have reported increased IER5L mRNA expression in multiple cancers, including NSCLC using multiple databases. Second, we performed comprehensive and in-depth bioinformatic analyses to explore the prognostic and therapeutic values of IER5L in NSCLC. Finally, qPCR validated the increased levels of IER5L mRNA expression in NSCLC. However, some limitations cannot be ignored in the current study. The current study is mainly based on bioinformatics analysis, and more experiments need to be performed to validate the function of IER5L in NSCLC. Moreover, the specific role of IER5L in cancer cell proliferation, migration, and invasion was not investigated in this study.

Conclusions

Pan-cancer analysis displayed that levels of IER5L expression were upregulated in multiple cancers and were associated with disease prognosis and progression, including NSCLC. An EGR1-hsa-miR-8075-IER5L network was constructed for NSCLC. DEGs are significantly enriched in the biological process of action potential, sodium ion transport, and regulation of membrane potential. Increased IER5L expression was primarily enriched in pathways of cell cycle, NOTCH signaling, and oxidative phosphorylation and was correlated with increased regulatory T cells and neutrophils, elevated levels of immune molecules, and high TMB. Thus, IER5L is a prognostic biomarker in multiple cancers and may correlate with immunotherapeutic response in NSCLC.