This study was approved by the Academic Committee of Cancer Hospital of Chinese Academy of Medical Sciences and carried out in strict accordance with the principles of the Helsinki Declaration. The data in this study were retrieved from online databases, all necessary written informed consent forms were obtained, and no human or animal experiments were performed.

In this study, R (version3.6.3) was used to analyze the expression level (https://xenabrowser.net/datapages/) of the KDM5 family in tumor and paracancerous tissues in the PAAD TCGA cohort and normal pancreatic tissue in the GTEx database, and visualized with ggplot2 (3.3.3). The statistical significance of the differential expression was evaluated by the Wilcoxon test. To analyze the correlation between the gene expression levels of KDM5 family members and clinical variables, we used the R (version 3.6.3) to analyze PAAD samples from the TCGA database (https://portal.gdc.cancer.gov/). The statistical significance of the differential expression was evaluated by Fisher’s test.

We used cBioPortal (https://www.cbioportal.org/) (Gao et al., 2013) to analyze the gene mutation of KDM5 family members in pancreatic cancer and further determined the correlation of these mutations with important clinicopathological factors. The statistical significance of the difference was evaluated by the chi-squared test. Furthermore, the log-rank test was used to evaluate the relationship between KDM5 family gene mutations and overall survival (OS) and disease-free survival (DFS) in patients with pancreatic cancer.

We used Kaplan–Meier Plotter (http://kmplot.com/analysis/) (Győrffy et al., 2013) to analyze the correlation between the expression of KDM5 family genes and OS and relapse-free survival (RFS) in pancreatic cancer. The cutoff for low expression versus high expression was set as the value automatically selected by the best cutoff model, and the array with deviation was excluded. The log rank test was used to calculate the p value (p < 0.05).

We used TIMER (http://timer.cistrome.org/) (Li et al., 2017) to analyze the correlations of KDM5 family gene expression in pancreatic cancer, and the statistical significance was evaluated by Spearman’s test. LinkedOmics (http://www.linkedomics.org/) (Vasaikar et al., 2018) was used to assess and draw a volcano plot of the correlations between KDM5 family members and 800 coexpressed genes in pancreatic cancer and a heatmap of the top 50 positively/negatively related genes. The genes with the strongest interaction with KDM5 family proteins were determined by combined score analysis via STRING (https://string-db.org/) (Szklarczyk et al., 2015), and the corresponding protein-protein interaction (PPI) network was constructed by STRING.

The first 10 motifs with positive and negative correlations with the expression of each member of the KDM5 family were functionally annotated (https://www.genecards.org/) by GeneCards.

We used Metascape (https://metascape.org) (Zhou et al., 2019) to visualize the enriched biological process (BP), cellular composition (CC), molecular function (MF) and KEGG pathway terms of the KDM5 family and its coexpressed genes. Furthermore, the R (version 3.6.3) and clusterProfiler package (version 3.14.3) were used to visualize the enriched BP, CC, MF, and KEGG pathway terms of the 9 genes with the strongest interaction with KDM5 family proteins.

We used TIMER to analyze the correlation between the expression level of the KDM5 family and the degree of lymphocyte infiltration. The statistical significance of the difference was evaluated by Spearman’s test. In addition, TISIDB (http://cis.hku.hk/TISIDB/) (Ru et al., 2019) was also used to analyze the correlation between KDM5 family expression and the expression of immune molecules in pancreatic cancer. The difference was evaluated by Spearman’s test. We used the GSCA Lite online tool (http://bioinfo.life.hust.edu.cn/web/GSCALite/) (Liu et al., 2018) to analyze the correlation between KDM5 family expression and sensitivity to current chemotherapeutic or targeted drugs for pancreatic cancer. The difference was evaluated by Spearman’s test.

We used the R (version 3.6.3) to analyze the expression level of the KDM5 family in tumor and paracancerous tissues in TCGA and normal tissues in GTEx. The analysis of gene expression levels showed that KDM5 family genes were upregulated or downregulated in tumors compared with normal or paracancerous tissues in each type of cancer. Compared with that in paracancerous tissues or normal tissues, the expression of KDM5A/B/C was significantly increased in 12 kinds of tumor tissues, such as PAAD, breast invasive carcinoma (BRCA) and stomach adenocarcinoma (STAD), but significantly downregulated in only adrenocortical carcinoma (ACC) and skin cutaneous melanoma (SKCM) (Figures 1A–C). In contrast, the expression of KDM5D was significantly decreased in 24 kinds of tumor tissues, such as PAAD, BRCA, and prostate adenocarcinoma (PRAD) (Figure 1D). In summary, our results show that KDM5A/B/C is significantly overexpressed in a variety of tumor tissues, including PAAD, while the expression of KDM5D is significantly downregulated, indicating that the expression of KDM5 family members is closely related to the occurrence and development of many kinds of human malignant tumors, including pancreatic cancer.

To further explore the mechanism of differential expression of the KDM5 family in pancreatic cancer, we used the cBioPortal online tool to analyze the gene mutations of the KD

M5 family. The KDM5 family had two or more types of gene variants in 17 samples (11%) from pancreatic cancer patients (Figure 2A). The most common variant was in the KDM5A/B gene (5%). The main types of variants were amplification and missense mutations. In addition, 2.7% of the variants occurred in KDM5C, and the most common type of variant was deep deletion, while no related mutations were detected in the KDM5D gene. Then, we comprehensively analyzed the clinical and pathological features of patients with KDM5 family gene mutations and those without mutations. The results showed that patients with a history of drinking and pancreatic cancer with mucinous adenocarcinoma were more likely to develop KDM5 family gene mutation (Figures 2B,C). Next, we analyzed the independent relationship between KDM5A/B/C gene mutation and clinicopathological features. The results showed that most of the patients with KDM5A gene mutation had a family history of breast cancer and melanoma or a pathological pancreatic cancer subtype of pancreatic mucinous adenocarcinoma (Figures 2D,E). Most of the patients with KDM5B gene mutation had a history of drinking (Figure 2F), and most of the patients with KDM5C gene mutation had a family history of breast cancer, a history of chronic pancreatitis and a higher pancreatic cancer tumor grade (Figures 2G–I). In the analysis of the prognostic characteristics of patients with KDM5 family gene mutation and those without mutation, it was not found that KDM5 family gene mutation had a significant effect on OS (Supplementary Figure S1A) and DFS (Supplementary Figure S1B) in patients with pancreatic cancer, which may have been related to the small sample size of patients. In summary, the frequency of KDM5A/B/C gene mutations was higher in patients with pancreatic cancer, while patients with a history of drinking, a family history of breast cancer, and a family history of melanoma or a history of chronic pancreatitis were more likely to develop KDM5A/B/C mutations. These mutations may lead to more malignant and higher-grade tumors. However, the effect of KDM5A/B/C gene mutation on the prognosis of patients with pancreatic cancer needs to be further clarified in studies with a large sample size.

To analyze the correlation between the gene expression levels of KDM5 family members and clinical variables, we used the basic R package (version 3.6.3) to analyze the samples of the PAAD cohort in the TCGA database. As shown in Table 1, high expression of KDM5A was significantly associated with poorer treatment efficacy (p = 0.037) and higher tumor grade (p = 0.019). High expression of KDM5B was significantly correlated with age ≤65 years old (p = 0.036), no history of diabetes (p = 0.038), higher tumor grade (p = 0.048), more advanced clinical stage (p = 0.004), and more advanced T stage (p = 0.020). High expression of KDM5C was significantly correlated with female sex (p < 0.001). High expression of KDM5D was significantly correlated with male sex (p < 0.001) and smoking history (p = 0.041). The results of this study suggest that the high expression levels of KDM5A and KDM5B can increase the degree of malignancy of PAAD and lead to a worse prognosis of patients, further confirming that members of KDM5 family may be oncogenes of PAAD. Next, we analyzed the correlation between the expression levels of KDM5 family genes and prognosis in patients with pancreatic cancer via Kaplan–-Meier Plotter. The results showed that upregulated expression of KDM5A was associated with shorter OS and RFS (Figures 3A,B). Downregulation of KDM5D expression was significantly correlated with shorter OS and RFS (Figures 3A,B). In addition, although there was a lack of consistent significance in the associations with OS and RFS, overexpression of KDM5B and KDM5C was significantly associated with shorter OS and shorter RFS, respectively (Figures 3A,B). The above data show that abnormal expression of the KDM5 family can be used as a biomarker to predict the prognosis of patients with pancreatic cancer, which is worthy of further experimental data verification.

To explore the interaction between KDM5 family genes and their co-expressed genes in pancreatic cancer, 800 co-expressed genes of the KDM5 family in pancreatic cancer were obtained from the LinkedOmics database (Supplementary Material S1), and a volcano plot of the KDM5 family and co-expressed genes in pancreatic cancer was drawn (Figure 4A). Then, we analyzed the GO and KEGG pathway terms of these 800 genes using the Metascape database. To study the functional mechanism of the KDM5 family in the occurrence and development of pancreatic cancer, GO analysis was performed. The results showed significant enrichment of the BP terms “peptide metabolic process”, “chromatin organization”, and “cellular response to DNA damage stimulus” (Figure 4B). The enriched CC terms mainly included “cytosolic ribosome”, “transferase complex” and “mitochondrial envelope” (Figure 4C). The enriched MF terms mainly included “structural constituent of ribosome”, “chromatin binding”, and “transcription factor binding” (Figure 4D). The KEGG pathway analysis showed that the target genes were mainly associated with the terms “ribosome”, “chemical carcinogenesis—reactive oxygen species”, and “transcriptional misregulation in cancer” (Figure 4E). To further analyze the mechanisms of the KDM5 family and its coexpressed genes, the top 50 genes with a positive correlation (Figure 5A) and the top 50 genes with a negative correlation (Figure 5B) with KDM5 family genes were visualized in a heatmap. GeneCards was used to annotate the top 10 genes positively related to the expression of KDM5A/B/C and the first 10 genes negatively related to the expression of KDM5D. The results showed that the functions of related genes included transcriptional regulation (DDI2, ASXL2, CCNT1, CNOT1, REST, TRIM44, IRF6, KDM6A, EIF1AX, DDX3X, TSIX, ZFX, TXLNG, KRAS, GPATCH2, and ZRSR2), protein modification (WNK1, UBXN7, RAPGEF6, CHML, WDR26, and SYAP1), DNA damage repair (SMC1A), cell cycle regulation (TP53BP2), and cell migration regulation (F11R). Then, the top 10 genes negatively related to the expression of KDM5A/B/C and the top 10 genes positively related to the expression of KDM5D were functionally annotated. The results showed that the functions of the genes included transcriptional regulation (DDX3Y, ZFY, and RPS4Y1), protein modification (USP9Y and UTY), translation regulation (EIF1AY and RPS4Y1), cell cycle regulation (GADD45GIP1), cell migration regulation (NLGN4Y), cell proliferation regulation (TMSB4Y, TP53I13, EGFL7, and TSPAN33) and apoptosis regulation (DPP7). Notably, our analysis showed that the KRAS gene was positively related to the expression of KDM5B (Figure 5A). KRAS can induce transcriptional silencing of tumor suppressor genes, and its mutations produce modified proteins that drive the occurrence and development of pancreatic cancer (Asimgil et al., 2022). This correlation further suggested that the KDM5 family plays a unique role in malignant tumors driven by KRAS mutations, including pancreatic cancer. At the same time, this important finding inspired us to analyze the correlation between the expression of KDM5 family genes and other pathogenic genes in pancreatic cancer. We used the R (version 3.6.3) to analyze the TCGA PAAD data. The results showed that the expression level of KDM5A/B/C was positively correlated with the expression of KRAS, BRCA1, BRCA2, ATM, SLC16A4, and RABL3. We speculated that although some of the genes known to be pathogenic in pancreatic cancer were not in the top 50 genes positively correlated with the expression of KDM5 family genes, they may still have a co-expression relationship with the KDM5 family. Furthermore, it is suggested that the KDM5 family may play an important role in the pathogenesis of pancreatic cancer. Finally, we used the TIMER database to analyze the relationships among the members of the KDM5 family, and the results showed that the expression of each member was significantly positively correlated: the correlation between KDM5A and KDM5B was the strongest (cor = 0.594) (Figure 5D). According to the above results, we speculated that members of the KDM5 family can cooperate with a variety of pathogenic genes, transcriptional regulatory factors, protein-modifying factors, DNA damage repair factors, cell cycle regulatory factors, and cell migration regulatory factors in pancreatic cancer, including factors related to the terms “chemical carcinogenesis—reactive oxygen species”, and “transcriptional misregulation in cancer”, to participate in the occurrence and development of pancreatic cancer.

To explore the PPIs between KDM5 family genes and related genes in pancreatic cancer, we analyzed the protein expression of members of the KDM5 family and related genes using the STRING database (Supplementary Material S2). The nine genes with the strongest protein-protein interaction (PPI) with KDM5 family members were determined. We analyzed the PPI network related to the KDM5 family expression in pancreatic cancer. The PPI network included 9 gene nodes and 24 edges (Figure 6A). The results showed that the protein expression of KDM5A/B/C was closely related to that of HDAC1. The protein encoded by HDAC1 can deacetylate part of the lysine residue of histone H3 and plays an important role in transcriptional regulation and cell proliferation (Cai et al., 2000). The proteins encoded by the KDM5 family play a role in the regulation of gene expression through the specific demethylation of histone H3 lysine 4 (Zhang et al., 2014). Both proteins act on histone H3 lysines, suggesting that they may have a synergistic effect in the regulation of gene expression. Next, we used the R (version 3.6.3) to analyze the GO and KEGG pathway terms of the 9 genes with the strongest interaction with KDM5 family proteins (Figure 6B). The results of the GO analysis showed significant enrichment of the BP terms “histone lysine demethylation”, “histone demethylation”, and “protein demethylation”. The significantly enriched CC terms mainly included “histone methyltransferase complex”, “methyltransferase complex”, and “Sin3 complex”. The significantly enriched MF terms mainly included “histone demethylase activity”, “demethylase activity”, and “transcription corepressor activity”. The KEGG pathway analysis showed that the target genes were mainly associated with the terms “longevity-regulating pathway”, “Notch signaling pathway”, and “amphetamine addiction”. The above results show that KDM5 family proteins may have strong interactions with HDAC1 and other proteins, modify the lysines of histone H3, and act on many pathways, including the “longevity-regulating pathway” and “Notch signaling pathway”, which play a regulatory role in the proliferation of pancreatic cancer cells.

The immune system is a complex system in which immune cells both act as the first line of defense against a variety of pathogens and provide surveillance by identifying and destroying latent cancer cells. However, in some cases, the immune system can help tumor cells escape immune control (Shalapour and Karin, 2015). Tumor-infiltrating lymphocytes are a unique kind of lymphocytes that infiltrate the tumor microenvironment by detecting cancer antigens and releasing proinflammatory cells. We used the TIMER database to further explore the relationship between the expression of KDM5 family genes and the level of infiltrating lymphocytes. The results showed that upregulation of KDM5A/B/C expression was associated with increased infiltration of B cells, CD8⁺ T cells, macrophages, neutrophils, dendritic cells and other infiltrating lymphocytes (Figures 7A–C) and the upregulation of KDM5D expression was not associated with increased infiltration of B cells, CD8⁺ T cells, macrophages, neutrophils, dendritic cells and other infiltrating lymphocytes (Figure 7D). Next, we analyzed the correlation between KDM5 family expression and immunostimulatory molecules (Figure 7E), immunosuppressive molecules (Figure 7F), MHC molecules (Figure 7G), chemokines (Figure 7H), and chemokine receptors (Figure 7I) in pancreatic cancer using the TISIDB database. The results showed that upregulated expression of KDM5A/B/C was associated with an increase in the expression of immunostimulatory molecules such as NT5E, TNFSF4, and TNFSF15, immunosuppressive molecules such as CD274 and TGFBR1, MHC molecules such as TAP2, chemokines such as CCL24, and chemokine receptors such as CCR8 and CCR9 (Figures 7E–I), which provides important information for predicting potential therapeutic targets. Finally, we used GSCALiteonlinetool to analyze the relationship between the expression of KDM5 family genes and sensitivity to current immune or targeted therapies for pancreatic cancer (Figure 7J). The results showed that the expression levels of KDM5A and KDM5C were negatively correlated with sensitivity to many pancreatic cancer-targeting or chemotherapeutic drugs, including axitinib and gemcitabine. Thus, the KDM5 family might represent a new target for predicting drug sensitivity and for developing multitarget combined therapy for pancreatic cancer.