Histone demethylase JMJD2D emerges as a novel prognostic biomarker and exhibits correlation with immune infiltration in lung adenocarcinoma

Article type: Research Article

Keywords: JMJD2D, Lung adenocarcinoma (LUAD), Overall survival, Immune cell infiltration, Tumor mutation burden (TMB)

Authors: Bona Liu, Hui Hao, Zhen Wang, Yingchun Li, Cheng Du, Jian Ming, Shuang Zhang, Lin Zhou, Dazhi Liu

Affiliations: Department of Oncology, General Hospital of Northern Theater Command, Shenyang, China; https://ror.org/027hqk105grid.477849.1Department of Medical Oncology, Hebei Cangzhou People’s Hospital, Cangzhou, China; https://ror.org/04c8eg608grid.411971.b0000 0000 9558 1426Department of Oncology, General Hospital of Northern Theater Command, Dalian Medical University, Shenyang, China; Department of Pathology, General Hospital of Northern Theater Command, Shenyang, China; Department of Military Patient Management, General Hospital of Northern Theater Command, Shenyang, China; https://ror.org/02gxych78grid.411679.c0000 0004 0605 3373Department of Thoracic Surgery, Yuebei People’s Hospital Affiliated to Shantou University Medical College, Shaoguan, China; Department of Thoracic Surgery, General Hospital of Northern Theater Command, Shenyang, China

License: © The Author(s) 2025 CC BY 4.0 Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article’s Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.

Article links: DOI: 10.1007/s12672-025-02871-y  |  PubMed: 40504421  |  PMC: PMC12162437

Relevance: Moderate: mentioned 3+ times in text

Full text: PDF (5.4 MB)

Introduction

Lung cancer is a serious malignancy that is estimated to cause approximately 2.2 million new cases and 1.8 million deaths in 2020 [ref. 1]. Lung cancer is the leading cause of cancer-related morbidity and mortality in China in the past decades, resulting in a significant health, financial, and societal burden [ref. 2, ref. 3]. The overall five-year survival rates of patients with lung cancer are only 75.2%, 45.25%, 24%, and 8%, respectively, for stages I to IV [ref. 2]. According to histological classification, lung cancer is divided into small-cell lung carcinoma (SCLC) and non-small-cell lung carcinoma (NSCLC), Small-cell lung carcinoma (SCLC) and non-small-cell lung carcinoma (NSCLC) [ref. 4]. Approximately 85% of cases are NSCLC, which is further classified into lung squamous cell carcinoma (LUSC), lung adenocarcinoma (LUAD), and large-cell carcinoma (LCC) [ref. 5]. Although the introduction novel agents and the application of predictive biomarkers for patients with advanced NSCLC, the long-term survival outcomes remain poor [ref. 6–ref. 8]. Given the poor clinical outcomes of patients with NSCLC, it is crucial to identify additional potential biomarkers for the diagnosis and treatment of patients with NSCLC.

Epigenetic modifications, including DNA or RNA methylation, non-coding RNAs, chromatin remodeling, and histone modifications, are involved in tumorigenesis, progression, metastasis, and drug resistance, all of which are closely associated with tumor heterogeneity in lung cancer [ref. 9–ref. 11]. Of these epigenetic modifications, histone modifications play important roles in the regulation of chromatin state, gene expression, and other nuclear events in the various cancer types [ref. 12]. Histone modifications represent a crucial epigenetic mechanism, encompassing acetylation, phosphorylation, methylation, ADP ribosylation, ubiquitination, and citrullination, and more [ref. 13–ref. 16]. Numerous studies have indicated that histone methylation and demethylation are involved in tumor progression and drug resistance [ref. 17–ref. 19]. Histone demethylases contain two types of demethylases, including histone lysine demethylases and histone arginine demethylases, which consist of the lysine-specific demethylase (LSD) family, the Jumonji C (JmjC)-domain-containing demethylase (JMJD) family and the histone arginine demethylases [ref. 20]. JMJD2D, also named KDM4D, is a demethylase that targets histone H3 on lysines 9 and 36 and histone H1.4 on lysine 26 [ref. 21]. JMJD2D is widely overexpressed in various tumors. For example, upregulated JMJD2D is observed in human liver cancers and is associated with poor survival [ref. 22]. Moreover, it also reported that upregulated JMJD2D is found in colorectal cancer (CRC) and promotes tumor progression [ref. 23–ref. 25]. JMJD2D is also overexpressed in prostate tumors and promotes prostate tumorigenesis [ref. 26]. Upregulation of JMDJ2D in gastrointestinal stromal tumors (GST) is observed and promotes tumor progression [ref. 27]. JMJD2D is also overexpressed in acute myeloid leukemia (AML) and is linked to poor overall survival in AML patients [ref. 28]. A study has reported that JMJD2D expression is associated with the metastatic spread of lung carcinomas [ref. 29]. However, studies on JMJD2D in lung cancer, particularly LUAD, remain limited. Given the scarcity of previous research, we conducted a comprehensive bioinformatics analysis to investigate the expression and molecular characteristics of JMJD2D in LUAD, aiming to identify potential diagnostic and therapeutic targets for the disease.

In the present study, we explored the correlation of JMJD2D with tumor development, immune cell infiltration, drug resistance, and prognosis values based on the data of LUAD from The Cancer Genome Atlas (TCGA) and Gene Expression Omnibus (GEO). Furthermore, we also identified the JMJD2D-related risk signature and then developed a predictive nomogram for patients with LUAD.

Methods

Human tissue microarray analysis of JMJD2D

Human LUAD tumor tissues microarray (AF-LucSur2202) containing 80 patients with LUAD and adjacent normal tissues were obtained from Yunnan Yantai Biotechnology Co. Ltd. Each sample plot with a diameter of 1.5 mm and a thinness of 4-µm was prepared based on a standard method. Immunohistochemistry (IHC) staining was performed using the anti-JMJD2D antibody (1:300, 22591-1-AP, Proteintech, Wuhan). The IHC score of each sample, which ranged from 0 to 12, was calculated according to the intensity of the nucleic staining (no staining = 0; weak brown staining = 1, moderate brown staining = 2 and strong brown staining = 3) and the extent of stained cells (0–5%=0, 5–25%=1, 26–50%=2, 51–75%=3 and 76–100%=4) as previously described [ref. 30].

Data sources and processing

Transcriptome profiles, along with corresponding clinical information and somatic mutation data of LUAD patients, were acquired from the Cancer Genome Atlas (TCGA, https://portal.gdc.cancer.gov/). The transcriptome profiles of 535 LUAD tumors and 59 normal samples were retrieved using the University of California, Santa Crus (UCSC Xena, https://xena.ucsc.edu/), while the corresponding clinical information was obtained from the TCGA-LUAD project within the TCGA database through the cBioPortal platform at (https://www.cbioportal.org/). After quality assessments, data integration, batch correction, and normalization, a total of 479 LUAD tumors and 59 normal samples were included for subsequent analyses. Additionally, the GSE68465 dataset, which comprises the transcriptome profiles and corresponding clinical information of 462 LUAD patients, was obtained from the Gene Expression Omnibus (GEO, https://www.ncbi.nlm.nih.gov/geo/) and generated by the GPL96 platform. Following quality control and normalization, 442 LUAD specimens with complete clinical information were incorporated into this study.

Differential expression of JMJD2D across multiple cancer types

The transcriptome profiles of sixteen tumors were downloaded from the TCGA database, including bladder cancer (BLCA), breast cancer (BRCA), cervical cancer (CESC), colon cancer (COAD), esophageal cancer (ESCA), head and neck cancer (HNSC), kidney clear cell carcinoma (KIRC), liver cancer (LIHC), lung adenocarcinoma (LUAD), lung squamous cell carcinoma (LUSC), pancreatic cancer (PAAD), prostate cancer (PRAD), rectal cancer (READ), stomach cancer (STAD), thyroid cancer (THCA), and glioblastoma (GBM). The Wilcoxon rank-sum test was used to compare JMJD2D expression between tumor and normal samples across all cancer types, and the visualization of results was performed using the “ggpubr” package (version 0.4.0) in R.

Survival analysis

The transcriptome profiles of LUAD patients, along with corresponding clinical information such as survival time and survival status, were integrated using R. The optimal cutoff value of JMJD2D expression was determined using the “survminer” (version 0.4.9) package in R based on the maximally selected rank statistics from the “maxstat” R package, which is an outcome-oriented method, provides a cutoff value that corresponds to the most significant relationship with survival [ref. 31]. Then, patients were stratified into high-JMJD2D and low-JMJD2D groups according to the optimal cutoff value (Fig. 3C). Kaplan-Meier (KM) survival analysis and log-rank test were performed using the “survival” (version 3.2-7) package to compare overall survival (OS) differences between the high-JMJD2D and low-JMJD2D groups (Fig. 3D). Additionally, the Kruskal-Wallis test was employed to assess differences in age, gender, survival status, pathological stage, American Joint Committee on Cancer (AJCC) TNM stage, and other clinical parameters between the high-JMJD2D and low-JMJD2D groups (Fig. 3E). This comprehensive analysis aims to elucidate the potential prognostic significance of JMJD2D expression in LUAD patients.

Identification of the differentially expressed genes (DEGs)

“Limma” R package (version 3.46.0) [ref. 32] was utilized to identify Differentially Expressed Genes (DEGs) between tumor and normal specimens, as well as between high and low JMJD2D expression groups. DEGs were screened based on cutoff values of absolute log-fold change (FC) > 0.5 and adjusted p-value < 0.05 (Table S5). “ggplot2” R package (version 3.3.2) (Fig. 8A) [ref. 33] and “Pheatmap” R package (version 1.0.12) [ref. 34] were employed for visualization purposes (Fig. 8B). The DEGs, as well as the top 50 DEGs ranked by adjusted p-value, were visualized using these packages, providing insights into gene expression alterations associated with tumor status and JMJD2D expression levels. This approach enhances our understanding of the molecular mechanisms underlying LUAD progression and JMJD2D-associated pathways.

Gene set enrichment analysis (GSEA)

GSEA was employed to explore the immune-related pathways associated with JMJD2D expression. The “org.Hs.eg.db” R package (version 3.12.0) [ref. 35] was used to retrieve the Entrez-ID of each Differentially Expressed Gene (DEG). Subsequently, the “clusterProfiler” R package (version 3.18.1) [ref. 36] was employed to conduct the Gene Ontology (GO) analysis. The top 10 significantly enriched pathways, ranked by adjusted P-value were visualized using the “enrichplot” R package (version 1.10.2) [ref. 37]. This visualization enhances the interpretation of the GO analysis results, providing insight into the immune-related pathways that are potentially regulated by JMJD2D expression. This comprehensive approach contributes to our understanding of the immunological implications of JMJD2D in the context of LUAD.

Immune cell infiltration analysis

The “estimate” R package (version 1.0.13) was employed to calculate the immune score, stromal score, and estimate score for each sample, providing valuable insights into the tumor microenvironment. Besides, the “immunedeconv” R package (version 2.0.4) [ref. 38] employed the quantiseq algorithm to evaluate the distribution of infiltrated immune cells based on the transcriptome profiles. Moreover, CIBERSORT (version 1.03), a deconvolution algorithm with a leukocyte signature matrix (LM22) consisting of 547 genes distinguishing 22 immune populations [ref. 39], was utilized to assess the distribution of 22 types of infiltrating immune cells based on transcriptome profiles. Differences in cell fractions between the two groups were determined via the Wilcoxon rank-sum test and visualized using the “vioplot” R package (version 0.3.7). Moreover, differences in cell fractions were assessed using the Wilcoxon rank-sum test and visualized using the “ggpubr” R package (version 0.4.0) [ref. 40]. A comprehensive analysis was performed on 43 immune checkpoint-related genes (ICRGs) obtained from a previous study [ref. 41], wherein differences in the expression of 25 ICRGs were evaluated using the Wilcoxon rank-sum test. These analyses shed light on the immune landscape and potential immunotherapeutic targets in LUAD.

Prediction of the therapeutic responses to immunotherapy

The Tumor Immune Dysfunction and Exclusion (TIDE) algorithm is a computational method designed to model two primary mechanisms of tumor immune evasion: high infiltration of cytotoxic T lymphocytes (CTL) inducing T cell dysfunction within tumors, and low infiltration of CTL resulting in reduced T cell infiltration [ref. 42]. To evaluate the TIDE score for each sample, data from the TIDE database (http://tide.dfci.harvard.edu/) were utilized. Subsequently, differences in TIDE scores between different groups were determined using the Wilcoxon rank-sum test.

Tumor mutation burden (TMB) analysis

The Mutation Annotation Format (MAF) profiles of simple nucleotide variation were downloaded from the TCGA database (https://portal.gdc.cancer.gov/cart). The Tumor Mutational Burden (TMB) score for each sample was assessed and analyzed using the “maftools” R package (version 2.6.05) [ref. 43]. Subsequently, differences in TMB between the high-JMJD2D and low-JMJD2D groups were determined via the Wilcoxon rank-sum test and visualized using the “ggpubr” R package (version 0.4.0). The relationship between TMB and JMJD2D expression was determined by Spearman correlation analysis and visualized using the “ggExtra” R package (version 0.10.0). Additionally, the top 20 frequent mutation genes between high-JMJD2D and low-JMJD2D groups were visualized using the “oncoplot” function, providing insights into the mutational landscape associated with JMJD2D expression levels.

Drug sensitivity analysis

The information on 138 common targeted drugs was obtained from the Genomics of Drug Sensitivity in Cancer (GDSC, https://www.cancerrxgene.org/) database. Then, the “pRRophetic” R package (version 0.5) [ref. 44] was performed to calculate the Half-maximal inhibitory concentration (IC50) value of each sample based on the gene expression profiles. The lower IC50 value indicates greater sensitivity to the targeted drug. The correlation between JMJD2D expression and the IC50 value of drugs was visualized using the “ggplot2” R package. Furthermore, differences in IC50 values between the high-JMJD2D and low-JMJD2D expression groups were detected by the Wilcoxon rank-sum test. This comprehensive analysis provides insights into the potential associations between JMJD2D expression levels and drug sensitivity, aiding in the identification of personalized therapeutic strategies for LUAD patients.

Identification and validation of the JMJD2D-related prognosis signature in LUAD

Differentially expressed JMJD2D-related genes (DE-JMJD2D-RGs) were obtained using the “venn” R package (version 1.11) by intersecting the DEGs identified from tumor compared with normal samples and from the high-JMJD2D compared with low-JMJD2D groups, respectively (Fig. 8C). Subsequently, univariate Cox analysis was performed using the “survival” R package (version 3.2-7) to identify the prognostic DE-JMJD2D-RGs, and the results were visualized using the “forestplot” R package (version 2.0.1). The DE-JMJD2D-RGs with p-value < 0.05 were then incorporated into a least absolute shrinkage and selection operator (LASSO) regression model to prevent overfitting using the “glmnet” package in R (version 4.1-1). A prognostic signature comprising eight DE-JMJD2D-RGs was constructed through LASSO regression analysis. The JMJD2D-related risk score was calculated according to the formula, risk score = \(\:\sum\:_{i=1}^{n}coef\:\left(genei\right)*expr\:\left(genei\right)\). Here, coef represented the coefficient of each gene, and expr represented the expression level of each gene. Then, all patients of the TCGA-LUAD cohort were divided into high-risk and low-risk groups based on the median value of the risk score (Fig. 8H). The overall survival (OS) between high-risk and low-risk groups was conducted using the “survminer” package in R (version 0.4.8) (Fig. 8I). Furthermore, Kaplan-Meier (KM) analysis based on the risk score and receiver operating characteristic (ROC) curve was drawn using the “survivalROC” R package (version 1.0.3) to evaluate the performance of the risk model (Fig. 8J). Additionally, the correlation between risk score and clinical information, such as age, gender, survival status, pathological stage, and AJCC TNM staging, was determined using the Wilcoxon rank-sum test or the Kruskal-Wallis test ((Fig. 8G). Moreover, the GSE68465 dataset served as an external validation set to validate the robustness of the risk model (Fig. 8K-M). This comprehensive analysis aims to develop and validate a prognostic model based on JMJD2D-related risk scores, providing valuable insights into LUAD patient prognosis and potential clinical applications.

Construction of a nomogram

The JMJD2D-related risk score and clinical information were incorporated into univariate and multivariate Cox analyses to identify independent prognostic factors for LUAD. Factors with a p-value < 0.05 were selected as independent factors. Subsequently, these factors were utilized to construct a nomogram using the “rms” package in R (version 6.2.0). The discrimination ability of the nomogram was evaluated using the Harrell concordance index (C-index), which assesses the model’s ability to distinguish between patients with different survival outcomes. Additionally, calibration curves were generated to assess the agreement between the predicted and observed probabilities of survival. This integrated approach aims to develop a reliable prognostic tool for LUAD patients, facilitating personalized treatment decisions and improving clinical outcomes.

Results

Identification of JMJD2D expression across various Cancer types

The study design is illustrated in Fig. 1. Differential expression analysis of JMJD2D across sixteen tumor types, as depicted in Fig. 2, revealed that JMJD2D expression was significantly elevated in certain cancers compared to normal tissues. Specifically, high JMJD2D expression was observed in BLCA, COAD, HNSC, LUAD, LUSC, and STAD. Conversely, whereas low JMJD2D expression was noted in BRCA, KIRC, and THCA. These findings suggest that JMJD2D may exhibit dual roles in tumorigenesis, with its expression varying across different cancer types. Notably, JMJD2D appears to function as an oncogene in lung cancer, particularly in LUAD and LUSC.

High expression of JMJD2D associated with poor overall survival of LUAD patients

To investigate the role of JMJD2D in LUAD, IHC staining was performed to detect the expression of JMJD2D in LUAD specimens. As shown in Fig. 3A-B, JMJD2D exhibited higher expression levels in LUAD tumor samples compared to adjacent normal tissues. Given this observation, we focused on the role of JMJD2D in the LUAD in the subsequent analyses. For our study, a total of 479 LUAD patients with complete clinical information were included, source from the TCGA database. The optimal cutoff value of JMJD2D expression was calculated based on the maximally selected rank statistics. Based on this analysis, the cohort of 479 patients was divided into high-JMJD2D (n = 207) and low-JMJD2D (n = 272) groups, with the optimal cutoff value set at 20.98 (Fig. 3C). Notably, high expression of JMJD2D was associated with poorer survival outcomes compared to low expression, as evidenced by a significant log-rank test result (P = 0.024) (Fig. 3D). Furthermore, correlation analyses between JMJD2D expression and clinical features revealed a significant association with the T stage (Fig. 3E). These findings underscore the potential prognostic significance of JMJD2D in LUAD and its potential role in tumor progression.

Correlation of low JMJD2D expression with activation of Immune-Related pathways

Given the significant differences in prognosis associated with varying JMJD2D expression levels, we conducted differential gene expression analysis between the high-JMJD2D and low-JMJD2D groups. As shown in Fig. 4A, our analysis identified a total of 109 DEGs, comprising 46 upregulated and 63 downregulated DEGs, meeting the criteria of |log FC| > 0.5 and P < 0.05 (Table S1). Furthermore, the top 50 DEGs were ranked based on their P-values, as illustrated in Fig. 4B. Additionally, GSEA results revealed that the DEGs predominantly enriched in immune-related processes. Specifically, pathways associated with adaptive and humoral immune responses, positive regulation of immune responses, activation of the immune response, and immune response-regulating signaling pathways were significantly enriched in the low-JMJD2D group (Fig. 4C, Table S2). These findings suggest a potential link between JMJD2D expression levels and immune-related processes in LUAD, highlighting its potential role in tumor immunity.

Correlation of JMJD2D expression with immune cell infiltration

Given the significant association of JMJD2D with immune-related biological processes, we further investigated its relationship with tumor immune cell infiltration. As shown in Fig. 5A, various immune cell proportions, including B cells, macrophages (M1/M2), monocytes, neutrophils, NK cells, CD4 + and CD8 + T cells, regulatory T cells (Tregs), and myeloid dendritic cells, exhibited significantly higher levels compared to other cell types. Further analysis of immune cell distribution between the high-JMJD2D and low-JMJD2D groups revealed notable differences in monocytes, neutrophils, and Tregs (Fig. 5B). Specifically, Tregs were significantly increased, while monocytes and neutrophils were decreased in the high-JMJD2D group compared to the low-JMJD2D group. Moreover, we explored the expression of 25 Immune Checkpoint-Related Genes (ICRGs) between the high-JMJD2D and low-JMJD2D groups. The results unveiled significantly higher CD40 expression and lower ICOSLG expression in the high-JMJD2D group compared to the low-JMJD2D group (Fig. 5C). These findings provide further insight into the immune landscape of LUAD and underscore the potential role of JMJD2D in modulating immune cell infiltration.

High JMJD2D expression correlates with high tumor mutation burden (TMB)

To investigate the relationship between JMJD2D expression and somatic mutations, we analyzed mutation data from 479 LUAD patients obtained from the TCGA database. As shown in Fig. 6A, the top 20 mutated genes between high-JMJD2D and low-JMJD2D groups, highlighting frequent genetic alterations. These genes included TP53 (42%), TTN (37%), MUC16 (34%), CSMD3 (34%), RYR2 (29%), LRP1B (28%), ZFHX4 (27%), USH2A (25%), KRAS (23%), FLG (20%), XIRP2 (19%), SPTA1 (19%), COL11A1 (17%), NAV3 (16%), FAT3 (16%), ANK2 (16%), ZNF536 (16%), CSMD1 (15%), MUC17 (16%), and APOB (15%). Moreover, our analysis revealed that TMB was significantly higher in the high-JMJD2D group compared to the low-JMJD2D group (Fig. 6B), Importantly, a weak positive correlation was observed between JMJD2D expression and TMB (Fig. 6C). These findings suggest that JMJD2D may be associated with increased genomic instability in LUAD, further underscoring its potential role in tumor progression.

Correlation of JMJD2D with therapeutic responses

To further elucidate the relationship between JMJD2D expression and therapeutic responses in LUAD, we assessed drug sensitivity by calculating IC50 values for various chemotherapeutic agents using data from the GDSC database (Table S3). Correlation analysis between JMJD2D expression and IC50 values revealed distinct drug-specific associations. Specifically, JMJD2D expression was positively correlated with the sensitivity nine drugs (BMS.708163, Roscovitine, Pyrimethamine, AZ628, Bicalutamide, PD.0325901, RDEA119, JNK.9 L, BMS.536924), whereas a negative correlation was observed for eleven drugs (KU.55933, AP.24534, IPA.3, Lenalidomide, ABT.263, CMK, SL.0101.1, Pazopanib, ATRA, Bosutinib, JNK. Inhibitor. VIII) (Fig. 7A, Table S4). Further analysis focused on six drugs with a correlation coefficient exceeding 0.2, revealing significant differences in IC50 values between the high-JMJD2D and low-JMJD2D groups. Specifically, the high-JMJD2D group exhibited higher IC50 values for BMS.708,163, Roscovitine, and Pyrimethamine compared to the low-JMJD2D group (Fig. 7B-D). Conversely, the high-JMJD2D group displayed lower IC50 values for ATRA, Bosutinib, and JNK. Inhibitor. VIII compared to the low-JMJD2D group (Fig. 7E-G). These findings suggest that JMJD2D may play a crucial role in modulating drug response in LUAD, providing potential insights for precision therapy.

Development of a JMJD2D-Related risk model for LUAD

Considering the multifaceted influences of JMJD2D on prognosis, tumor immune-related processes, TMB, and drug sensitivity in LUAD, we endeavored to develop a comprehensive JMJD2D-related prognostic signature. As shown in Fig. 8A-B, a total of 5374 DEGs (2951 upregulated and 2423 downregulated) between tumor and normal were identified with |log FC| > 0.5 and adjusted P-value < 0.05 (Table S5). Subsequently, we overlapped these DEGs with the 109 JMJD2D-related DEGs, resulting in the selection of 77 DE-JMJD2D-RGs for constructing the prognostic signature (Fig. 8C). Utilizing Kaplan-Meier survival analysis, we screened nine prognostic DE-JMJD2D-RGs with a P-value < 0.05 (Fig. 8D). Next, using LASSO regression analysis, we further refined the prognostic genes, ultimately identifying eight genes (NR0B2, IGF2BP3, TFF1, CFTR, ELAVL2, GSTA1, CTNND2, HGD) for inclusion in the prognostic signature (Fig. 8E-F). The risk score was calculated using the following formula: RiskScore = NR0B2*-0.005066 + IGF2BP3*0.0039992 + TFF1*0.0033018 + CFTR*(-0.0013022) + ELAVL2*0.015465 + GSTA1*(-0.0021139) + CTNND2*(-0.0065844) + HGD*0.002528 (Table 1). Subsequently, 479 patients from the TCGA-LUAD cohort were stratified into high- and low-risk groups using the median risk score of 0.944 (Fig. 8H). Kaplan-Meier survival analysis revealed poorer survival outcomes in the high-risk group compared to the low-risk group (Fig. 8I). Additionally, the predictive accuracy of the prognostic signature was evaluated using ROC curves for 1-, 3-, and 5-year OS, yielding AUC values of 0.662, 0.650, and 0.623, respectively (Fig. 8J). Furthermore, external validation was performed using the GSE68465 dataset, where 442 patients were divided into high- and low-risk groups based on the median risk score of 1.25 (Fig. 8K). Consistent with the TCGA-LUAD cohort, the high-risk group exhibited worse survival outcomes than the low-risk group (Fig. 8L). The AUC values of the ROC curves for 1-, 3-, and 5-year OS were 0.678, 0.632, and 0.604, respectively (Fig. 8M). Additionally, we investigated the correlation between the risk score and clinical features, revealing differences in gender, survival status, overall survival time, pathologic stages, and AJCC stages between the high-risk and low-risk groups (Fig. 8G; Table 2). Overall, the JMJD2D-related prognostic signature demonstrates promising predictive performance and clinical utility in prognostic assessment for LUAD patients.

Construction of a predictive nomogram

The risk score, age, gender, pathologic stages, and AJCC stages, were subjected to univariate and multivariate Cox regression analyses to identify the independent factors for LUAD. The results indicated that risk score was only independent factor for prognosis in LUAD (Fig. 9A-B). Subsequently, a nomogram was constructed by integrating the risk score and relevant risk factors (Fig. 9C). The C-index was determined to be 0.6251, indicating moderate discriminative ability. Additionally, calibration curves for survival probability at 1-, 3-, and 5-years demonstrated the accurate predictive performance of the nomogram across both short- and long-term intervals (Fig. 9D). Overall, these findings underscore the utility of the risk score as an independent prognostic factor for LUAD, and the nomogram provides a valuable tool for individualized prognostic assessment in clinical practice.

Correlation of risk score with immune characteristics

Given the association of the risk score with poor clinical outcomes of LUAD patients, we investigated risk score correlation with immune characteristics. Initially, we utilized the CIBERSORT algorithm to determine the distribution of 22 types of immune cells in tumor samples (Fig. 10A). Subsequently, we compared the immune cell distribution between high- and low-risk groups, revealing a significant decrease in the fraction of activated CD4 memory T cells, Macrophages M0, and activated dendritic cells (DCs), alongside increased resting DCs, resting mast cells, and neutrophils in the high-risk group compared to the low-risk group (Fig. 10B). Moreover, we observed higher TIDE score and exclusion score, but lower dysfunction scores in the high-risk group relative to low-risk group (Fig. 10C-D). Besides, higher expression of CD274 was detected in the high-risk group compared to the low-risk group (Fig. 10E). These findings suggest that high-risk scores may be linked to immune suppression in LUAD, potentially contributing to a less favorable prognosis.

Discussion

With the development of next-generation sequencing (NGS), numerous crucial molecular biomarkers have been identified in the diagnosis and treatment of NSCLC, including alterations in EGFR, BRAF, MET, ALK, ROS1, RET, and NTRK, along with the corresponding kinase inhibitors [ref. 8, ref. 29, ref. 45]. Despite significant progress has been made in the diagnosis and treatment of NSCLC over recent decades, the overall survival rates for patients with NSCLC remain suboptimal. Consequently, the search for more sensitive biomarkers to improve the diagnosis and treatment of NSCLC continues to be a critical challenge in enhancing patient survival outcomes.

In recent years, accumulating evidence has highlighted the pivotal role of epigenetic modifications in cancer initiation and progression [ref. 46]. Reversible histone methylation, is a crucial process involved in tumor onset and advancement [ref. 47–ref. 49]. Numerous studies have elucidated the significant involvement of JMJD2D in regulating tumor initiation, progression, and drug resistance, positioning it as a promising therapeutic target across various cancers [ref. 50–ref. 52]. In our present study, we observed significant overexpression of JMJD2D in BLCA, COAD, HNSC, LUAD, LUSC, and STAD. Elevated JMJD2D expression was closely associated with poorer survival outcomes in LUAD patients, inhibition of immune-related pathways, heightened TMB, and resistance to certain anti-tumor drugs. Previous research has indicated that JMJD2D facilitates CRC immune evasion by inducing PD-L1 overexpression [ref. 24]. Here, we not only found the high-JMJD2D-induced inhibition of immune-related pathways but also noted an upregulation of CD40 expression, suggesting that JMJD2D may contribute to immune escape in LUAD, thereby accelerating tumor progression. Additionally, JMJD2D has been implicated in DNA damage response (DDR), linking it to cancer-relevant epigenetic regulation and genome stability [ref. 53–ref. 56]. Similarly, we also found some mutated genes associated with JMJD2D and observed a weak correlation between high JMJD2D expression and elevated TMB. Furthermore, we found that JMJD2D conferred resistance to BMS.7.8163, Roscovitine, and Pyrimethamine, while displaying sensitivity to ATRA, Bosutinib, and JNK. Inhibitor. VIII. These findings collectively suggest that JMJD2D may serve as a potential biomarker for guiding anti-tumor drug selection.

Based on the connection between JMJD2D and tumorigenesis, particularly in LUAD progression, we developed a risk score model using JMJD2D-related genes, including NR0B2, IGF2BP3, TFF1, CFTR, ELAVL2, GSTA1, CTNND2, and HGD. Among these genes, NR0B2, CFTR, GSTA1, and CTNND2 were found to exert a protective role, whereas IGF2BP3, TFF1, ELAVL2, and HGD were associated with higher-risk scores. NR0B2 has been identified as a favorite prognosis factor in liver cancer [ref. 57]. CFTR, a cAMP-activated chloride channel, regulates fluid homeostasis and transports relevant substrates [ref. 58]. Alterations in CFTR are implicated in various tumors, with dysregulation linked to more aggressive lung cancers [ref. 59, ref. 60]. GSTA1, a component of glutathione S-transferases (GSTs), participates in carcinogen detoxification and plays a crucial role in carcinogenesis [ref. 61]. Clinical trials have shown that genetic polymorphisms in GSTA1 are linked to lung cancer occurrence [ref. 62]. CTNND2 is overexpressed in lung cancer, promoting tumorigenesis and metastasis [ref. 63, ref. 64].

In contrast, IGF2BP3, an N6-Methyladenosine (m6A) reader, acts as an oncogene gene in lung cancer [ref. 65, ref. 66]. TFF1 overexpression has been observed in various cancers, including breast, colonic, pancreatic, ovary, prostate, thyroid, and lung cancers, is associated with chemoresistance, metastasis, and poor prognosis [ref. 67–ref. 70]. ELAVL2, although its role in lung cancer is less explored, has been reported as an oncogene promotes ovarian cancer cells resistant to paclitaxel [ref. 71]. HGD, another relatively unexplored gene in cancers, is involved in cholangiocarcinoma genesis [ref. 72], and has been identified as a candidate co-expressed with B-type Raf kinase (BRAF) V600E mutation in papillary thyroid carcinoma [ref. 73]. Collectively, these eight genes hold promise as potential biomarkers for molecular diagnosis and treatment strategies in LUAD. Further research on their roles and mechanisms in LUAD pathogenesis may provide valuable insights into the disease’s management.

Based on the risk score, we also explored the association of immune characteristics with the risk score in LUAD. The findings revealed a decrease in the proportions of activated memory CD4 + T cells, Macrophages M0, and activated DCs, along with an increase in resting DCs, resting mast cells, and neutrophils in samples with a high-risk score. Moreover, higher TIDE score and exclusion score, as well as overexpression of CD274, were associated with a high-risk score, suggesting potential involvement of high-risk scores in immunosuppression.

Despite the comprehensive exploration of the role and regulatory mechanisms of JMJD2D in LUAD in our study, there are still some limitations that need to be addressed. First, a larger sample cohort is needed to further validate our findings and enhance the statistical power of the analysis. Second, experimental studies are required to confirm the regulatory mechanisms proposed in this study, providing deeper insights into the functional roles of JMJD2D in LUAD.

Conclusion

In summary, our study elucidated the oncogenic role of JMJD2D in LUAD, highlighting its involvement in tumorigenesis, progression, antitumor resistance, and poor prognosis. Additionally, we identified a JMJD2D-related prognostic signature and developed both a risk model and predictive model for LUAD, offering potential tools for personalized prognosis and therapeutic strategies in LUAD patients.

References

1. H Sung. Global cancer statistics 2020: GLOBOCAN estimates of incidence and mortality worldwide for 36 cancers in 185 countries. CA Cancer J Clin, 2021. [DOI | PubMed]
2. YH Luo. Lung cancer in Republic of China. J Thorac Oncol, 2021. [DOI | PubMed]
3. S Gao. Lung cancer in People’s Republic of China. J Thorac Oncol, 2020. [DOI | PubMed]
4. T Sher, GK Dy, AA Adjei. Small cell lung cancer. Mayo Clin Proc, 2008. [DOI | PubMed]
5. TV Denisenko, IN Budkevich, B Zhivotovsky. Cell death-based treatment of lung adenocarcinoma. Cell Death Dis, 2018. [DOI | PubMed]
6. D Planchard. Metastatic non-small cell lung cancer: ESMO clinical practice guidelines for diagnosis, treatment and follow-up. Ann Oncol, 2018. [DOI | PubMed]
7. EB Garon. Five-Year overall survival for patients with advanced non–small-cell lung cancer treated with pembrolizumab: results from the phase I KEYNOTE-001 study. J Clin Oncol, 2019. [DOI | PubMed]
8. RS Herbst, D Morgensztern, C Boshoff. The biology and management of non-small cell lung cancer. Nature, 2018. [DOI | PubMed]
9. 9.Chao YL, Pecot CV. Targeting epigenetics in lung Cancer. Cold Spring Harb Perspect Med, 2021. 11(6).
10. H Luo. Targeting the epigenetic processes to enhance antitumor immunity in small cell lung cancer. Semin Cancer Biol, 2022. [DOI | PubMed]
11. HW Tan. Epigenetic regulation of angiogenesis in lung cancer. J Cell Physiol, 2021. [DOI | PubMed]
12. ME Neganova. Histone modifications in epigenetic regulation of cancer: perspectives and achieved progress. Semin Cancer Biol, 2022. [DOI | PubMed]
13. S Messner, MO Hottiger. Histone ADP-ribosylation in DNA repair, replication and transcription. Trends Cell Biol, 2011. [DOI | PubMed]
14. D Zhu, Y Zhang, S Wang. Histone citrullination: a new target for tumors. Mol Cancer, 2021. [DOI | PubMed]
15. TJ Jarome. Ubiquitination of histone H2B by proteasome subunit RPT6 controls histone methylation chromatin dynamics during memory formation. Biol Psychiatry, 2021. [DOI | PubMed]
16. N Wang, T Ma, B Yu. Targeting epigenetic regulators to overcome drug resistance in cancers. Signal Transduct Target Ther, 2023. [DOI | PubMed]
17. 17.Liu Q et al. TdIF1-LSD1 Axis regulates Epithelial-Mesenchymal transition and metastasis via histone demethylation of E-Cadherin promoter in lung Cancer. Int J Mol Sci, 2021. 23(1).
18. S Chang, S Yim, H Park. The cancer driver genes IDH1/2, JARID1C/ KDM5C, and UTX/ KDM6A: crosstalk between histone demethylation and hypoxic reprogramming in cancer metabolism. Exp Mol Med, 2019. [DOI | PubMed]
19. Y Chen. Lung Cancer therapy targeting histone methylation: opportunities and challenges. Comput Struct Biotechnol J, 2018. [DOI | PubMed]
20. 20.Shi YG, Tsukada Y. The discovery of histone demethylases. Cold Spring Harb Perspect Biol, 2013. 5(9).
21. WL Berry, R Janknecht. KDM4/JMJD2 histone demethylases: epigenetic regulators in cancer cells. Cancer Res, 2013. [DOI | PubMed]
22. M Li. Demethylase-independent function of JMJD2D as a novel antagonist of p53 to promote liver Cancer initiation and progression. Theranostics, 2020. [DOI | PubMed]
23. K Peng. Histone demethylase JMJD2D interacts with β-Catenin to induce transcription and activate colorectal Cancer cell proliferation and tumor growth in mice. Gastroenterology, 2019. [DOI | PubMed]
24. Q Chen. Demethylase JMJD2D induces PD-L1 expression to promote colorectal cancer immune escape by enhancing IFNGR1-STAT3-IRF1 signaling. Oncogene, 2022. [DOI | PubMed]
25. K Peng. Histone demethylase JMJD2D activates HIF1 signaling pathway via multiple mechanisms to promote colorectal cancer glycolysis and progression. Oncogene, 2020. [DOI | PubMed]
26. R Gu. Methylation of the epigenetic JMJD2D protein by SET7/9 promotes prostate tumorigenesis. Front Oncol, 2023. [DOI | PubMed]
27. F Hu. Histone demethylase KDM4D promotes gastrointestinal stromal tumor progression through HIF1β/VEGFA signalling. Mol Cancer, 2018. [DOI | PubMed]
28. W Wu, X Cao, L Mo. Overexpression of KDM4D promotes acute myeloid leukemia cell development by activating MCL-1. Am J Transl Res, 2021. [PubMed]
29. Y Soini, VM Kosma, R Pirinen. KDM4A, KDM4B and KDM4C in non-small cell lung cancer. Int J Clin Exp Pathol, 2015. [PubMed]
30. P Xie. The covalent modifier Nedd8 is critical for the activation of Smurf1 ubiquitin ligase in tumorigenesis. Nat Commun, 2014. [DOI | PubMed]
31. 31.Kassambara A. surv_cutpoint: Determine the Optimal Cutpoint for Continuous Variables in survminer. https://search.r-project.org/CRAN/refmans/survminer/html/surv_cutpoint.html.
32. ME Ritchie. Limma powers differential expression analyses for RNA-sequencing and microarray studies. Nucleic Acids Res, 2015. [DOI | PubMed]
33. 33.Wickham HJPa. ggplot2: elegant graphics for data analysis Springer-Verlag New York; 2009. 2016.
34. 34.Kolde R. pheatmap: Pretty Heatmaps. R package version 1.0. 12. 2019.
35. 35.Carlson MJBL. United Kingdom: Genome Biology, org. Mm. eg. db: Genome wide annotation for Mouse. R package version 3.2. 3. 2019.
36. 36.Yu G et al. ClusterProfiler: an R package for comparing biological themes among gene clusters. 2012. 16(5): pp. 284–7.
37. 37.Yu GJMTNA. enrichplot: Visualization of functional enrichment result. R package version 1.10. 2. 2021.
38. G Sturm, F Finotello, M List. Immunedeconv: an R package for unified access to computational methods for estimating immune cell fractions from bulk RNA-Sequencing data. Methods Mol Biol, 2020. [DOI | PubMed]
39. B Chen. Profiling tumor infiltrating immune cells with CIBERSORT. Methods Mol Biol, 2018. [DOI | PubMed]
40. 40.Kassambara A. ggpubr:‘ggplot2’based publication ready plots. R package version 0.4. 0. 2020.
41. Y Xie. Identification of Lactate-Related gene signature for prediction of progression and immunotherapeutic response in skin cutaneous melanoma. Front Oncol, 2022. [DOI | PubMed]
42. P Jiang. Signatures of T cell dysfunction and exclusion predict cancer immunotherapy response. Nat Med, 2018. [DOI | PubMed]
43. A Mayakonda. Maftools: efficient and comprehensive analysis of somatic variants in cancer. Genome Res, 2018. [DOI | PubMed]
44. P Geeleher, N Cox, RS Huang. pRRophetic: an R package for prediction of clinical chemotherapeutic response from tumor gene expression levels. PLoS ONE, 2014. [DOI | PubMed]
45. EN Imyanitov, AG Iyevleva, EV Levchenko. Molecular testing and targeted therapy for non-small cell lung cancer: current status and perspectives. Crit Rev Oncol Hematol, 2021. [DOI | PubMed]
46. L Morera, M Lübbert, M Jung. Targeting histone methyltransferases and demethylases in clinical trials for cancer therapy. Clin Epigenetics, 2016. [DOI | PubMed]
47. J Marsolier. H3K27me3 conditions chemotolerance in triple-negative breast cancer. Nat Genet, 2022. [DOI | PubMed]
48. CA Zahnow. Inhibitors of DNA methylation, histone deacetylation, and histone demethylation: a perfect combination for Cancer therapy. Adv Cancer Res, 2016. [DOI | PubMed]
49. X Kong. Defining UHRF1 domains that support maintenance of human Colon cancer DNA methylation and oncogenic properties. Cancer Cell, 2019. [DOI | PubMed]
50. Z Fang. Discovery of a potent and selective inhibitor of histone lysine demethylase KDM4D. Eur J Med Chem, 2021. [DOI | PubMed]
51. DH Lee. Advances in histone demethylase KDM4 as cancer therapeutic targets. Faseb J, 2020. [DOI | PubMed]
52. YW Chin, SY Han. KDM4 histone demethylase inhibitors for anti-cancer agents: a patent review. Expert Opin Ther Pat, 2015. [DOI | PubMed]
53. LC Young, DW McDonald, MJ Hendzel. Kdm4b histone demethylase is a DNA damage response protein and confers a survival advantage following γ-irradiation. J Biol Chem, 2013. [DOI | PubMed]
54. H Khoury-Haddad. The emerging role of lysine demethylases in DNA damage response: dissecting the recruitment mode of KDM4D/JMJD2D to DNA damage sites. Cell Cycle, 2015. [DOI | PubMed]
55. H Khoury-Haddad. PARP1-dependent recruitment of KDM4D histone demethylase to DNA damage sites promotes double-strand break repair. Proc Natl Acad Sci U S A, 2014. [DOI | PubMed]
56. R Wu. H3K9me3 demethylase Kdm4d facilitates the formation of pre-initiative complex and regulates DNA replication. Nucleic Acids Res, 2017. [DOI | PubMed]
57. R Zhu. The orphan nuclear receptor gene NR0B2 is a favorite prognosis factor modulated by multiple cellular signal pathways in human liver cancers. Front Oncol, 2021. [DOI | PubMed]
58. Q Wang. Analysis of CFTR gene expression as an immunological and prognostic biomarker in pan-cancers. Comput Biol Med, 2022. [DOI | PubMed]
59. X Shi. Integrative pan cancer analysis reveals the importance of CFTR in lung adenocarcinoma prognosis. Genomics, 2022. [DOI | PubMed]
60. K Vekens. Lung cancer in a CF patient: combination of bad luck or is there more to say?. Acta Clin Belg, 2021. [DOI | PubMed]
61. Q Deng. Polymorphisms of GSTA1 contribute to elevated cancer risk: evidence from 15 studies. J Buon, 2015. [PubMed]
62. M Sørensen. Genetic polymorphisms in CYP1B1, GSTA1, NQO1 and NAT2 and the risk of lung cancer. Cancer Lett, 2005. [DOI | PubMed]
63. F Huang. δ-Catenin promotes tumorigenesis and metastasis of lung adenocarcinoma. Oncol Rep, 2018. [PubMed]
64. Q Lu. δ-Catenin dysregulation in cancer: interactions with E-cadherin and beyond. J Pathol, 2010. [DOI | PubMed]
65. A Yan. IGF2BP3 worsens lung Cancer through modifying long Non-coding RNA CERS6-AS1/microRNA-1202 Axis. Curr Med Chem, 2023. [DOI | PubMed]
66. Z Liu. N(6)-methyladenosine-modified circIGF2BP3 inhibits CD8(+) T-cell responses to facilitate tumor immune evasion by promoting the deubiquitination of PD-L1 in non-small cell lung cancer. Mol Cancer, 2021. [DOI | PubMed]
67. S Zhao, Y Ma, X Huang. Trefoil factor 1 elevates the malignant phenotype of mucinous ovarian cancer cell through Wnt/β-catenin signaling. Int J Clin Exp Pathol, 2015. [PubMed]
68. T Arumugam. Trefoil factor 1 stimulates both pancreatic cancer and stellate cells and increases metastasis. Pancreas, 2011. [DOI | PubMed]
69. A Yusufu. TFF3 and TFF1 expression levels are elevated in colorectal cancer and promote the malignant behavior of colon cancer by activating the EMT process. Int J Oncol, 2019. [PubMed]
70. K Minegishi. TFF-1 functions to suppress multiple phenotypes associated with lung cancer progression. Onco Targets Ther, 2021. [DOI | PubMed]
71. GB Park. Glucose deprivation enhances resistance to paclitaxel via ELAVL2/4-mediated modification of glycolysis in ovarian cancer cells. Anticancer Drugs, 2022. [DOI | PubMed]
72. R Aukkanimart. Altered expression of oxidative metabolism related genes in cholangiocarcinomas. Asian Pac J Cancer Prev, 2015. [DOI | PubMed]
73. X Yu. Key candidate genes associated with BRAF(V600E) in papillary thyroid carcinoma on microarray analysis. J Cell Physiol, 2019. [DOI | PubMed]