Hao Li
Lei Miao2†
Jiayuan Wu- 1Department of Orthopedics, Beijing Jishuitan Hospital, Capital Medical University, Beijing, China
- 2School of Computer Science, South China Business College of Guangdong University of Foreign Studies, Guangzhou, China
Objective: To identify cell senescence-related genes and immune cell infiltration characteristics in intervertebral disc degeneration (IVDD) using bioinformatics and investigate their biological functions and signaling pathways.
Methods: GSE56081 and GSE23130 datasets were downloaded from the Gene Expression Omnibus (GEO) database. The former was used for the analysis of differentially expressed genes (DEGs) and the latter was used as the validation set. Intersection analysis of DEGs and cell senescence-related genes was performed to screen for senescence-related differentially expressed genes (SRDEGs). SRDEGs were analyzed using Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG). A protein–protein interaction (PPI) network was also drawn, and the hub SRDEGs were obtained using 11 different algorithms. Immune infiltration analysis was then performed. Receiver operating characteristic (ROC) curve and quantitative polymerase chain reaction (qPCR) were used to evaluate the diagnostic value of the hub SRDEGs.
Results: Four hub SRDEGs, SP1, FOXO1, ESR1, and MAPK1, were identified. SP1 was downregulated in IVDD, while the other three hub genes were upregulated. ROC curve verification showed that the AUC for SP1, FOXO1, ESR1, and MAPK1 were 0.483, 0.683, 0.683, and 0.725, respectively. RT-qPCR confirmed that MAPK1 expression was higher in the degeneration group (t = 3.229, P <0.001). Immune infiltration analysis demonstrated elevated proportions of CD8 T cells and activated NK cells in IVDD samples, with MAPK1 showing positive correlations with CD8 T cells, activated NK cells, and neutrophils but negative correlations with naive CD4 T cells, B memory cells, and resting NK cells.
Conclusion: These findings highlight MAPK1 as a pivotal regulator of cellular senescence and immune cell infiltration in IVDD pathogenesis, offering novel therapeutic targets for intervention.
1 Introduction
Intervertebral disc degeneration (IVDD), which is highly age-related and can give rise to a host of spinal surgical diseases, such as lower back pain, disc herniation, spinal stenosis, and degenerative scoliosis, is the pathological basis of many degenerative spinal diseases. These conditions affect patients’ quality of life and impose a heavy burden on social healthcare (1–4). Moreover, with an increasing aging population, this will inevitably become a problem that needs to be urgently solved in today’s society. IVDD mainly includes the death of nucleus pulposus cells and degradation of the extracellular matrix, which is closely related to the senescence of nucleus pulposus cells (5–7).
Cellular senescence, which is characterized by irreversible cessation of the cell cycle, is a process in which cells gradually lose their ability to proliferate and differentiate and experience a decline in various physiological functions during the course of life (8). Although the cell cycle is arrested and proliferation ceases, senescent cells still possess metabolic activity and an enhanced secretory function. They influence the microenvironment of adjacent cells and tissues through autocrine and paracrine pathways, leading to and accelerating the occurrence of senescence-related diseases such as IVDD (4, 9, 10). Inflammatory responses, oxidative stress, and mitochondrial dysfunction can accelerate cellular senescence, thereby promoting intervertebral disc degeneration. However, the specific gene targets and molecular mechanisms related to cellular senescence in IVDD remain unclear. Immune cell infiltration is the accumulation and aggregation of immune cells in tissues, such as B cells, T cells, and macrophages. A previous study showed that when the physical barrier between the intervertebral disc and systemic circulation is damaged, various immune cells are involved in the inflammatory response or IVDD progression (11). The present study, through bioinformatics analysis methods, comprehensively explores the key roles of cellular senescence-related genes and immune cell infiltration in the onset and progression of IVDD, aiming to provide new diagnostic and therapeutic targets for IVDD.
2 Materials and methods
2.1 Data collection
The datasets GSE56081 and GSE23130 used in this study were obtained from the Gene Expression Omnibus (GEO) database (https://www.ncbi.nlm.nih.gov/geo/). Dataset GSE56081 (12), containing five degenerated nucleus pulposus samples and five non-degenerated ones, was used to identify differentially expressed genes (DEGs). Dataset GSE23130 (13), with eight degenerated and 15 non-degenerated nucleus pulposus samples, was the validation analysis dataset.
Cellular senescence-related genes were obtained from the CellAge online database (https://www.genomics.senescence.info/cells), total 866 in number (14).
2.2 Determination of DEGs and senescence-related DEGs
The DEGs in IVDD were acquired using the R package “limma,” and the criteria for identifying the DEGs were set at |log2 FC| >1 and adjusted p-value <0.05. A Venn diagram was used to intersect the DEGs from GSE56081 with the cellular senescence-related genes from the CellAge database to obtain senescence-related DEGs (SRDEGs). Volcano plots of DEGs and hierarchical clustering heatmaps of SRDEGs were generated using the R package “ggplot2.”
2.3 Functional enrichment analysis
The “clusterProfiler” package in R was used for GO and KEGG enrichment analyses of SRDEGs. GO enrichment analysis classified gene functions into three categories: Biological Process (BP), Cellular Component (CC), and Molecular Function (MF). KEGG enrichment analysis covers a wide range of biochemical processes, divided into seven main categories: metabolism, genetic information processing, environmental information processing, cellular processes, organismal systems, human diseases, and drug development. These processes are represented as pathway maps of molecular interactions and reactions. The threshold for significance was set at p <0.05. Enrichment analysis plots were visualized using the R package “ggplot2.”
2.4 Protein–protein interaction network analysis
PPI network analysis of SRDEGs was conducted using STRING (https://cn.string-db.org/), a common tool for assessing protein interactions. The results were then imported into Cytoscape v3.10.3 (https://cytoscape.org/) to construct a PPI network. In this network, the nodes represent SRDEGs and lines (connections between nodes) represent the interactions between different SRDEGs.
2.5 Screening and correlation analysis of the hub SRDEGs
Based on the constructed PPI network, 11 algorithms (CycloHubba MCC, DMNC, MNC, Degree, EPC, BottleNeck, EcCentricity, Closeness, Radiality, Betweenness, Stress) were used to calculate the top 30 SRDEGs. Different algorithms represent different topological features based on network analysis, and the selection of multiple algorithms can provide a more robust and comprehensive evaluation of potential hub genes from multiple network perspectives (15). The intersection of these results was considered the hub SRDEGs. Correlations between hub SRDEGs were assessed using Pearson correlation analysis. In Pearson correlation analysis, the r value refers to the correlation coefficient and is used to assess the effect size. A scatter plot of the correlation was drawn using the R package “ggplot2.”
2.6 Functional enrichment analysis of the hub SRDEGs
GO and KEGG enrichment analyses of the hub SRDEGs were performed using the R software package “clusterProfiler,” and the threshold p-value was set at <0.05. The R package “ggplot 2” was used for visualization of enrichment analysis plots.
2.7 Immune infiltration analysis by CIBERSORT algorithm
Immune cell infiltration in IVDD and control samples was evaluated using the CIBERSORT bioinformatic algorithm. This algorithm, along with the R package “CIBERSORT” and the leukocyte gene signature matrix LM22, was employed to simulate and calculate the transcription feature matrix of 22 immune cells with 1,000 permutations and a significance threshold of p <0.05. Differences in the percentage of immune cells between IVDD and control samples were assessed using the Wilcoxon test, with P <0.05 to be statistically significant. The correlations between immune cells and potential biomarkers, as well as the correlations between different immune cells, were calculated using the R package “corrplot” and Pearson correlation analysis.
2.8 ROC analysis of the hub SRDEGs
The R package “pROC” was used for Receiver operating characteristics (ROC) analysis. Gene expression of the hub SRDEGs from dataset GSE23130 was analyzed using the “roc function”, and the “ci function” was used to evaluate the area under the curve (AUC) and confidence interval. In the ROC curve, sensitivity is on the Y-axis and false-positive rate (1−specificity) on the X-axis. An ROC curve with an AUC of at least 0.70 indicates sufficient predictive value.
2.9 RNA extraction and RT-qPCR
Eight nucleus pulposus tissue samples were collected at our hospital for PCR verification. According to the Pfirrmann grading system for disc degeneration, eight patients were divided into two groups: IVDD group (grades III–IV) and control group (grades I–II) (16). Written informed consent was obtained from all the patients. Total RNA was extracted from the nucleus pulposus tissue using the TRIzol reagent (Thermo Fisher Scientific). Prime-ScriptTM Master Mix (code no. RR036A) was used to reverse transcribe the total RNA into cDNA. Quantitative PCR was performed using SYBR® Premix Ex Taq™ (Tli RNaseH Plus) on a LightCycler 480 II real-time fluorescence quantitative PCR instrument (ROCHE) with GAPDH as the internal control. Relative mRNA content was calculated using the 2−ΔΔCt method. The primers used were as follows:
GAPDH forward primer (F): 5’-AGACACCATGGGGAAGGTGAA-3’, reverse primer (R): 5’-ATTGCTGATGATCTTGAGGCTG-3’;
MAPK1 forward primer (F): 5’-CTGCTGCTCAACACCACCT-3’, reverse primer (R): 5’-GCCACATATTCTGTCAGGAACC-3’;
FoxO1 forward primer (F): 5’-GAGGAGCCTCGATGTGGATG-3’, reverse primer (R): 5’-CCGAGATTTGGGGGAACGAA-3’;
ESR1 forward primer (F): 5’-TGGCTAGCGCTAGCTAGGGGGCTACA-3’, reverse primer (R): 5’-TGGGGCTAGGCTGCTCTAGCC-3’.
2.10 Statistical analysis
In the present study, statistical analyses were conducted using R 4.4.2 and GraphPad Prism 10.1.2. The unpaired Student’s t-test was applied to assess the difference between the IVDD and control groups, and a p-value of less than 0.05 was considered to be statistically significant.
3 Results
3.1 Determination of the senescence-related differentially expressed genes
Differentially expressed genes (DEGs) between the IVDD and control groups were identified using the GSE56081 dataset. A volcano plot (Figure 1A) showed 2,801 DEGs, with 2002 upregulated and 799 downregulated. We then intersected 866 senescence-related genes from the CellAge database with 2,801 DEGs and identified 106 senescence-related differentially expressed genes (SRDEGs) (Figure 1B), including 74 upregulated and 32 downregulated genes (Table 1). A hierarchical clustering heatmap demonstrated significant differences in SRDEG expression between the IVDD and control groups (Figure 2).
Figure 1. Identification of SRDEGs between IVDD and control groups. (A) Volcano plot of DEGs based on dataset GSE56081, where Y-axis is −log 10 (adjusted p value) and X-axis is log2 FC. Green dots represent upregulated genes and purple dots represent down-regulated genes. Gray dots represent genes that were not significantly differentially expressed. (B) Venn diagram based on DEGs and senescence-related genes based on CellAge database.
Table 1. List of the senescence-related differentially expressed genes.
Figure 2. The heatmap of the 106 senescence-related differentially expressed genes (SRDEGs). Color scale represents the relative gene expression of each sample, the yellow represents high expression abundance and blue represents low expression abundance; Group 1 is the five samples from the control group and Group 2 is the five samples from the IVDD group.
3.2 GO and KEGG enrichment analysis
GO enrichment analysis indicated the following. In biological processes, SRDEGs were significantly enriched in terms of regulation of protein stability, regulation of miRNA transcription, and regulation of miRNA metabolism (Figure 3A). In the cellular components, SRDEGs were highly enriched in the vesicle lumen, transcription repressor complex, and secretory granule lumen (Figure 3B). In terms of the molecular function, SRDEGs were notably enriched in ubiquitin and ubiquitin-like ligase binding, transcriptional co-repressor activity, and transcriptional co-regulator binding (Figure 3C). KEGG enrichment analysis showed that SRDEGs were mainly associated with the FoxO, MAPK, and AGE-RAGE signaling pathways and cellular senescence (Figure 3D).
Figure 3. Functional enrichment analysis of SRDEGs. (A) The top 10 significantly enriched GO terms in the category of biological process (BP) for the SRDEGs; (B) The top 10 significantly enriched GO terms in the category of cellular component (CC); (C) The top 10 significantly enriched GO terms in the category of molecular function (MF); (D) Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analysis for the SRDEGs.
3.3 PPI network analysis and correlation analysis of the SRDEGs
PPI network analysis based on the STRING database was performed and visualized using the Cytoscape software (Figure 4). The top 30 SRDEGs were then calculated using 11 algorithms (Figure 5A) in the CycloHubba plug-in, and the intersection was considered as the hub SRDEGs. Four hub SRDEGs were identified: SP1, FoxO1, ESR1, and MAPK1 (Figure 5B). SP1 was downregulated in IVDD, while the other three hub genes were upregulated (Table 2). In addition, Pearson correlation analysis was performed on the four hub SRDEGs (Figure 5C). It can be seen that the most negatively correlated is SP1-ESR1, and the most positively correlated is MAPK1-ESR1. Their correlation analysis is shown in Figure 6A, which shows a negative correlation between SP1 and ESR1 expression with an r value of −0.9 and p <0.01. Figure 6B shows that FoxO1 expression is positively correlated with MAPK1 expression with an r value of 0.77 and p <0.01.
Figure 4. PPI network analysis of SRDEGs.
Figure 5. (A) Upset graph of eleven algorithms (CycloHubba MCC, DMNC, MNC, Degree, EPC, BottleNeck, EcCentricity, Closeness, Radiality, Betweenness, Stress); (B) Four hub SRDEGs: SP1, FoxO1, ESR1 and MAPK1; (C) A Pearson correlation analysis of four hub SRDEGs.
Table 2. List of four hub SRDEGs.
Figure 6. (A) Pearson correlation analysis showed that ESR1 expression was negatively correlated with SP1 expression. (B) Pearson correlation analysis showed that ESR1 expression was positively correlated with MAPK1 expression.
3.4 Function analysis of the hub SRDEGs
GO enrichment analysis revealed that the hub SRDEGs were mainly enriched in cellular response to starvation, cellular response to insulin stimulus, cellular response to estrogen stimulus, and cellular response to nutrient levels in biological processes (Figure 7A), chromatin and pseudopodia in cellular components (Figure 7B), and DNA-binding transcriptional activator activity, transcriptional coregulator binding, and MAP kinase activity in molecular functions (Figure 7C). KEGG enrichment analysis indicated that the key SRDEGs were mainly associated with the AGE-RAGE, TNF-β, and prolactin signaling pathways (Figure 7D).
Figure 7. Functional enrichment analysis of hub SRDEGs. (A) The top 10 significantly enriched GO terms in the category of biological process for these hub SRDEGs; (B) The top 10 significantly enriched GO terms in the category of cellular component; (C) The top 10 significantly enriched GO terms in the category of molecular function; (D) KEGG pathway analysis for the hub SRDEGs.
3.5 Immune infiltration analysis
Figure 8 shows the proportion of 22 immune cells in the samples using a bar chart. As indicated by the CIBERSORT algorithm, neutrophils and M0 macrophages showed the most significant positive correlation. In contrast, CD8 T cells and naive CD4 T cells displayed the strongest negative correlation (Figure 9). A box plot of immune cells (Figure 10) revealed that the control group had relatively higher levels of naive CD4 T cells and resting NK cells. Conversely, the IVDD groups exhibited higher quantities of CD8 T cells and activated NK cells, which may be related to the IVDD process.
Figure 8. The proportion of 22 immune cells in samples in GSE56081. The x-axis represents the proportion of 22 immune cells, and the y-axis represents the 10 samples from GSE56081.
Figure 9. The correlations of 22 immune cells. The redder the color, the stronger the negative correlation; the bluer the color, the stronger the positive correlation.
Figure 10. Boxplot of the difference in the proportion of immune cells between the IVDD group and the control group. The x-axis represents the 22 types of immune cells, and the y-axis represents the proportion of cellular composition in the samples. The levels of naive CD4 T cells and resting NK cells were relatively high in the control group, and the levels of CD8 T cells and activated NK cells were higher in the IVDD group. *p <0.05, **p <0.01.
In addition, the correlations between the hub SRDEGs and immune cells were as follows: SP1 was positively correlated with naive CD4 T cells and B memory cells and negatively correlated with CD8 T cells, activated NK cells, and neutrophils. MAPK1 is negatively correlated with naive CD4 T cells, B memory cells, and resting NK cells, whereas it is positively correlated with CD8 T cells, activated NK cells, and neutrophils. FoxO1 is negatively associated with naive CD4 T cells, resting dendritic cells, and resting NK cells and positively associated with CD8 T cells, activated NK cells, and neutrophils. ESR1 was negatively correlated with naive CD4 T cells and B memory cells and positively correlated with CD8 T cells and activated NK cells (Figure 11).
Figure 11. The correlations between the hub SRDEGs and immune cells. The x-axis represents the 22 types of immune cells, and the y-axis represents the four hub SRDEGs. Red indicates positive correlation (darker red signifies stronger positive correlation), while blue indicates negative correlation (darker blue signifies stronger negative correlation). *p <0.05, **p <0.01.
3.6 The validation of hub SRDEGs in GSE23130 dataset
Based on the dataset GSE23130, the hub SRDEGs were verified using the ROC curve, as shown in Figure 12. The AUC for SP1 is 0.483, FoxO1 is 0.683, ESR1 is 0.683, and MAPK1 is 0.725, which is above the threshold of 0.7.
Figure 12. ROC curve of hub SRDEGs based on GSE23130 dataset, including SP1, FOXO1, ESR1, and MAPK1. The X-axis represents the (1 − Specificity), and the Y-axis represents the Sensitivity.
3.7 Expression validation of hub SRDEGs
RT-qPCR confirmed the difference in MAPK1 expression in the nucleus pulposus between the degeneration and control groups, showing higher MAPK1 expression in the IVDD group (t = 3.229, P <0.001). FoxO1 and ESR1 expression was also confirmed to be higher in the IVDD group (t = 5.675, P = 0.0013; t = 5.341, P = 0.0018, respectively) (Figure 13).
Figure 13. The expression differences of MAPK1, FoxO1 and ESR1 genes between IVDD group and control group. **p <0.01, ***p <0.001.
4 Discussion
Intervertebral disc degeneration (IVDD), marked by nucleus pulposus cell reduction and extracellular matrix degradation, is a major cause of lower back pain and many degenerative diseases, posing a heavy global healthcare burden (17, 18). However, the pathogenesis of this condition remains controversial. Oxidative stress, inflammatory conditions, abnormal mechanical stress, and nutritional metabolic imbalance are contributing factors. Cell senescence plays a key role in these mechanisms and is increasingly being found to be closely related to the occurrence and development of degenerative disc diseases (9, 19, 20).
The present study analyzed GSE56081 datasets using bioinformatics and found 2,801 DEGs, with 2,002 upregulated and 799 downregulated genes. Intersection with senescence-related genes yielded 106 SRDEGs (74 upregulated and 32 downregulated). GO function analysis showed that SRDEGs were mainly involved in regulating protein stability, miRNA transcription, and miRNA metabolism in biological process. Intracellular protein degradation is primarily maintained through the ubiquitin-proteasome system (UPS) and autophagy. In senescent cells, the decline in UPS function and autophagy efficiency leads to the accumulation of misfolded or oxidatively damaged proteins, activating the unfolded protein response and oxidative stress, and ultimately triggering cell cycle arrest (21, 22). Additionally, the significant enrichment of ubiquitin and ubiquitin-like ligase binding in molecular functions further supports the above viewpoints. KEGG pathway analysis showed that SRDEGs were significantly enriched in the FoxO, MAPK, and AGE-RAGE signaling pathways, which are crucial for cellular senescence.
The FoxO transcription factor family consists of four members, among which FoxO1/3/4 are all involved in senescence-related signaling pathways that influence cellular development, metabolism, and lifespan.
The FoxO signaling pathway is primarily regulated by protein kinase B (PKB/Akt) and mitogen-activated protein kinases (MAPK). Nogueira et al. (23) demonstrated that in mouse embryonic fibroblasts, Akt acts as a senescence-inducing agent by inactivating FoxO1 and FoxO3, including in contexts of replicative senescence, oxidative stress-induced senescence, and oncogenic Ras-induced senescence. The MAPK signaling pathway comprises four major branches, ERK, JNK, p38/MAPK, and ERK5 (24), all of which play critical roles in cellular senescence. The MAPK/ERK pathway (25), a classical pro-survival signaling cascade, is closely linked to the differentiation of mesenchymal stem cells (MSCs). Inhibition of ERK blocks spontaneous MSC differentiation, whereas ERK1/2 suppresses apoptosis. Furthermore, the MAPK/ERK signaling pathway regulates the expression of downstream target genes, such as p16 and p53, with p16 being a key senescence-associated gene (26). Xie et al. (27) demonstrated that ERK, JNK, and p38/MAPK collectively exert significant effects during induced senescence in testicular stem Leydig cells. The AGE-RAGE signaling pathway refers to a cascade triggered by the binding of advanced glycation end products (AGEs), glycated proteins resulting from hyperglycemia, to their receptor RAGE, ultimately activating downstream cellular responses (28). Ferrucci et al. (29) identified the most prominent aging-associated metabolic pathways, including IGF signaling, MAPK signaling, HIF-1 signaling, cytokine signaling, FoxO metabolic pathways, folate metabolism, and AGE-RAGE-related metabolic pathways.
To further identify the hub genes associated with cellular senescence in intervertebral disc degeneration, PPI network analysis combined with 11 computational algorithms identified four hub genes: SP1, FoxO1, ESR1, and MAPK1. Among these, SP1 was downregulated in IVDD, whereas the other three genes (FoxO1, ESR1, and MAPK1) exhibited upregulated expression. External dataset validation using ROC curves revealed the following AUC values: SP1, 0.483; FoxO1, 0.683; ESR1, 0.683; and MAPK1, 0.725. Interestingly, while SP1 was identified as a hub gene based on its high centrality in the SRDEG PPI network, suggesting a key structural role in potential senescence-related interactions, its diagnostic value, assessed by ROC analysis on the validation dataset, was low (AUC = 0.483). This discrepancy highlights the fact that topological importance within a network (hub status) does not always translate directly into strong discriminatory power based on expression levels alone (diagnostic value). The functional impact of SP1 may depend more on changes in activity or post-translational modifications rather than significant alterations in mRNA abundance. RT-qPCR confirmed that MAPK1 expression was higher in the IVDD group. These results showed that MAPK1 may be a key driver of cellular senescence, contributing to IVDD progression. The intervertebral disc is an immune-privileged tissue in the body, covered by the annulus fibrosus and upper and lower cartilaginous endplates (3, 30). However, when nucleus pulposus cells (NPCs) degenerate, they disrupt the NPC-blood barrier, directly exposing NPCs to the host’s immune system. This ultimately triggers immune responses and massive infiltration of immune cells, including macrophages, neutrophils, mast cells, CD4 T cells, CD8 T cells, Treg cells, and B cells (7). In the present study, immune cell infiltration analysis showed that the expression levels of naive CD4 T cells and resting NK cells were higher in the control group, while CD8 T cells and activated NK cells were higher in the IVDD groups, which play pivotal roles in immune-related inflammatory responses associated with IVDD. The correlations between MAPK1 and immune cells showed that MAPK1 was negatively correlated with naive CD4 T cells, B memory cells, and resting NK cells, whereas it was positively correlated with CD8 T cells, activated NK cells, and neutrophils.
The MAPK1 gene, also termed ERK or ERK2, is located on human chromosome 22q11.22, and encodes the ERK2 protein, a core component of the MAPK/ERK signaling pathway. As mentioned previously, the MAPK/ERK pathway (25) is a classical pro-survival signaling cascade that is closely related to the differentiation of MSCs. ERK2 plays an essential role in fundamental biological processes, including cell proliferation, differentiation, survival, senescence, and stress responses. Nojima et al. (31) found that in mouse embryonic fibroblasts, downregulation of IGFBP5 promoted replicative senescence by activating ERK2. Furthermore, in IGFBP5-knockdown cells, silencing of ERK2, but not ERK1, prevented an increase in SA-β-gal-positive cells. Shin et al. (32) also found that the depletion of ERK2, but not ERK1, abrogates oncogenic Ras-induced senescence. Furthermore, ERK2 regulates the expression of the downstream target gene p16, a key gene involved in cell senescence (26, 31). MAPK1/ERK2 exhibits dual regulatory effects in cellular senescence: in normal cells, ERK2 delays senescence by promoting proliferation and inhibiting apoptosis, whereas in senescent cells, its sustained activation may negatively regulate the MAPK pathway through feedback mechanisms and crosstalk with other pathways, such as PI3K/AKT and p53 signaling, ultimately inducing cell cycle arrest (33–35). Furthermore, ERK2 modulates cellular senescence by regulating telomerase activity and telomere length (36). The MAPK cascade plays diverse roles in immune responses, including initiation of innate immunity, activation of adaptive immunity, and regulation of immunogenic cell death (37). ERK2 is a critical upstream regulator of neuroinflammatory pathologies. In demyelinating mouse models, ERK2 activation exacerbates inflammatory demyelination by inducing proinflammatory mediators and gliosis, whereas ERK2-deficient mice exhibit reduced levels of these inflammatory mediators (38). The role of ERK2 in immune responses during intervertebral disc degeneration remains poorly understood, and future studies should focus on elucidating its regulatory mechanisms in immune cell infiltration, pro-inflammatory cytokine production, and extracellular matrix degradation within the disc microenvironment.
The present study had some limitations. First, this study predominantly relied on bioinformatic analyses, with only preliminary RT-qPCR validation performed for hub genes. Protein-level validation, such as western blotting, was not performed to assess total and phosphorylated MAPK1/ERK2. Additionally, the lack of functional assays, such as siRNA knockdown of MAPK1 followed by senescence marker quantification, SASP-related gene expression profiling, or co-culture assays with immune cells, makes it difficult to draw definitive conclusions regarding the causal role of MAPK1 in cellular senescence and immune infiltration. Future studies that incorporate these missing experimental validations are essential to confirm the proposed mechanism.
5 Conclusion
In conclusion, the present study employed bioinformatics techniques to investigate the mechanisms of cellular senescence and immune infiltration in intervertebral disc degeneration. By integrating multiple algorithms to screen for senescence-related genes in IVDD, our findings highlight MAPK1 as a critical factor that plays a pivotal role in regulating both senescence processes in nucleus pulposus cells and immune infiltration mechanisms.
Data availability statement
Publicly available datasets were analyzed in this study. This data can be found here: GSE56081 and GSE23130.
Ethics statement
The studies involving humans were approved by the Ethics Committee of the Affiliated Beijing Jishuitan Hospital. The studies were conducted in accordance with the local legislation and institutional requirements. The participants provided their written informed consent to participate in this study.
Author contributions
HL: Conceptualization, Formal Analysis, Investigation, Methodology, Software, Writing – original draft. LM: Conceptualization, Data curation, Methodology, Software, Visualization, Writing – original draft. JW: Conceptualization, Methodology, Supervision, Writing – original draft. YW: Methodology, Writing – original draft. JX: Methodology, Writing – original draft. DH: Conceptualization, Methodology, Project administration, Supervision, Writing – review & editing.
Funding
The author(s) declare that financial support was received for the research and/or publication of this article. This work was supported in part by the Beijing Hospitals Authority Ascent Plan (DFL20240401 to DH).
Conflict of interest
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
Generative AI statement
The author(s) declare that no Generative AI was used in the creation of this manuscript.
Publisher’s note
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.
Supplementary material
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fimmu.2025.1589849/full#supplementary-material
References
1. Maher C, Underwood M, and Buchbinder R. Non-specific low back pain. Lancet (London England). (2017) 389:736–47. doi: 10.1016/S0140-6736(16)30970-9
2. Luoma K, Riihimäki H, Luukkonen R, Raininko R, Viikari-Juntura E, and Lamminen A. Low back pain in relation to lumbar disc degeneration. Spine. (2000) 25:487–92. doi: 10.1097/00007632-200002150-00016
3. Risbud MV and Shapiro IM. Role of cytokines in intervertebral disc degeneration: pain and disc content. Nat Rev Rheumatol. (2014) 10:44–56. doi: 10.1038/nrrheum.2013.160
4. Zhang X, Zhang Z, Zou X, Wang Y, Qi J, Han S, et al. Unraveling the mechanisms of intervertebral disc degeneration: an exploration of the p38 MAPK signaling pathway. Front Cell Dev Biol. (2023) 11:1324561. doi: 10.3389/fcell.2023.1324561
5. Feng C, Yang M, Zhang Y, Lan M, Huang B, Liu H, et al. Cyclic mechanical tension reinforces DNA damage and activates the p53-p21-Rb pathway to induce premature senescence of nucleus pulposus cells. Int J Mol medicine. (2018) 41:3316–26. doi: 10.3892/ijmm.2018.3522
6. Silwal P, Nguyen-Thai AM, Mohammad HA, Wang Y, Robbins PD, Lee JY, et al. Cellular senescence in intervertebral disc aging and degeneration: molecular mechanisms and potential therapeutic opportunities. Biomolecules. (2023) 13:686. doi: 10.3390/biom13040686
7. Wang F, Cai F, Shi R, Wang XH, and Wu XT. Aging and age related stresses: a senescence mechanism of intervertebral disc degeneration. Osteoarthritis cartilage. (2016) 24:398–408. doi: 10.1016/j.joca.2015.09.019
8. He S and Sharpless NE. Senescence in health and disease. Cell. (2017) 169:1000–11. doi: 10.1016/j.cell.2017.05.015
9. Childs BG, Durik M, Baker DJ, and van Deursen JM. Cellular senescence in aging and age-related disease: from mechanisms to therapy. Nat medicine. (2015) 21:1424–35. doi: 10.1038/nm.4000
10. Gruber HE, Ingram JA, Norton HJ, and Hanley EN Jr. Senescence in cells of the aging and degenerating intervertebral disc: immunolocalization of senescence-associated beta-galactosidase in human and sand rat discs. Spine. (2007) 32:321–7. doi: 10.1097/01.brs.0000253960.57051.de
11. Ye F, Lyu FJ, Wang H, and Zheng Z. The involvement of immune system in intervertebral disc herniation and degeneration. JOR spine. (2022) 5:e1196. doi: 10.1002/jsp2.1196
12. Wan ZY, Song F, Sun Z, Chen YF, Zhang WL, Samartzis D, et al. Aberrantly expressed long noncoding RNAs in human intervertebral disc degeneration: a microarray related study. Arthritis Res Ther. (2014) 16:465. doi: 10.1186/s13075-014-0465-5
13. Gruber HE, Hoelscher GL, Ingram JA, and Hanley EN Jr. Genome-wide analysis of pain-, nerve- and neurotrophin -related gene expression in the degenerating human annulus. Mol pain. (2012) 8:63. doi: 10.1186/1744-8069-8-63
14. de Magalhães JP, Abidi Z, Dos Santos GA, Avelar RA, Barardo D, Chatsirisupachai K, et al. Human Ageing Genomic Resources: updates on key databases in ageing research. Nucleic Acids Res. (2024) 52:D900–d8. doi: 10.1093/nar/gkad927
15. Zhang X, Yu H, Bai R, and Ma C. Identification and characterization of biomarkers and their role in opioid addiction by integrated bioinformatics analysis. Front Neurosci. (2020) 14:608349. doi: 10.3389/fnins.2020.608349
16. Wang YXJ. Several concerns on grading lumbar disc degeneration on MR image with Pfirrmann criteria. J orthopaedic translation. (2022) 32:101–2. doi: 10.1016/j.jot.2021.12.003
17. Liang H, Luo R, Li G, Zhang W, Song Y, and Yang C. The proteolysis of ECM in intervertebral disc degeneration. Int J Mol Sci. (2022) 23:1715. doi: 10.3390/ijms23031715
18. Zhou Z, Hui ES, Kranz GS, Chang JR, de Luca K, Pinto SM, et al. Potential mechanisms underlying the accelerated cognitive decline in people with chronic low back pain: A scoping review. Ageing Res Rev. (2022) 82:101767. doi: 10.1016/j.arr.2022.101767
19. Dimozi A, Mavrogonatou E, Sklirou A, and Kletsas D. Oxidative stress inhibits the proliferation, induces premature senescence and promotes a catabolic phenotype in human nucleus pulposus intervertebral disc cells. Eur Cells materials. (2015) 30:89–102. doi: 10.22203/eCM.v030a07
20. Li P, Hou G, Zhang R, Gan Y, Xu Y, Song L, et al. High-magnitude compression accelerates the premature senescence of nucleus pulposus cells via the p38 MAPK-ROS pathway. Arthritis Res Ther. (2017) 19:209. doi: 10.1186/s13075-017-1384-z
21. Tsakiri EN and Trougakos IP. The amazing ubiquitin-proteasome system: structural components and implication in aging. Int Rev Cell Mol Biol. (2015) 314:171–237. doi: 10.1016/bs.ircmb.2014.09.002
22. Palmer JE, Wilson N, Son SM, Obrocki P, Wrobel L, Rob M, et al. Autophagy, aging, and age-related neurodegeneration. Neuron. (2025) 113:29–48. doi: 10.1016/j.neuron.2024.09.015
23. Nogueira V, Park Y, Chen CC, Xu PZ, Chen ML, Tonic I, et al. Akt determines replicative senescence and oxidative or oncogenic premature senescence and sensitizes cells to oxidative apoptosis. Cancer Cell. (2008) 14:458–70. doi: 10.1016/j.ccr.2008.11.003
24. Dent P, Yacoub A, Fisher PB, Hagan MP, and Grant S. MAPK pathways in radiation responses. Oncogene. (2003) 22:5885–96. doi: 10.1038/sj.onc.1206701
25. Sun Y, Liu WZ, Liu T, Feng X, Yang N, and Zhou HF. Signaling pathway of MAPK/ERK in cell proliferation, differentiation, migration, senescence and apoptosis. J receptor Signal transduction Res. (2015) 35:600–4. doi: 10.3109/10799893.2015.1030412
26. Kim J and Wong PK. Oxidative stress is linked to ERK1/2-p16 signaling-mediated growth defect in ATM-deficient astrocytes. J Biol Chem. (2009) 284:14396–404. doi: 10.1074/jbc.M808116200
27. Qigen X, Haiming C, Kai X, Yong G, and Chunhua D. Prenatal DEHP exposure induces premature testicular aging by promoting leydig cell senescence through the MAPK signaling pathways. Advanced Biol. (2023) 7:e2300130. doi: 10.1002/adbi.202300130
28. Vitorakis N and Piperi C. Pivotal role of AGE-RAGE axis in brain aging with current interventions. Ageing Res Rev. (2024) 100:102429. doi: 10.1016/j.arr.2024.102429
29. Moaddel R, Ubaida-Mohien C, Tanaka T, Lyashkov A, Basisty N, Schilling B, et al. Proteomics in aging research: A roadmap to clinical, translational research. Aging Cell. (2021) 20:e13325. doi: 10.1111/acel.13325
30. Wang L, He T, Liu J, Tai J, Wang B, Zhang L, et al. Revealing the immune infiltration landscape and identifying diagnostic biomarkers for lumbar disc herniation. Front Immunol. (2021) 12:666355. doi: 10.3389/fimmu.2021.666355
31. Nojima I, Hosoda R, Toda Y, Saito Y, Ueda N, Horimoto K, et al. Downregulation of IGFBP5 contributes to replicative senescence via ERK2 activation in mouse embryonic fibroblasts. Aging. (2022) 14:2966–88. doi: 10.18632/aging.203999
32. Shin J, Yang J, Lee JC, and Baek KH. Depletion of ERK2 but not ERK1 abrogates oncogenic Ras-induced senescence. Cell signalling. (2013) 25:2540–7. doi: 10.1016/j.cellsig.2013.08.014
33. Gong T, Si K, Liu H, and Zhang X. Research advances in the role of MAPK cascade in regulation of cell growth, immunity, inflammation, and cancer. Zhong nan da xue xue bao Yi xue ban. (2022) 47:1721–8. doi: 10.11817/j.issn.1672-7347
34. Zou J, Lei T, Guo P, Yu J, Xu Q, Luo Y, et al. Mechanisms shaping the role of ERK1/2 in cellular senescence (Review). Mol Med reports. (2019) 19:759–70. doi: 10.3892/mmr.2018.9712
35. Hernandez-Segura A, Nehme J, and Demaria M. Hallmarks of cellular senescence. Trends Cell Biol. (2018) 28:436–53. doi: 10.1016/j.tcb.2018.02.001
36. Deschênes-Simard X, Gaumont-Leclerc MF, Bourdeau V, Lessard F, Moiseeva O, Forest V, et al. Tumor suppressor activity of the ERK/MAPK pathway by promoting selective protein degradation. Genes Dev. (2013) 27:900–15. doi: 10.1101/gad.203984.112
37. Dong C, Davis RJ, and Flavell RA. MAP kinases in the immune response. Annu Rev Immunol. (2002) 20:55–72. doi: 10.1146/annurev.immunol.20.091301.131133
Keywords: GEO database, intervertebral disc degeneration, cellular senescence, immune infiltration, MAPK1
Citation: Li H, Miao L, Wu J, Wang Y, Xia J and He D (2025) Identification and validation of cellular senescence-related genes and immune cell infiltration characteristics in intervertebral disc degeneration. Front. Immunol. 16:1589849. doi: 10.3389/fimmu.2025.1589849
Received: 08 March 2025; Accepted: 07 May 2025;
Published: 29 May 2025.
Edited by:
Zhe Chen, Shanghai Jiao Tong University, ChinaReviewed by:
Guoqing Tang, Kunshan Traditional Chinese Medicine Hospital, ChinaLiang Kang, First Affiliated Hospital of Anhui Medical University, China
Ding-Qiang Chen, First Affiliated Hospital of Xiamen University, China
Xiu Yang, Fuzong Clinical Medical College of Fujian Medical University, China
Copyright © 2025 Li, Miao, Wu, Wang, Xia and He. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Da He, aGVkYUBqc3QtaG9zcC5jb20uY24=
†These authors have contributed equally to this work