Skip to main content

Identification and validation of biomarkers related to Th1 cell infiltration in neuropathic pain

Abstract

Neuropathic pain (NP) is a widespread chronic pain with a prevalence of 6.9–10% in the general population, severely affecting patients’ physical and mental health. Accumulating evidence indicated that the immune environment is an essential factor causing NP. However, the mechanism is unclear. This study attempted to analyze NP-related immune infiltration patterns. We downloaded the expression profiles from the Gene Expression Omnibus (GEO) database. The novel method of single-sample gene set enrichment analysis (ssGSEA) algorithm and weighted gene co-expression network analysis (WGCNA) was applied to identify immune-related genes and verified in vitro and in vivo experiments. The spared nerve injury (SNI) group was closely related to type1 T helper cells (Th1 cells), and two key genes (Abca1 and Fyb) positively correlated with Th1 cell infiltration. At the single-cell level, Abca1 and Fyb were significantly expressed in macrophages. In addition, we verified that Abca1 could affect the function of macrophages. Finally, we hypothesized that Abca1 is involved in the infiltration of Th1 cells into dorsal root ganglion (DRG) tissues and induces NP via immunoinflammatory response. Hence, the present study aimed to elucidate the correlation between NP and neuroinflammation and identify a new therapeutic target for treating NP.

Introduction

Neuropathic pain (NP) was redefined as pain caused by a lesion or disease of the somatosensory nervous system by The International Association for the Study of Pain (IASP) in 2011 [1]. The main signs and symptoms include spontaneous pain, allodynia, and hyperalgesia [2]. The prevalence of NP is about 6.9–10% worldwide [3], but the diverse etiology, complex pathogenesis, and variable clinical presentation have no effect on the clinical treatment for NP. Thus, understanding the pathogenesis of NP is an urgent requirement to explore the therapeutic targets and develop the targeted drugs [4].

Previous studies on NP have focused on the neuronal changes and the role of peripheral and central sensitization after nerve injury. Accumulating evidence suggested that immune cells and related cytokines are closely associated to NP [5,6,7]. After peripheral nerve injury, the activation and migration of immune cells and the release of immunologically active substances stimulate the immune response or inflammatory damage of the nervous system, affect the excitability of neurons and the synaptic connection between neurons, and promote the occurrence and development of NP [5, 8].

Some studies illustrated that the interaction between neurons and non-neurons in the spinal dorsal horn is a vital factor that induces and maintains NP after peripheral nerve injury [8, 9]. Furthermore, Costigan et al. [10] demonstrated that T cells infiltration and activation in the spinal dorsal horn after peripheral nerve injury contributed to neuropathic hypersensitivity, especially high-level Th1 cells in adults. However, the contribution of non-neuronal cells in the DRG has been studied only slightly. Recently, Liu et al. [11]. demonstrated that macrophages in mice DRG induce and maintain NP. The study also confirmed the interaction between sensory neurons and macrophages, revealing the potential peripheral DRG targets for NP therapy. Therefore, exploring the infiltration pattern of immune cells in the DRG provides a new vision for studying the NP caused by peripheral nerve injury.

Various tools have been developed to search for immune-related biomarkers since the rapid development of bioinformatics [12] but have been rarely used in NP research. In this study, we first applied the single sample gene set enrichment analysis (ssGSEA) algorithm to evaluate the difference in immune cell content in DRG tissues and found higher levels of Th1 cells infiltration in the spinal nerve injury (SNI) group compared to the control (Con) group. Subsequently, the module genes associated with Th1 cells infiltration were identified using weighted gene co-expression network analysis (WGCNA) [13], the module genes most associated with Th1 cells infiltration were identified. Next, the expression characteristics of crucial genes were analyzed at the single-cell level. Finally, the immune-related key genes were validated, and the potential mechanism of Abca1 was explored in the occurrence and development of NP to provide critical guidance in developing effective immunotherapy strategies.

Results

Selecting and preprocessing data

The flowchart of analysis steps is shown in Fig. 1. All data were downloaded from the GEO database and sorted in Table 1. First, we merged GSE102721 and GSE149770 datasets to obtain 9 Con and 10 SNI samples. Then, the batch effect between the data was eliminated. The clustering distribution was uniform after the normalization of the dataset, indicating data reliability (Fig. 2A and B).

Fig. 1
figure 1

Flow chart of the analysis step in this study. The datasets were downloaded from the GEO database and the differentially expressed genes (DEGs) that were significantly correlated with NP were screened. Enrichment analysis was performed to describe the function of DEGs. ssGSEA algorithm was used to identify significantly different immune cells. WGCNA package was applied to screen Hub genes related to immune cells. Then, Hub genes were verified through external datasets, and spearman was used to analyze the correlation between Hub genes and immune cells. Next, the expression patterns of these Hub genes were detected at the single-cell level. Finally, the expression level of the Hub gene and its related functions will be verified by in vitro and in vivo experiments

Table 1 The enrolled datasets in the current study
Fig. 2
figure 2

The differentially expressed genes (DEGs) in Dorsal root ganglion (DRG) tissues and enriched analysis. A UMAP diagram, the sample distribution of each dataset before the batch effect is removed, and each point represents a sample; (B) the sample distribution of each dataset after removing the batch effect; (C) Heat map of DEGs, the color represents expression level, the higher expression level, the darker the color (red for up-regulation and blue for down-regulation). The tree cluster on the left represents significant clustering results in different samples; (D) Volcano map of DEGs, red represents significantly up-regulated genes, blue represents significantly down-regulated genes, and gray represents genes with no significant differences. E Gene Ontology (GO) enrichment analysis bubble diagram, the ordinate is the name of GO and the abscissa is the proportion of genes. The circle size represents the count of DEGs enriched in the GO (the larger the circle, the more DEGs are enriched), and the color represents enrichment significance (the redder the color is, the more significant the DEGs are enriched on this GO). GO analysis includes biological processes (BP), cellular components (CC) and molecular functions (MF); (F) Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analysis bubble diagram: the circle size represents the count of DEGs enriched to the pathway, the color represents significance

Screening for DEGs and enrichment analysis

A total of 163 DEGs (114 upregulated and 49 downregulated) were identified between the Con and SNI groups from the merged gene expression matrix, as shown in the volcano map and heat map (Figs. 2C and D). Figures 2E and F illustrated the results of the functional correlation analysis of GO and KEGG. The results of GO analysis suggested that DEGs are involved in three cellular functions: biological processes (BP), cell components (CC), and molecular function (MF) (Fig. 2E). The main changes in BP were the metabolism and activation of amino acids, such as icosanoid metabolic process, response to glucocorticoid, and response to corticosteroid. The primary variations in CC involved synaptic membrane, receptor complex, and postsynaptic membrane. The significant changes in MF included G protein-coupled peptide receptor activity, monooxygenase activity, peptide receptor activity, and fatty acid binding. KEGG pathway analysis results (Fig. 2F) indicated that DEGs are significantly enriched in neuroactive ligand-receptor interaction and calcium, PI3K-Akt, MAPK, and TGF-beta signaling pathways.

Screening of significantly different immune cells

The ssGSEA algorithm was used to analyze the specific immune cell types infiltrated into DRG tissues and elucidate their immune microenvironment. The boxplot showed differences in immune cell infiltration in DRG tissues (Fig. 3A). Compared to the Con group, the infiltration level of type 1 T helper (Th1) cells (P < 0.05), CD56bright natural killer (NK) cells (P < 0.05), and neutrophils (P < 0.01) were higher in the SNI group. Therefore, we defined these three types of significantly different immune cells as the trait data analyzed by WGCNA. In addition, Fig. 3B shows the correlation between 28 immune cells. Correlation heat maps showed that Th1 cells were positively correlated with gamma delta T cells and regulatory T cells. CD56bright NK cells were positively associated with central memory CD8 T cells, gamma delta T cells, and Tfh cells. A negative correlation was established between the activated B cells and effector memory CD8 T cells.

Fig. 3
figure 3

Immune infiltration landscape in DRG tissue and screened Hub genes related to immune traits by WGCNA. A Box chart, the difference of immune cell infiltration between the control group and the SNI group (symbol "*", "**", and "ns" respectively represent p-value < 0.05, p-value < 0.01, and non-significance.); (B) Heat map of the correlation of 28 types immune cells, the pie chart will be displayed only when the correlation P < 0.05, color represents correlation, color depth represents correlation significance. The darker the red represents the more significant the positive correlation, the darker the blue represents the more significant the negative correlation; (C) Network topology analysis under various soft threshold power. When the power value was 4, the independence degree was ≥ 0.85. (On the left: The X-axis shows the soft threshold power. The Y-axis reflects the fit index of the scaled-free topological model. On the right: The X-axis shows the soft threshold power. The Y-axis represents average connectivity (degree)); (D) The clustering tree and co-expression network of co-representation network modules are constructed based on 1-TOM matrix; (E) Correlation of modules (column) with immune traits (row). Red represents positive correlation, blue, negative correlation. Each unit includes correlation coefficient and p value. F Scatter plot showed the correlation between gene significance (GS > 0.5) and module membership (MM > 0.8). G Venn diagram, the intersection genes between module genes and DEGs

Identifying immune-cell infiltration-related genes

During WGCNA analysis, we chose β = 4 (scale-free R2 = 0.85) to construct a scale-free network (Fig. 3C). The modules with similarities > 0.75 were merged, and 12 modules were retrieved (Fig. 3D). Figure 3E illustrates the correlation between each module and immune traits. We found that Th1 cells were highly correlated with blue modules (Cor = 0.942, P = 1.57E-06). CD56bright NK cells were highly correlated with turquoise module (Cor = 0.917, P = 1.05E-05). Neutrophils were highly correlated with the turquoise module (Cor = 0.693, P = 0.008). Numerous studies have shown that the immune response mediated by Th1 cells is closely related to NP. To clarify the potential role and mechanism of Th1 cells in NP, we selected the blue module with the highest correlation with Th1 cells as the core module for subsequent analysis. According to the screening criteria, GS > 0.5 and MM > 0.8, 305 genes were identified in the blue module (Fig. 3F). Next, we constructed a Venn plot to show the intersection genes of the core module and DEGs (Fig. 3G). These intersection genes, including Abca1, Bach1, Cckbr, Col5a3, and Fyb, were selected as the key candidates.

Identifying and verifying key genes

The expression of the candidate key genes was verified using an external dataset (GSE24982). Figures 4A and B show that Abca1 and Fyb were significantly upregulated in the SNI group of the GSE24982 dataset, and the receiver operating characteristic (ROC) curves were drawn to evaluate the verification reliability (Figs. 4C and D). However, it was found that Bach1 did not exist in the GSE24982 dataset, while Cckbr and Col5a3 did not differ significantly between the experimental and the control group. Therefore, follow-up studies mainly focused on Abca1 and Fyb. Eventually, Abca1 [area under the curve (AUC) = 0.885, 95% confidence interval (CI): 0.768–0.975] and Fyb (AUC = 0.998, 95% CI: 0.985–1.000) were identified as the key genes for Th1 cell infiltration in DRG. Subsequently, Spearman’s correlation analysis was performed to verify the correlation between Abca1, Fyb and immune cell infiltration. The results showed a correlation between Abca1 and Th1 T cells, regulatory T cells, CD56bright NK cells, type 17 Th cells, gamma delta T cells, immature B cells, neutrophils, and activated dendritic cells was statistically significant (Fig. 4E). Notably, Abca1 was significantly correlated with Th1 cells (R = 0.81, P = 1E-05) (Fig. 4F). Moreover, the correlation between Fyb and Th1 cells, regulatory T cells, gamma delta T cells, myeloid-derived suppressor cells (MDSCs), type 17 Th cells, neutrophils, activated dendritic cells, and CD56bright NK cells was statistically significant (Fig. 4G). Especially, Fyb was significantly correlated with Th1 cells (R = 0.8, P = 1.7E-05) (Fig. 4H).

Fig. 4
figure 4

Verification the correlation between Hub gene expression and immunity. A, B Hub gene expression level was verified by external dataset; (C, D) ROC curve shows validation efficiency; (E) Lollipop plot showed the correlation between Abca1 and infiltrating immune cells. (The larger the circle size, the stronger the correlation; p < 0.05 was considered statistically significant.); (F) The correlation between Abca1 and Th1 cells. G Lollipop plot showed the correlation between Fyb and infiltrating immune cells. (The larger the circle size, the stronger the correlation; p < 0.05 was considered statistically significant.); (H) The correlation between Fyb and Th1 cells

scRNA-Seq data revealed high cellular heterogeneity in DRG tissues

We downloaded mouse scRNA sequencing data from Sham and SNI models to determine the single-cell level pattern in DRG tissues. First, we performed quality control on the gene expression matrix (Fig. 5A). Then, the RNA-seq data were normalized, and 15 PCs (P < 0.05) were screened for subsequent analysis (Fig. 5B). Then t-distributed random neighborhood embedding (t-SNE) was used for unsupervised cell clustering analysis (Fig. 5C). Figure 5D shows that DRG tissue is isolated into nine distinct clusters, including schwann, SGC, endothelial, macrophage, pericyte, fibroblast, neuron, and myeloid. Next, we tested the expression patterns of candidate key genes in these cell clusters and verified their expression patterns at the single-cell level. As expected, they were significantly expressed in immune-related cells; the key gene Abca1 (Figs. 5E and F) and Fyb (Figs. 5G and H) was highly expressed in macrophages.

Fig. 5
figure 5

Flow analysis of single cell data and expression pattern of Hub gene at single cell level. A, B Gene filtering and PCA clustering of the gene expression matrix. C, D t-SNE projections and cell subset annotation of SNI model DRG tissue. E, F Expression pattern of Abca1 at the single cell level, shown in t-SNE figures and violin maps. (G, H) Expression pattern of Fyb at the single cell level, shown in t-SNE figures and violin maps

Abca1 was highly expressed in mouse DRG tissue and significantly affected M1 polarization of macrophages

The expression of Abca1 was further validated in animal and cell experiments. RAW264.7 cells were stimulated with 1µg/mL LPS at different time points to detect Abca1 and the corresponding inflammatory indicators (Fig. 6A). As expected, the results showed that Abca1 expression was significant at 6 h after LPS stimulation, and the inflammatory factors, TNF-α and IL-1β of macrophage M1 polarization were significantly increased (Figs. 6B and C). Subsequently, qRT-PCR results (Fig. 6D) showed that Abca1 expression level increased significantly in SNI model DRG tissues. However, there was no significant difference in the expression level of Fyb (Fig. 6D). Therefore, our subsequent studies mainly focus on Abca1. To further explore the effect of Abca1 on macrophage function, we knocked down Abca1 in macrophages. Figure 6E shows the siRNA knockdown efficiency; finally, siRNA4 was selected for subsequent assays. The results showed that when Abca1 expression was downregulated, the expression of macrophage inflammatory factors, TNF-α, and IL-1β was decreased (Fig. 6F, G, and H). Therefore, we speculated that Abca1 affects the polarization function of macrophage M1.

Fig. 6
figure 6

In vivo and in vitro experimental verification. A, B, C Expression levels of Abca1 and inflammatory factors in Raw264.7 cells stimulated by (LPS 1µg/mL) at different times. D Expression levels of Abca1 and Fyb in SNI model DRG tissue. E Abca1 silencing efficiency. F, G, H Expression levels of inflammatory factors after Abca1 was down-regulated

Discussion

The pathogenesis of NP is complex, but effective clinical treatment is lacking [14]. The interference strategies of immune cell molecular function and signaling have shown great therapeutic potential in many diseases, yet their mechanism of action in NP is poorly understood [15, 16]. The present study aimed to analyze the patterns of immune cell infiltration in DRG tissues in NP and identify the key biomarkers. T cells act on intracellular signal transduction pathways by secreting various cytokines and expressing multiple transcription factors to participate in the occurrence and development of NP [17,18,19]. Subsequently, six hub genes (Abca1, Bach1, Cckbr, Col5a3, and Fyb) were closely related to Th1 cells infiltration, indicating their role in the occurrence and progression of NP. Among the candidate Th1 infiltration-related genes, Abca1 was validated as a potential key biomarker and target of NP. Novel hub genes and pattern recognition of immune infiltration have expanded the understanding of NP pathogenesis.

In this study, GO analysis (Figs. 2E) showed that the SNI group was significantly enriched in glucocorticoid metabolism, synaptic component composition, G protein-coupled peptide receptor activity, and fatty acid linkage compared to the control group. Although glucocorticoids are classic anti-inflammatory drugs, activation of central glucocorticoid receptors might aggravate hippocampal neuronal death and increase the neurotoxicity of CNS inflammation [20, 21]. In addition, previous studies have shown that inflammatory cytokines regulate synaptic structure and function after peripheral nerve injury, which might underlie NP and memory deficits caused by peripheral nerve injury [18, 19]. Jiang et al. [22] demonstrated that GPR151, a Gβγ-coupled receptor, induces ERK-dependent neuroinflammatory response and is involved in the maintenance of trigeminal neuralgia. Therefore, GPR151 may be a potential drug target for treating trigeminal neuralgia. Furthermore, a recent study showed that ω-6 polyunsaturated fatty acids exacerbate preclinical inflammation and induce reversible peripheral nerve dysfunction resulting in pain [23]. Moreover, KEGG analysis (Figs. 2F) showed that neuroactive ligand-receptor effects [24], calcium signaling pathway [25], PI3K-Akt signaling pathway [26], MAPK signaling pathway [27], and TGF-β signaling pathway differ significantly between the SNI and control groups. These signaling pathways have been studied extensively with respect to inflammatory effects and NP. In conclusion, the GO and KEGG results support the role of inflammation in NP.

In this study, the infiltration level of Th1 cells was significantly higher in the SNI group than in the Con group in DRG tissue. Reportedly, CD4 T lymphocytes (mainly Th1 cells) infiltrating the spinal cord develop nerve injury-induced behavioral hypersensitivity responses [28]. WGCNA indicated that the blue module was enriched for genes expressed during Th1 cell infiltration. Furthermore, the intersection genes of the blue module and DEGs were considered hub genes related to Th1 infiltration (Abca1, Bach1, Cckbr, Col5a3, and Fyb). Finally, we identified Abca1 and Fyb as the key genes closely associated with Th1 cell infiltration in NP pain through external data set validation. Correlation analysis between Abca1 and immune cells showed that Abca1 was significantly correlated with Th1 cells (R = 0.81, P = 1E-05). Interestingly, single-cell analysis revealed that Th1 infiltration-related genes (Abca1 and Fyb) were expressed in macrophages. Reportedly, the non-neuronal interactions in DRG tissue influence the development of NP [29]. Davoli-Ferreira et al. [30] demonstrated that regulatory T cells relieve NP by inhibiting the Th1 response at the site of peripheral nerve injury. Although the results of immune infiltration showed no significant difference in the infiltration level of macrophages between the SNI and Con groups, macrophages have extremely high plasticity and complex heterogeneity, so it is essential to understand the molecular characteristics of macrophages in different times and spaces and different lineages. Studies have shown that there are two sources of macrophages: macrophages (MDMs) derived from monocytes; the other are macrophages (TRMs) that reside in tissues [31]. Miriam Merad [32] study found that TRMs highly expressed Ccl17 and Tgfb1 and proved that these genes could promote the recruitment, differentiation, and amplification of Treg cells, resulting in a significant increase in the number of regulatory T (Treg) cells in tissues. Undoubtedly, it is substantial to trace the origin of macrophages with high Abca1 expression and explore how it affects Th1 cell infiltration. Therefore, we hypothesized that Abca1 promotes Th1 cells infiltration by affecting macrophage function in DRG tissues, thereby inducing NP.

Abca1 is a member of the ATP-binding cassette transporter (ABC) superfamily and a crucial cell-surface protein promoting cholesterol efflux. Glial cells express TLR4, mediating the secretion of inflammatory cytokines, chemokines, and bioactive lipids [33, 34]. Neuroinflammation mediated by infiltrating immune cells in the spinal cord and DRGs is an essential component of neuropathy [35]. However, the mechanism by which Th1 cells infiltration in DRG tissue induces abnormal pain after peripheral nerve injury has not been studied extensively. The current results showed that the expression level of Abca1 was significantly increased in DRG tissues of the SNI model. Subsequently, after stimulation of RAW264.7 cells with 1 g/mL LPS, the levels of TNF-α and IL-1β were increased significantly; however, this effect was reversed when Abca1 was knocked down in macrophages. Therefore, we proposed that Abca1 affects the function of M1 polarization in macrophages. Niehaus et al. [36] found that macrophages in DRG are critical in initiating and maintaining mechanical hypersensitivity in NP in mice. Furthermore, interferon-gamma (IFN-γ) is a potent macrophage activator [37], and thus, the mechanism underlying the interaction between macrophages, Th1 cells and inflammatory effect in DRG leading to NP needs an in-depth exploration. Therefore, the underlying mechanism may be that Abca1 induces M1 activation of macrophages to secrete proinflammatory factors, which in turn promotes naive CD4 T cells to differentiate into Th1 cells and enhances Th1 cell-derived IFN-γ-mediated inflammatory effects.

Reportedly, T cells infiltrate the damaged sciatic nerve after CCI, and passive transfer of Th1 T cells can restore the nerve sensitivity in athymic nude mice [38]. Costigan et al. [10] demonstrated a functional role of IFN-γ signaling in generating pain-like hypersensitivity responses after nerve injury in adults. Intrathecal injection of IFN-γ into the spinal cord induces hyperalgesia in naive animals [39], stimulates microglia via IFN-γR, interferes with the signaling, and inhibits neural mechanical hypersensitivity [40]. Despite the close association of Th1 cytokines with nerve injury-induced pain, the exact contribution of Th1 cells remains to be determined. An in-depth study of immune cell crosstalk would clarify the specific role of infiltrating CD4 T lymphocytes in developing NP and guide the designing therapies for the cytokines/mechanisms underlying NP. Therefore, intervention strategies targeting Th1 cells infiltration into the DRG might contribute to the establishment of anti-NP hypersensitivity responses in humans.

Neuroimmune signaling might underlie the development of abnormal pain after nerve injury [41, 42]. Our study showed that Abca1 siRNA treatment reduced the mRNA levels of proinflammatory mediators, including TNF-α, and IL-1β, in cultured RAW264.7 cells activated by LPS. Interestingly, Abca1 was closely associated with Th1 cell infiltration in the model of peripheral nerve injury. Therefore, the effect of Abca1 on macrophage function and mediating Th1 cell infiltration will further expand Neuroimmune signaling communication in neuropathic pain; the underlying mechanism is shown in Fig. 7. Developing targeted therapies for peripheral NP is essential by using immune infiltration-related genes as starting points.

Fig. 7
figure 7

Hypothesis mechanism diagram (Created with BioRender.com)

Materials and Methods

Data acquisition and preprocessing

The bulk datasets (GSE102721, GSE149770, and GSE24982) and single-cell transcriptome data (GSE174430) were obtained from Gene Expression Omnibus (GEO) database (https://www.ncbi.nlm.nih.gov/geo/). First, the practical extraction and reporting language (Perl) (https://www.perl.org/get.html) was used for quick and accurate processing of the file paths that required R package analysis. Then, R language (Rv3.6.3 and Rv4.1.1) (https://www.r-project.org/) was used for data analysis. Next, two gene expression matrices were merged using the “limma” package [43], and batch effects were removed for subsequent analysis. Finally, the “ggplot2” package [44] was used to draw UMAP to visualize the results of batch effect removal.

Difference and enrichment analysis

The “limma” package was applied to screen the differentially expressed genes (DEGs) in the Sham and SNI groups. P-value < 0.05 and |log2 FC|> 0.5 defined the significantly different genes. The “ggplot2” package drew the heat map and the volcano plot. The “clusterProfiler” [45], “org.Mm.eg.db”, and “ggplot2” packages were applied for the analysis of Gene Ontology (GO) [46] and Kyoto Encyclopedia of Genes and Genomes (KEGG) [47, 48]. The GO analysis included biological process (BP), molecular function (MF), and cellular component (CC). KEGG analysis was used to determine the pathways of biological molecular interaction. A false-discovery rate (FDR) < 0.25 and P < 0.05 were used to screen out significant functional enrichment.

Evaluation of immune cell infiltration

The ssGESA algorithm was applied to explore the level of immune cell infiltration using packages “GSVA” [49] and “GSEABase.” In addition, the “ggplot2” package was applied to draw a boxplot and illustrate the differences in infiltrating immune cells. Furthermore, 28 types of infiltrating immune cells were plotted using the “corrplot” package (https://github.com/taiyun/corrplot).

WGCNA

WGCNA [13] was employed to identify the immune cell infiltration-related module genes and explore the correlation between the immunophenotypes and module genes in the network. The top 5000 genes with median absolute deviation were screened out for subsequent analysis. First, the outliers of gene expression matrices were filtered by hierarchical cluster analysis. β = 4 (scale-free R2 = 0.85) was selected to construct a scale-free network. Then, a dynamic hybrid cutting technology was used to construct the gene modules, and a cluster dendrogram was drawn (minModuleSize = 50, mergeCutHeight = 0.25, the colors represent different modules). Next, the correlation heatmap between modules and immune traits was drawn; rows represented modules, columns represented traits, and the correlation and P-value were displayed in boxes. The correlation between gene expression and immune characteristics was determined by screening gene significance (GS > 0.5) and module membership (MM > 0.8). Finally, the intersection of genes of interesting modules and DEGs constituted the candidate key genes. The “VennDiagram” package [50] was used to draw the Venn plot and show the intersection results.

Identification and verification of key genes

The differential expression of candidate key genes was verified through an external dataset, and the key genes with the highest significance were determined. Next, the “ggpubr” and “ggExtra” packages were applied for Spearman’s correlation analysis on key genes and infiltrating immune cells, and a lollipop plot was used to visualize the analysis results.

Analysis at the single-cell level

DRG tissue-derived single-cell RNA sequencing dataset in the SNI model (GSE174430) was downloaded from the GEO database. “Seurat” [51] package was used for downstream analysis. The cells were filtered with nFeature_RNA > 500 and nFeature_RNA < 5000 and percent.mt < 5%. “LogNormalize” method was utilized to normalize and scale gene expression. Then, PCA was used to identify the major principal components (PCs), and the JackStraw and ScoreJackStraw functions were used to visualize the P-value distribution. “Harmony” package was used for batch correction to avoid batch effects, and the FindClusters function was used to classify the cells into eight clusters with a resolution of 0.1. The cell types were manually annotated by cellMarker and panglaoDB. Finally, the expression patterns of key genes were identified and visualized by t-SNE and VlnPlot.

Animals and models

A total of 12 healthy 8-week-old male mice were purchased from the Experimental Animal Center of Guangzhou university of Chinese Medicine and randomly and equally divided into two groups (Sham and SNI). Before any experiment, all animals were acclimatized for one week. First, SNI model [52] was established, and 10% chloral hydrate (0.3 mL/100 g) was injected intraperitoneally for anesthesia. The sciatic nerve and its branches were then exposed. Next, the tibial and the common peroneal nerve were ligated with 6–0 non-invasive sutures and closely connected with distal resection, leaving the sural nerve intact. In the Sham group, the sciatic nerve and its branches were exposed but not ligated.

Cell culture and inflammatory model

Mouse macrophage cell line RAW264.7 was cultured in DMEM/High Glucose with 10% fetal bovine serum (FBS), 100 U/mL of penicillin, and 100 µg/mL of streptomycin and maintained at 37 °C in the atmosphere of 95% humidity and 5% CO2. RAW 264.7 cells were plated at a density of 1 × l05 cells/mL/ well in 6-well plates and divided into four groups. Lipopolysaccharide (LPS) was solubilized in phosphate-buffered saline (PBS). The control group was treated only with DMEM maintenance medium. The LPS group was stimulated at 6, 12, and 24 h, respectively, and the optimal time points were selected for subsequent studies.

Identification of efficient Abca1 siRNA

The siRNA of Abca1 was designed and synthesized by Tsingke Biotechnology (Beijing, China), and the sequence was as follows: siAbca1-Mouse (#1): 5’-UUGAUGAGCCUGACUUCUGTT-3’ and 5’-CAGAAGUCAGGCUCAUCAATT-3’; siAbca1-Mouse (#2): 5’-UAGUUGUUAUCCUCGUACCTT-3’ and 5’-GGUACGAGGAUAACAACUATT-3’; siAbca1-Mouse (#3): 5’-AGUAGAUCUUGGAAGGGACTT-3’ and 5’-GUCCCUUCCAAGAUCUACUTT-3’; siAbca1-Mouse (#4): 5’-AUGACAACCUUGGAUCCACTT-3’ and 5’-GUGGAUCCAAGGUUGUCAUTT-3’. RAW 264.7 cells were plated at a density of 1 × l05 cells/mL/well in 6-well plates and divided into five groups. According to the manufacturer's instructions, when the cells had grown to 60–70%, the cells were transfected using the ZetaLife Advanced DNA/RNA Transfection reagent. Total RNA of cells was extracted to detect the expression of Abca1, and siRNA with the highest knockdown efficiency was screened out for subsequent studies.

RNA extraction and quantitative real-time polymerase chain reaction (qRT-PCR)

Trizol reagen (Beyotime, R0016) was used to extract total RNA from cells and DRG tissue. Reverse transcription was performed using Prime Script RT reagent kit (Perfect Real Time). β-actin was used as an internal reference. RT-qPCR assays were carried out using Light Cycle instrument, and the relative expression level was analyzed using the 2^-(∆∆Ct) formula. The primer sequences are listed in Table 2.

Table 2 Primers used for qRT-PCR

Statistical analysis

All statistical analyses were conducted in the R language (Rv3.6.3 and Rv4.1.1). All statistical tests were bilateral, and P < 0.05 indicated a statistically significant difference.

Conclusion

Herein, whole transcriptome and single-cell sequencing data were combined to investigate the pattern of immune infiltration in DRG tissues and the key genes in NP. The current study showed that Th1 infiltration in DRG tissue was involved in the pathological process of NP induced by peripheral nerve injury. The underlying mechanism may be that Abca1 affects the M1 polarization of macrophages and promotes Th1 cell infiltration in DRG tissues. However, the mechanism of neuroimmune inflammation is rather complex and needs further investigation. This study provided novel ideas and potential targets for the pathogenesis and treatment of NP.

Availability of data and materials

Original contributions are included in the article/supplementary materials. Publicly available datasets can be found here:

https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE102721

https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE149770

https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE24982

https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE174430

References

  1. Jensen TS, Baron R, Haanpaa M, et al. A new definition of neuropathic pain. Pain. 2011;152(10):2204–5.

    Article  PubMed  Google Scholar 

  2. von Hehn CA, Baron R, Woolf CJ. Deconstructing the neuropathic pain phenotype to reveal neural mechanisms. Neuron. 2012;73(4):638–52.

    Article  Google Scholar 

  3. van Hecke O, Austin SK, Khan RA, Smith BH, Torrance N. Neuropathic pain in the general population: a systematic review of epidemiological studies. Pain. 2014;155(4):654–62.

    Article  PubMed  Google Scholar 

  4. Widerstrom-Noga E. Neuropathic pain and spinal cord injury: phenotypes and pharmacological management. Drugs. 2017;77(9):967–84.

    Article  PubMed  Google Scholar 

  5. Malcangio M. Role of the immune system in neuropathic pain. Scand J Pain. 2019;20(1):33–7.

    Article  PubMed  Google Scholar 

  6. Grace PM, Hutchinson MR, Maier SF, Watkins LR. Pathological pain and the neuroimmune interface. Nat Rev Immunol. 2014;14(4):217–31.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  7. McMahon SB, La Russa F, Bennett DL. Crosstalk between the nociceptive and immune systems in host defence and disease. Nat Rev Neurosci. 2015;16(7):389–402.

    Article  CAS  PubMed  Google Scholar 

  8. Ji RR, Chamessian A, Zhang YQ. Pain regulation by non-neuronal cells and inflammation. Science. 2016;354(6312):572–7.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  9. Inoue K, Tsuda M. Microglia in neuropathic pain: cellular and molecular mechanisms and therapeutic potential. Nat Rev Neurosci. 2018;19(3):138–52.

    Article  CAS  PubMed  Google Scholar 

  10. Costigan M, Moss A, Latremoliere A, et al. T-cell infiltration and signaling in the adult dorsal spinal cord is a major contributor to neuropathic pain-like hypersensitivity. J Neurosci. 2009;29(46):14415–22.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  11. Yu X, Liu H, Hamel KA, et al. Dorsal root ganglion macrophages contribute to both the initiation and persistence of neuropathic pain. Nat Commun. 2020;11(1):264.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  12. Newman AM, Liu CL, Green MR, et al. Robust enumeration of cell subsets from tissue expression profiles. Nat Methods. 2015;12(5):453–7.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  13. Langfelder P, Horvath S. WGCNA: an R package for weighted correlation network analysis. BMC Bioinformatics. 2008;9:559.

    Article  PubMed  PubMed Central  Google Scholar 

  14. Alles SRA, Smith PA. Etiology and pharmacology of neuropathic pain. Pharmacol Rev. 2018;70(2):315–47.

    Article  CAS  PubMed  Google Scholar 

  15. Gentles AJ, Newman AM, Liu CL, et al. The prognostic landscape of genes and infiltrating immune cells across human cancers. Nat Med. 2015;21(8):938–45.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  16. Jiang T, Shi J, Dong Z, et al. Genomic landscape and its correlations with tumor mutational burden, PD-L1 expression, and immune cells infiltration in Chinese lung squamous cell carcinoma. J Hematol Oncol. 2019;12(1):75.

    Article  PubMed  PubMed Central  Google Scholar 

  17. Sorge RE, Mapplebeck JC, Rosen S, et al. Different immune cells mediate mechanical pain hypersensitivity in male and female mice. Nat Neurosci. 2015;18(8):1081–3.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  18. Liu Y, Zhou LJ, Wang J, et al. TNF-alpha differentially regulates synaptic plasticity in the hippocampus and spinal cord by microglia-dependent mechanisms after peripheral nerve injury. J Neurosci. 2017;37(4):871–81.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  19. Chirila AM, Brown TE, Bishop RA, Bellono NW, Pucci FG, Kauer JA. Long-term potentiation of glycinergic synapses triggered by interleukin 1beta. Proc Natl Acad Sci U S A. 2014;111(22):8263–8.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  20. Takeda K, Sawamura S, Sekiyama H, Tamai H, Hanaoka K. Effect of methylprednisolone on neuropathic pain and spinal glial activation in rats. Anesthesiology. 2004;100(5):1249–57.

    Article  CAS  PubMed  Google Scholar 

  21. Yin X, Chen S, Eisenbarth SC. Dendritic cell regulation of T helper cells. Annu Rev Immunol. 2021;39:759–90.

    Article  CAS  PubMed  Google Scholar 

  22. Jiang BC, Zhang J, Wu B, et al. G protein-coupled receptor GPR151 is involved in trigeminal neuropathic pain through the induction of Gbetagamma/extracellular signal-regulated kinase-mediated neuroinflammation in the trigeminal ganglion. Pain. 2021;162(5):1434–48.

    Article  CAS  PubMed  Google Scholar 

  23. Boyd JT, LoCoco PM, Furr AR, et al. Elevated dietary omega-6 polyunsaturated fatty acids induce reversible peripheral nerve dysfunction that exacerbates comorbid pain conditions. Nat Metab. 2021;3(6):762–73.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  24. Pathirathna S, Todorovic SM, Covey DF, Jevtovic-Todorovic V. 5alpha-reduced neuroactive steroids alleviate thermal and mechanical hyperalgesia in rats with neuropathic pain. Pain. 2005;117(3):326–39.

    Article  CAS  PubMed  Google Scholar 

  25. Liu Q, Chen W, Fan X, et al. Upregulation of interleukin-6 on Cav3.2 T-type calcium channels in dorsal root ganglion neurons contributes to neuropathic pain in rats with spinal nerve ligation. Exp Neurol. 2019;317:226–43.

    Article  CAS  PubMed  Google Scholar 

  26. Jiang SP, Zhang ZD, Kang LM, Wang QH, Zhang L, Chen HP. Celecoxib reverts oxaliplatin-induced neuropathic pain through inhibiting PI3K/Akt2 pathway in the mouse dorsal root ganglion. Exp Neurol. 2016;275(Pt 1):11–6.

    Article  CAS  PubMed  Google Scholar 

  27. Wang L, Yin C, Liu T, et al. Pellino1 regulates neuropathic pain as well as microglial activation through the regulation of MAPK/NF-kappaB signaling in the spinal cord. J Neuroinflammation. 2020;17(1):83.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  28. Zhang X, Wu Z, Hayashi Y, Okada R, Nakanishi H. Peripheral role of cathepsin S in Th1 cell-dependent transition of nerve injury-induced acute pain to a chronic pain state. J Neurosci. 2014;34(8):3013–22.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  29. Scholz J, Woolf CJ. The neuropathic pain triad: neurons, immune cells and glia. Nat Neurosci. 2007;10(11):1361–8.

    Article  CAS  PubMed  Google Scholar 

  30. Davoli-Ferreira M, de Lima KA, Fonseca MM, et al. Regulatory T cells counteract neuropathic pain through inhibition of the Th1 response at the site of peripheral nerve injury. Pain. 2020;161(8):1730–43.

    Article  CAS  PubMed  Google Scholar 

  31. Hashimoto D, Chow A, Noizat C, et al. Tissue-resident macrophages self-maintain locally throughout adult life with minimal contribution from circulating monocytes. Immunity. 2013;38(4):792–804.

    Article  CAS  PubMed  Google Scholar 

  32. Casanova-Acebes M, Dalla E, Leader AM, et al. Tissue-resident macrophages provide a pro-tumorigenic niche to early NSCLC cells. Nature. 2021;595(7868):578–84.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  33. Bruno K, Woller SA, Miller YI, et al. Targeting toll-like receptor-4 (TLR4)-an emerging therapeutic target for persistent pain states. Pain. 2018;159(10):1908–15.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  34. Papageorgiou IE, Lewen A, Galow LV, et al. TLR4-activated microglia require IFN-gamma to induce severe neuronal dysfunction and death in situ. Proc Natl Acad Sci U S A. 2016;113(1):212–7.

    Article  CAS  PubMed  Google Scholar 

  35. Lees JG, Makker PG, Tonkin RS, et al. Immune-mediated processes implicated in chemotherapy-induced peripheral neuropathy. Eur J Cancer. 2017;73:22–9.

    Article  CAS  PubMed  Google Scholar 

  36. Niehaus JK, Taylor-Blake B, Loo L, Simon JM, Zylka MJ. Spinal macrophages resolve nociceptive hypersensitivity after peripheral injury. Neuron. 2021;109(8):1274-1282 e1276.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  37. Schroder K, Hertzog PJ, Ravasi T, Hume DA. Interferon-gamma: an overview of signals, mechanisms and functions. J Leukoc Biol. 2004;75(2):163–89.

    Article  CAS  PubMed  Google Scholar 

  38. Moalem G, Xu K, Yu L. T lymphocytes play a role in neuropathic pain following peripheral nerve injury in rats. Neuroscience. 2004;129(3):767–77.

    Article  CAS  PubMed  Google Scholar 

  39. Vikman KS, Siddall PJ, Duggan AW. Increased responsiveness of rat dorsal horn neurons in vivo following prolonged intrathecal exposure to interferon-gamma. Neuroscience. 2005;135(3):969–77.

    Article  CAS  PubMed  Google Scholar 

  40. Tsuda M, Masuda T, Kitano J, Shimoyama H, Tozaki-Saitoh H, Inoue K. IFN-gamma receptor signaling mediates spinal microglia activation driving neuropathic pain. Proc Natl Acad Sci U S A. 2009;106(19):8032–7.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  41. Grace PM, Tawfik VL, Svensson CI, Burton MD, Loggia ML, Hutchinson MR. The Neuroimmunology of chronic pain: from rodents to humans. J Neurosci. 2021;41(5):855–65.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  42. Goncalves Dos Santos G, Delay L, Yaksh TL, Corr M. Neuraxial cytokines in pain States. Front Immunol. 2019;10:3061.

    Article  PubMed  Google Scholar 

  43. Ritchie ME, Phipson B, Wu D, et al. limma powers differential expression analyses for RNA-sequencing and microarray studies. Nucleic Acids Res. 2015;43(7):e47.

    Article  PubMed  PubMed Central  Google Scholar 

  44. Steenwyk JL, Rokas A. ggpubfigs: colorblind-friendly color palettes and ggplot2 graphic system extensions for publication-quality scientific figures. Microbiol Resour Announc. 2021;10(44): e0087121.

    Article  PubMed  Google Scholar 

  45. Yu G, Wang LG, Han Y, He QY. clusterProfiler: an R package for comparing biological themes among gene clusters. OMICS. 2012;16(5):284–7.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  46. Gene Ontology C. Gene Ontology Consortium: going forward. Nucleic Acids Res. 2015;43(Database issue):D1049-1056.

    Article  Google Scholar 

  47. Kanehisa M, Goto S. KEGG: kyoto encyclopedia of genes and genomes. Nucleic Acids Res. 2000;28(1):27–30.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  48. Kanehisa M, Furumichi M, Tanabe M, Sato Y, Morishima K. KEGG: new perspectives on genomes, pathways, diseases and drugs. Nucleic Acids Res. 2017;45(D1):D353–61.

    Article  CAS  PubMed  Google Scholar 

  49. Hanzelmann S, Castelo R, Guinney J. GSVA: gene set variation analysis for microarray and RNA-seq data. BMC Bioinformatics. 2013;14:7.

    Article  PubMed  PubMed Central  Google Scholar 

  50. Chen H, Boutros PC. VennDiagram: a package for the generation of highly-customizable Venn and Euler diagrams in R. BMC Bioinform. 2011;12:35.

    Article  Google Scholar 

  51. Butler A, Hoffman P, Smibert P, Papalexi E, Satija R. Integrating single-cell transcriptomic data across different conditions, technologies, and species. Nat Biotechnol. 2018;36(5):411–20.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  52. Decosterd I, Woolf CJ. Spared nerve injury: an animal model of persistent peripheral neuropathic pain. Pain. 2000;87(2):149–58.

    Article  PubMed  Google Scholar 

Download references

Institutional review board statement

Not applicable.

Informed consent statement

Not applicable.

Funding

This work was supported partly by National Natural Science Foundation of China (81870879 to JZ), Guangzhou Key Laboratory of Spinal Cord Level Neuropathic Pain Research (202102100005 to JZ); Guangdong Science and Technology Department (A2022086) and President's Foundation (YQ2021012) to SMT. The funders had no role in the design of the study, the collection, analysis, and interpretation of the data, the writing of the manuscript, and the decision to submit the manuscript for publication.

Author information

Authors and Affiliations

Authors

Contributions

XSZ and JZ contributed to the study conception and design, conducted data analysis, and drafted and critically revised the manuscript. JRC, CYG and YDD assisted data analysis. YC, SW and CXZ contributed to data collection. ZQL and SMT contributed to literature search. All authors contributed to the article and approved the submitted version.

Corresponding author

Correspondence to Jun Zhou.

Ethics declarations

Competing interests

The authors declare no competing interests.

Additional information

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Rights and permissions

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/. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated in a credit line to the data.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Zhang, X., Cheng, J., Deng, Y. et al. Identification and validation of biomarkers related to Th1 cell infiltration in neuropathic pain. J Inflamm 20, 19 (2023). https://doi.org/10.1186/s12950-023-00343-0

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doi.org/10.1186/s12950-023-00343-0

Keywords