Skip to main content
Download PDF

Identification of hub genes in the subacute spinal cord injury in rats

BMC Neuroscience volume 23, Article number: 51 (2022) Cite this article

Abstract

Background

Spinal cord injury (SCI) is a common trauma in clinical practices. Subacute SCI is mainly characterized by neuronal apoptosis, axonal demyelination, Wallerian degeneration, axonal remodeling, and glial scar formation. It has been discovered in recent years that inflammatory responses are particularly important in subacute SCI. However, the mechanisms mediating inflammation are not completely clear.

Methods

The gene expression profiles of GSE20907, GSE45006, and GSE45550 were downloaded from the GEO database. The models of the three gene expression profiles were all for SCI to the thoracic segment of the rat. The differentially expressed genes (DEGs) and weighted correlation network analysis (WGCNA) were performed using R software, and functional enrichment analysis and protein–protein interaction (PPI) network were performed using Metascape. Module analysis was performed using Cytoscape. Finally, the relative mRNA expression level of central genes was verified by RT-PCR.

Results

A total of 206 candidate genes were identified, including 164 up-regulated genes and 42 down-regulated genes. The PPI network was evaluated, and the candidate genes enrichment results were mainly related to the production of tumor necrosis factors and innate immune regulatory response. Twelve core genes were identified, including 10 up-regulated genes and 2 down-regulated genes. Finally, seven hub genes with statistical significance in both the RT-PCR results and expression matrix were identified, namely Itgb1, Ptprc, Cd63, Lgals3, Vav1, Shc1, and Casp4. They are all related to the activation process of microglia.

Conclusion

In this study, we identified the hub genes and signaling pathways involved in subacute SCI using bioinformatics methods, which may provide a molecular basis for the future treatment of SCI.

Peer Review reports

Introduction

Spinal cord injury (SCI) is a serious complication of spinal fracture, in which the spinal cord or cauda equina is damaged to different degrees due to the displacement of the vertebral body or the intrusion of bone fragments into the spinal canal. It is a common trauma in clinics and has the characteristics of high incidence and disability rate. SCI is the main cause of long-term physical impairment and disability, which is a huge burden on the quality of life of patients and the medical system [1,2,3,4]. Despite significant advances in the state of the art of medical care in spinal surgery, there are currently no effective treatment options for this neurological disorder, mostly limited to supportive measures [5,6,7]. Traumatic SCI is divided into primary injury and secondary injury. Primary injury refers to mechanical direct injury to the spinal cord [8]; Secondary injury refers to the subsequent pathological reaction caused by direct phases. Secondary injury includes three stages: acute, subacute, and chronic injury. Acute secondary injury occurs within 0–48 h of injury and is initially characterized by increased calcium influx, ion imbalance, lipid peroxidation, free radical production, inflammation, and edema [9]. The stage of subacute secondary injury occurs two days to two weeks after injury and is characterized by neuronal apoptosis, axonal demyelination, Waller degeneration, axonal remodeling, and glial scar formation [9]. Over time, subacute lesions developed into chronic secondary lesions characterized by glial scar maturation, capfsular formation, and axonal dieback [10].

The most important stage in the pathophysiological process of SCI involves a secondary injury caused by neuroinflammation in the lesion area, accompanied by abnormal molecular signals, vascular changes, and secondary cell dysfunction, which is uncontrolled and destructive cascade [11,12,13,14]. Gene expression and signal pathways of that inflammatory cascade in subacute injury are complex [15]. During this period, the inflammatory cascade has dual effects, not only aggravating the tissue cell damage but also serving as an important promoting factor for the plastic change after SCI [16,17,18]. How to balance these dual effects is very important for the effect of intervention and treatment. Therefore, the relevant molecular mechanisms of inflammation in the subacute phase after SCI are worthy of in-depth study.

Bioinformatics is an emerging interdisciplinary subject that combines molecular biology and information technology, which opens up a new way for the diagnosis and treatment of human diseases [19]. Gene chip, as an emerging technology, has been used for efficient and large-scale access to biological information and can widely collect disease expression profile data. Some scholars have unveiled the activation pathways and molecular targets during acute or chronic SCI by identifying differentially expressed genes (DEGs) using microarray or RNA sequencing analysis [20,21,22,23,24,25]. However, the sequencing analysis for subacute SCI is relatively rare. In this paper, the bioinformatics tools are used to analyze the data from the sham operation group and SCI group of GSE20907, GSE45006, GSE45550 in the common gene chip databases, aiming to identify the key biomarkers of abnormally expressed genes in subacute SCI and provide targets for the diagnosis and treatment of subacute SCI.

Materials and methods

Download expression matrix data

Expression matrices of GSE20907, GSE45006, and GPL45550 were downloaded from the GEO database (https://www.ncbi.nlm.nih.gov/geo/). The ages of rats in the three data sets were about 9–11 weeks. GSE20907 [26] microarray data were based on the GPL6247 platform by performing a T9-T10 laser axonotomy in female Long Evans rats, followed by a moderate spinal contusion injury by dropping 10 g of the rod from a height of 25 mm using NYU Impactor. GSE45006 [27] microarray data were based on the GPL1355 platform and operated on female Wistar rats for T6-T8 laminectomy and moderate to severe impact compression damage using a 35 g Walsh clip for 1 min. GSE45550 [28] microarray data were based on the GPL1355 platform by performing dorsal laminectomy in the thoracic vertebrae T7-T9 of female Sprague Dawley rats, and the contusion was generated by dropping 10 g cylinders onto the T8 segment of the spinal cord from a height of 25 mm. Four groups of data from sham, on the 3rd, 7th, and 14th day after operation were extracted from each expression matrix for analysis. Information on these data sets is shown in Table 1. The experimental strategy in this paper is shown in Fig. 1-a.

Table 1 Number of samples per group contained in each dataset
Full size table
Fig. 1
figure 1

Research roadmap and data preprocessing. a Research roadmap. b PCA chart of three data sets without removing batch effect. c PCA chart of three data sets after removing batch effect. d Venn diagram of DEGs in sham operation group and 3 d, 7 d, and 14 d groups after subacute SCI. In each annotation circle, Red represented the number of up-regulated genes, Blue represented the number of down-regulated genes, and Yellow represented the number of genes with the opposite trend in the intersection set. Principal component analysis, PCA

Full size image

Data preprocessing and identification of DEGs

The expression matrix was subjected to batch effect elimination and batch normalization using R software (version 4.1.0; https://www.r-project.org/) and R-package SVA [29]. The Limma package in R (Limma;http://www.bioconductor.org/packages/release/bioc/html/limma.html) was used to identify DEGs by comparing expression values between sham and subacute SCI. The corresponding P-value of the gene symbols after the t-test was used. The adjusted P < 0.05 and |logFC|> 1 were used as the selection criteria.

Weighted correlation network analysis (WGCNA)

WGCNA (v1.61; https://cran.r-project.org/web/packages/WGCNA/index.html) is a tool [30] for constructing gene co-expression networks and identifying gene clusters or modules. Thus, the WGCNA integration algorithm R (v3.4.1) was used to analyze highly relevant native gene clusters or modules for subacute SCI. There were ≥ 10 cutoff genes, cutting height = 0.85, Z‐ score ≥ 5, and stability-related stability correlation P ≤ 0.05 in this study. The connection of nodes (genes) between the two was used to calculate the data set, and genes with the correlation coefficient < 0.5 were excluded. The conservation status of WGCNA module and the traits related characteristics were analyzed.

Time dynamic clustering analysis

To characterize the dynamics of subacute SCI gene expression, R-package Mfuzz[31] (v2.1; http://cran.r-project.org/web/packages/Mfuzz/index.html) was used for temporal dynamic cluster analysis. Fuzzy C-means clustering analysis, the core algorithm of Mfuzz, can aggregate genes with similar expression patterns. The differential genes were divided into four different clusters according to the expression pattern of subacute SCI. The cluster-score of the gene indicated the similarity of each cluster.

Functional enrichment analysis.

The KEGG/GO analysis was performed using the Metascape website (http://www.metascape.org/) [32].

Protein–protein interaction network (PPI)

The PPI analysis was performed using the Metascape website (http://www.metascape.org/) [32]. The connectivity (degree) and hub nodes (genes) in PPI [33] were obtained using scale-free property to obtain. And MCODE algorithm was applied to this network, and GO enrichment analysis was applied to each MCODE network, each MCODE network being assigned a unique shape. For the hub nodes, the size of the shape represented the value of MCODE-degree. The results of PPI were imported into Cytoscape software (v3.9.0; http://www.cytoscape.org/) and further analyzed in combination with the results of temporal dynamic clustering analysis.

Comparison of expression of candidate genes and real-time polymerase chain reaction (RT-PCR)

The heat map was made using R-packages ComplexHeatmap (v3.1; https://cran.r-project.org/web/packages/ComplexHeatmap/index.html) and GGplots (v3.0; https://cran.r-project.org/web/packages/ggplots/index.html) to compare the expression levels of candidate genes.

The hub genes were selected for RT-PCR. Female Sprague Dawley rats (6–8 weeks of age, average weight 210 g) were acquired from the Laboratory Animal Center of Nantong University (Nantong, China). Rats were provided with normal food and water and housed at 20–26 ℃ under 55%-65% humidity with a 12 h/12 h artificial diurnal cycle. Rats underwent general anesthesia (20 ml/kg) by intraperitoneal injection of avertin (2, 2, 2-tribromoethanol, Sigma-Aldrich) in 0.9% saline solution and were injured by impact compression using a 35 g Walsh clip for 1 min at thoracic level T7-T9 and sacrificed (CO2 asphyxia) for spinal cord tissue removal 3, 7, and 14 days after SCI. Five rats in each group. Total RNA was isolated using TRIZOL reagent (Invitrogen, Thermo Fisher Scientific, Inc.) and reverse transcribed into cDNA according to the manufacturer's instructions. RT-PCR was performed using SYBR green dye(Takara) in a thermal cycler under the following parameters: initial denaturation step at 95 ℃ for 30 min; 40 cycles at 95 ℃ for 5 s; And 60 ℃ for 30 s. The complete experimental procedure was performed on each sample in duplicate. The mRNA primers are shown in Table 2.

Table 2 Hub genes primers used in this study
Full size table

Data analysis

SPSS 22.0 software (IBM, Armonk, NY, USA) was used for all statistical analyses. If the data both showed normality and homogeneity of variance, they were expressed as mean and standard deviation. The Unpaired Student’s t-test was used for inter-group comparison. If the data were not showed normality or homogeneity of variance, the quantitative variables were expressed as the median of the ranges and compared between groups using the Mann–Whitney Wilcoxon test, respectively. All significance levels were set to P < 0.05 and plotted using GraphPad Prism 6 (GraphPad Software Inc., CA, USA).

Results

Data preprocessing and identification of DEGs

To eliminate the batch effect of merging different datasets, microarray results from GSE20907, GSE45006, and GSE45550 were batch corrected and normalized using PCA (Fig. 1b, c), and DEGs were selected by comparing the surgical group and the subacute SCI groups (3 days /7 days /14 days) using R-package Limma. The Venn diagram was drawn using R software (Fig. 1d). There were 2414 genes in the union of DEGs detected, and 402 genes in the intersection of DEGs detected. In the sham operation group, there were 771 down-regulated genes and 972 up-regulated genes compared with the 3-day group, 477 down-regulated genes and 677 up-regulated genes compared with the 7-day group, and 398 down-regulated genes and 571 up-regulated genes compared with the 14-day group.

WGCNA of the union of DEGs

The union of 2414 DEGs was taken as the expression matrix and used as the input data for network construction. According to the rules of scale-free networks, the larger the correlation coefficient is, the more significant the scale-free property of the network is. According to the prerequisite of the approximate scale-free topology processed by WGCNA, the soft threshold power of the adjacency matrix is 9, and the standard that the square of the intrinsic gene correlation coefficient is greater than 0.85 is taken as the standard for module identification (Fig. 2a). At this time, the average connectivity was 1, indicating that the gene module was constructed according to the approximate scale-free topology standard (Fig. 2b). After the soft threshold was determined, the expression matrix of the differential genes was converted into an adjacency matrix, a topology matrix, and a dissimilarity matrix between genes. On this basis, the hierarchical clustering method is used for gene clustering, and the dynamic cutting algorithm is used for module identification of the system clustering tree. Fifteen different co-expression modules were obtained and expressed in different colors (Fig. 2c). These modules were correlated to clinical features and modules of continuous correlation were found (Fig. 2d). From the results, we could see that the correlation coefficients of the three subgroups of subacute SCI in the Turquoise module were gradually decreased (Rho3 = 0.24, Rho7 = 0.21, Rho14 = 0.19), suggesting that the genes in this module had a continuous correlation with subacute SCI. Figure 2e shows the number of genes in each module, in which the Turquoise module contained 530 genes.

Fig. 2
figure 2

WGCNA of the union of DEGs. a scale independence, b mean connectivity. The network topology analysis for adjacency matrix with different soft threshold power. Red numbers in the boxes indicate the soft thresholding power corresponding to the correlation coefficient square value(y-axis). c consensus module dendrogram was produced by clustering of 2414 genes with a variation coefficient of expression > 0.1, based on the criteria of correlation coefficient square of eigengenes above 0.85, soft threshold power of 9, the number of genes > 10, and cut height = 0.95. d Module-trait associations. Each row corresponds to a module trait gene, and each column corresponds to a trait. Red indicated a positive correlation between modular trait genes and traits, and blue indicated a negative correlation. Each cell contains the correlation coefficient Rho and the P-value in parentheses. e Pie chart of the number of genes in modules, each color representing each module. WGCNA, weighted correlation network analysis

Full size image

Functional enrichment analysis of candidate genes

398 genes in the Intersection (402 genes exclude 4 genes with different trends) and 530 genes in the Turquoise module were crossed again to obtain 206 candidate genes (Fig. 3a). The GO function and KEGG pathway of candidate genes were annotated with Metascape, and the top 20 results were shown in Fig. 3b and Table 3. We found that these genes were mainly associated with 12 GO Biological Processes, including Tumor Necrosis Factor Production, Myeloid Leukocyte Activation, Leukocyte Migration, Positive Regulation of Cell Migration, Lymphocyte Proliferation, Innate Immune Response, Wound Healing, Negative Regulation of Cytokine Production, Immune Effector Process, Regulation of Cell Adhesion, Positive Regulation of Smooth Muscle Cell Proliferation, and Regulation of Tumor Necrosis Factor-mediated Signaling Pathway; It was associated with five KEGG Pathways, including Leishmaniasis, Pertussis, Staphylococcus Aureus Infection, Proteoglycans in Cancer, and Legionellosis; It was associated with two Reactome Gene Sets, including Innate Immune System and Hemostasis; Related to one WikiPathways for IL-5 signaling pathway.

Fig. 3
figure 3

Functional enrichment analysis and PPI analysis of candidate genes. a Venn diagram analysis of genes between the Turquoise module and the Intersection. b Functional enrichment analysis for 206 candidate genes. c Pie chart of the number of genes of four clusters in time-dynamic clustering analysis of 206 candidate genes. Four different colors represent different clusters. d Change trend of the relative expression level of four clusters in time-dynamic clustering analysis of 206 candidate genes. e PPI analysis of candidate genes. Four different colors represent different clusters in the temporal dynamic clustering analysis MCODE algorithm was applied to this network, and GO enrichment analysis was applied to each MCODE network, each MCODE network being assigned a unique shape. For the hub nodes, the size of the shape represented the value of MCODE-degree. PPI, protein–protein interaction

Full size image
Table 3 The top 20 results of the GO function and KEGG pathway of candidate genes
Full size table

Temporal dynamic clustering analysis and PPI analysis of candidate genes

Temporal dynamic clustering analysis was performed on the candidate genes using R-package Mfuzz. As shown in Fig. 3c, d, the candidate genes were divided into four clusters with different expression modes. Cluster-1 represented 79 genes with continuous high expression over time after SCI, Cluster-2 represented 16 genes with continuous low expression over time after SCI, Cluster-3 represented 85 genes with high expression after SCI, and Cluster-4 represented 26 genes with low expression after SCI. PPI analysis of candidate genes was performed using the Metascape, and the results showed that 117 genes among the candidate genes served as hub nodes, as shown in Fig. 3e. There were 52 hub nodes in Cluster-1, 5 hub nodes in Cluster-2, 49 hub nodes in Cluster-3, and 11 hub nodes in Cluster-4. The top three network shapes for Log(q-value) values are Diamond, Rectangle, and Triangle. As shown in Table 4, the hub nodes are mainly related to the innate immune regulation of GO function. MCODE-1 is mainly related to the regulation of the GPCR signaling pathway, MCODE-2 is mainly related to the cross-endothelial migration of leukocytes, and MCODE-3 is mainly related to the Kit receptor signaling pathway and cytokine.

Table 4 GO enrichment analysis of MCODE network in PPI.MCODE algorithm was applied to this network, and GO enrichment analysis was applied to each MCODE network
Full size table

Comparison of hub genes expression levels and RT-PCR verification

The 117 genes were produced into heat map as shown in Fig. 4, and the 12 hub genes selected according to the MCODE-degrees and the cluster-scores are shown in Table 5, which are Itgb1, Fcgr2b, Ptprc, S100a4, Cd63, Lgals3, Lamc1, Vav1, Shc1, Casp4, Mapk12, and Vegfa. Figure 5 shows the different expression levels of the 12 hub genes in the control and experimental groups in the combined data set. Figure 6 shows the results of hub genes validation in rats by RT-PCR. The relative expression levels of the seven central genes, Itgb1, Ptprc, Cd63, Lgals3, Vav1, Shc1, and Casp4, were consistent with the expression matrix. The relative expression level of Fcgr2b increased on the 3rd day after SCI but decreased on the 7th and 14th days. The relative expression levels of S100a4 and Lamc1 were statistically significant on the 3rd and 14th days after SCI, but not on the 7th day. The relative expression levels of Mapk12 and Vegfa decreased on the 3rd day after SCI, which was consistent with the expression matrix, but increased on the 7th and 14th days after SCI, which was inconsistent with the expression matrix.

Fig. 4
figure 4

Heat map of 117 hub nodes in PPI. The numbers on the left represent clusters of time dynamic analysis. PPI, protein–protein interaction. PPI, protein–protein interaction

Full size image
Table 5 The MCODE-degrees and the cluster-scores of 12 hub genes
Full size table
Fig. 5
figure 5

Violin plot of the normalized expression levels of 12 hub genes based on the combined expression matrix

Full size image
Fig. 6
figure 6

RT-PCR validation. a–l Relative expression levels of 12 hub genes. 5 samples per group in duplicate were analyzed using RT-PCR and summarized as mean average ± SE with P < 0.05. A pairwise comparison was made between the sham operation group and 3, 7, and 14 days after subacute SCI, a Mann–Whitney Wilcoxon's test was performed. *P < 0.05, **P < 0.01. RT-PCR, real-time polymerase chain reaction; SCI, spinal cord injury; SE, standard error

Full size image

Discussion

SCI leads to motor and sensory dysfunction, and the cascade of primary injury leads to the complex cascade of secondary injury events. Studies have shown that triggering immunoreaction in the subacute phase after spinal cord injury combined with rehabilitation training is more conducive to the recovery of spinal cord function [34]. As an animal with high homology to humans, low price, and easy feeding, rats have been widely used to make spinal cord injury models. In this study, gene expression profiles of GSE20907, GSE45006, and GSE45550 were combined for the first time to conduct DEGs analysis of the subacute SCI in rats. Although the differences among the SCI models with three gene expression profiles may lead to differences in signaling pathways and functional molecules, we assume that common and collective molecular mechanisms of injury may occur even in different SCI models, and identifying these mechanisms can provide new ideas for the treatment and prognosis judgment of SCI.

In this study, we first performed WGCNA using the union of DEGs. WGCNA has been well applied in biomedical research, and its analysis has mainly focused on specific phenotypes and co-expression modules. Genes in the same module are considered to be functionally related, with higher reliability and biological significance [35, 36]. Therefore, this analysis allows the identification of biologically relevant modules and core genes that can ultimately be used as biomarkers for SCI detection or treatment. The analysis results of WGCNA showed that the Turquoise module was considered to be the module most related to the subacute SCI (3 days to 14 days), and then intersected with the intersection of DEGs to obtain 206 candidate genes. At this time, the functional enrichment analysis and PPI analysis of these candidate genes were conducted to analyze the possible potential interactions and potentially significant molecular regulatory network mechanisms between DEGs-encoded proteins, and the results showed that the candidate genes were mainly related to the production of cellular inflammatory tumor necrosis factor and the innate immune regulatory response. Studies have shown that identification of various immune responses, including activation of the complement system, induction of innate and adaptive immune responses, and antibody production [37, 38] based on GO functional analysis is the most significantly upregulated biological process during acute SCI. This study also showed the same results during subacute SCI, indicating that inflammatory injury is still dominant in subacute SCI. The results of temporal dynamic clustering analysis showed that the expression patterns of 206 candidate genes were mainly Cluster 1 which continuously showed high expression over time and Cluster 3 which always maintained a high expression level and the expression level did not significantly change over time. Besides, the main hub nodes in PPI were Cluster 1 and Cluster 3 with up-regulated DEGs, while only five genes in Cluster 2 were the hub nodes. The Log(q-value) value of the hub nodes in MCODE-1 was the largest, indicating that the GPCR signaling pathway played a very important role in subacute SCI. The peptide energy system was the most abundant network for human ligand receptor-mediated signaling, which has been widely studied in inflammation. It has been shown that the GPCR ligand promotes neuroinflammation in central nervous system diseases [39], which supports our findings. It is well known that the cross-endothelial migration of leukocytes is an indispensable link in tissue inflammation, so the results in MCODE-2 also indicate this point. In addition, the Kit receptor and cytokine-related signaling pathways in MCODE-3 have been widely studied in inflammatory mechanisms, and it has been shown to play a role in the peripheral nervous inflammation mechanism [40].

According to MCODE-degree and cluster-score, twelve hub genes were screened and verified by RT-PCR. The results showed that the relative expression levels of seven hub genes, namely Itgb1, Ptprc, Cd63, Lgals3, Vav1, Shc1, and Casp4, were consistent with microarray hybridization. It is worth mentioning that all of them were related to the activation process of microglia. As is known to all, activated microglia are the key factors for the occurrence and development of SCI [41,42,43]. Itgb1 belongs to the family of adhesion molecules, and its signal transmission is bidirectional, including cell polarity change, regulation of movement, gene expression, and other mechanisms. The "internal–external" and "external-internal" signal regulation are closely related and affect each other [44]. Studies have shown that Itgb1 is widely expressed in the nervous system and may play a central role in physiological processes [45]. Moreover, abnormal expression of Itgb1 is related to neuropathic pain, inflammation, and malignant diseases due to peripheral nerve injury [46]. Meanwhile, a recent study has pointed out that the administration of anti-Itgb1 antibody (β1-Ab) in the subacute SCI successfully prevented glial scar formation and enhanced axonal regeneration [47]. Ptprc can be used as a specific marker for the activation of microglia [48]. After T-cell activation, Ptprc recruits and dephosphorylates SKAP1 and FYN, and dephosphorylates LYN, thereby regulating LYN activity. Cd63 acts as a cell surface protein in the activation regulation of the signaling cascade of cell development, activation, growth, and movement. Cd63 plays a role in the activation of Itgb1 and integrin signaling, leading to the activation of AKT, FAK/PTK2, and MAP kinases [49]. A recent study has shown that the expression of Cd63 in both plasma and spinal cord tissues is increased after SCI in mice [50]. Lgals3 is involved in central nervous system inflammation, including neutrophil activation and adhesion, chemical traction of mononuclear macrophages, regulation of apoptotic neutrophils, and activation of mast cells [51]. Some scholars have pointed out that Lgals3 plays an important role in regulating microglial activation and neuroinflammation, serving as a biomarker of neurodegenerative diseases [52]. Vav1 can mediate the Rho activation of JAK as a guanine nucleotide exchange factor of Rho family GTP enzymes and plays an important role in the development and activation of T cells and B cells [53]. Vav1 is highly expressed in the early reaction and regeneration stage of sciatic nerve injury and activates the Rac1 GTP enzyme to promote axonal regeneration of DRG neurons [54]. Shc1, a downstream target of tumor suppressor p53, is essential for the ability of stress-activated p53 to induce increased intracellular oxidants, cytochrome c release, and apoptosis. Although studies have suggested that Shc1 plays an important role in regulating microglia polarization [55], research on Shc1 in the nervous system is relatively rare. Casp4 is involved in the activation of inflammatory bodies and has been widely studied. An Alzheimer's disease study pointed out that Casp4 regulates microglia [56] in a way that increases the pro-inflammatory process. Our study showed that these seven hub genes could be used as potential targets for exploring the molecular mechanisms related to the development of subacute SCI over time and for balancing the injury and repair induced by SCI.

In addition, we also verified the mRNA relative expression levels of five other hub genes (Fcgr2b, S100a4, Lamc1, Mapk12, and Vegfa) in spinal cord tissues in the subacute SCI, although their results do not conform to the expression matrix. Fcgr2b is a low-affinity receptor in the Fc region of immunoglobulin γ complex, and its cytoplasmic domain contains an inhibitory motif (ITIM), which is the only inhibitory Fcg receptor. The results of a mouse experimental model with cerebral ischemia showed that the expression level of Fcgr2b of microglia/macrophage activated by inflammatory reaction remained unchanged [57]. but the research on Fcgr2b in the nervous system was relatively rare, the reason why the relative mRNA expression level of Fcgr2b increased in the early stage of subacute SCI but decreased in the advanced stage needs to be further explored. S100a4 is a calcium-binding protein. It plays a role in a variety of cellular processes, including motility, angiogenesis, cell differentiation, apoptosis, and autophagy [58]. Lamc1 belongs to the extracellular matrix glycoprotein family and is the major non-collagen component of the basement membrane. They are related to various biological processes, including cell adhesion, differentiation, migration, signal transduction, neurite outgrowth, and metastasis [59]. In this study, there was no statistical difference between the 7-day after subacute SCI group and the sham operation group possibly due to sample size. Mapk12 plays a vital role in cellular processes, involving inflammation, cell growth, cell differentiation and cell cycle action [60]. Vegfa is able to inhibit apoptosis and induce endothelial cell proliferation, promote cell migration, and axonal regeneration [61]. RT-PCR verification showed that their expression levels decreased when stimulated by inflammatory injury, but increased when the disease progressed due to regeneration and repair of the injury site. Further experimental studies are needed to explain the mechanism of why the expression matrix is inconsistent with RT-PCR verification. This is also a limitation of our study. In our study, the predictive results were mainly based on bioinformatics analysis, and the research sample was from rats but not humans, so clinical and relevant biological experiments are needed to further explore and prove the action mechanism of central genes. To sum up, our study provides some meaningful insights into the disease progression mechanisms and targeted therapeutic strategies for subacute SCI patients, and helps to develop new therapeutic strategies for SCI.

In this study, comprehensive bioinformatics analysis was performed on the gene expression profiles of rats at three different time points after subacute SCI. A total of 12 hub genes were identified and seven hub genes with statistical significance in both the RT-PCR results and the expression matrix were identified. The biological functions and pathways of the identified genes provide more detailed molecular mechanisms for understanding the disease progression of subacute SCI. In conclusion, we have identified the hub genes and signaling pathways involved in subacute SCI by using bioinformatics methods, which may provide a molecular basis for the future treatment of SCI.

Availability of data and materials

All data generated or analyzed during this study are included in this published article and its supplementary information files.

References

  1. Jain NB, Ayers GD, Peterson EN, Harris MB, Morse L, O’Connor KC, Garshick E. Traumatic spinal cord injury in the United States, 1993–2012. JAMA. 2015;313(22):2236–43.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  2. James ND, McMahon SB, Field-Fote EC, Bradbury EJ. Neuromodulation in the restoration of function after spinal cord injury. Lancet Neurol. 2018;17(10):905–17.

    Article  PubMed  Google Scholar 

  3. O’Shea TM, Burda JE, Sofroniew MV. Cell biology of spinal cord injury and repair. J Clin Invest. 2017;127(9):3259–70.

    Article  PubMed  PubMed Central  Google Scholar 

  4. Vismara I, Papa S, Rossi F, Forloni G, Veglianese P. Current options for cell therapy in spinal cord injury. Trends Mol Med. 2017;23(9):831–49.

    Article  PubMed  CAS  Google Scholar 

  5. Calvert JS, Grahn PJ, Zhao KD, Lee KH. Emergence of epidural electrical stimulation to facilitate sensorimotor network functionality after spinal cord injury. Neuromodulation. 2019;22(3):244–52.

    Article  PubMed  Google Scholar 

  6. Lavis T, Goetz LL. Comprehensive care for persons with spinal cord injury. Phys Med Rehabil Clin N Am. 2019;30(1):55–72.

    Article  PubMed  Google Scholar 

  7. Thomaz SR, Cipriano G Jr, Formiga MF, Fachin-Martins E, Cipriano GFB, Martins WR, Cahalin LP. Effect of electrical stimulation on muscle atrophy and spasticity in patients with spinal cord injury - a systematic review with meta-analysis. Spinal Cord. 2019;57(4):258–66.

    Article  PubMed  Google Scholar 

  8. Stahel PF, VanderHeiden T, Finn MA. Management strategies for acute spinal cord injury: current options and future perspectives. Curr Opin Crit Care. 2012;18(6):651–60.

    Article  PubMed  Google Scholar 

  9. Alizadeh A, Dyck SM, Karimi-Abdolrezaee S. Traumatic spinal cord injury: an overview of pathophysiology, models and acute injury mechanisms. Front Neurol. 2019;10:282.

    Article  PubMed  PubMed Central  Google Scholar 

  10. Tran AP, Warren PM, Silver J. The biology of regeneration failure and success after spinal cord injury. Physiol Rev. 2018;98(2):881–917.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  11. Beattie MS, Li Q, Bresnahan JC. Cell death and plasticity after experimental spinal cord injury. Prog Brain Res. 2000;128:9–21.

    Article  PubMed  CAS  Google Scholar 

  12. Blight AR. Spinal cord injury models: neurophysiology. J Neurotrauma. 1992;9(2):147–9.

    Article  PubMed  CAS  Google Scholar 

  13. Grossman SD, Rosenberg LJ, Wrathall JR. Relationship of altered glutamate receptor subunit mRNA expression to acute cell loss after spinal cord contusion. Exp Neurol. 2001;168(2):283–9.

    Article  PubMed  CAS  Google Scholar 

  14. Silva NA, Sousa N, Reis RL, Salgado AJ. From basics to clinical: a comprehensive review on spinal cord injury. Prog Neurobiol. 2014;114:25–57.

    Article  PubMed  Google Scholar 

  15. Sykova E, Homola A, Mazanec R, Lachmann H, Konradova SL, Kobylka P, Padr R, Neuwirth J, Komrska V, Vavra V, et al. Autologous bone marrow transplantation in patients with subacute and chronic spinal cord injury. Cell Transplant. 2006;15(8–9):675–87.

    Article  PubMed  Google Scholar 

  16. Gensel JC, Zhang B. Macrophage activation and its role in repair and pathology after spinal cord injury. Brain Res. 2015;1619:1–11.

    Article  PubMed  CAS  Google Scholar 

  17. Jones TB, McDaniel EE, Popovich PG. Inflammatory-mediated injury and repair in the traumatically injured spinal cord. Curr Pharm Des. 2005;11(10):1223–36.

    Article  PubMed  CAS  Google Scholar 

  18. Rust R, Kaiser J. Insights into the dual role of inflammation after spinal cord injury. J Neurosci. 2017;37(18):4658–60.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  19. Ortuno FM, Torres C, Glosekotter P, Rojas I. New trends in biomedical engineering and bioinformatics applied to biomedicine - special issue of IWBBIO 2014. Biomed Eng Online. 2015;14(Suppl 2):I1.

    Article  PubMed  PubMed Central  Google Scholar 

  20. Duran RC, Yan H, Zheng Y, Huang X, Grill R, Kim DH, Cao Q, Wu JQ. The systematic analysis of coding and long non-coding RNAs in the sub-chronic and chronic stages of spinal cord injury. Sci Rep. 2017;7:41008.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  21. Du H, Shi J, Wang M, An S, Guo X, Wang Z. Analyses of gene expression profiles in the rat dorsal horn of the spinal cord using RNA sequencing in chronic constriction injury rats. J Neuroinflammation. 2018;15(1):280.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  22. Guo L, Lv J, Huang YF, Hao DJ, Liu JJ. Bioinformatics analyses of differentially expressed genes associated with spinal cord injury: a microarray-based analysis in a mouse model. Neural Regen Res. 2019;14(7):1262–70.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  23. Liu ZG, Li Y, Jiao JH, Long H, Xin ZY, Yang XY. MicroRNA regulatory pattern in spinal cord ischemia-reperfusion injury. Neural Regen Res. 2020;15(11):2123–30.

    Article  PubMed  PubMed Central  Google Scholar 

  24. Niu SP, Zhang YJ, Han N, Yin XF, Zhang DY, Kou YH. Identification of four differentially expressed genes associated with acute and chronic spinal cord injury based on bioinformatics data. Neural Regen Res. 2021;16(5):865–70.

    Article  PubMed  Google Scholar 

  25. Fang S, Zhong L, Wang AQ, Zhang H, Yin ZS. Identification of regeneration and hub genes and pathways at different time points after spinal cord injury. Mol Neurobiol. 2021;58(6):2643–62.

    Article  PubMed  CAS  Google Scholar 

  26. Siebert JR, Middelton FA, Stelzner DJ. Intrinsic response of thoracic propriospinal neurons to axotomy. BMC Neurosci. 2010;11:69.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  27. Chamankhah M, Eftekharpour E, Karimi-Abdolrezaee S, Boutros PC, San-Marina S, Fehlings MG. Genome-wide gene expression profiling of stress response in a spinal cord clip compression injury model. BMC Genomics. 2013;14:583.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  28. Baligand C, Chen YW, Ye F, Pandey SN, Lai SH, Liu M, Vandenborne K. Transcriptional pathways associated with skeletal muscle changes after spinal cord injury and treadmill locomotor training. Biomed Res Int. 2015;2015:387090.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  29. Leek JT, Johnson WE, Parker HS, Jaffe AE, Storey JD. The sva package for removing batch effects and other unwanted variation in high-throughput experiments. Bioinformatics. 2012;28(6):882–3.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

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

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  31. Kumar L. M EF: Mfuzz: a software package for soft clustering of microarray data. Bioinformation. 2007;2(1):5–7.

    Article  PubMed  PubMed Central  Google Scholar 

  32. Zhou Y, Zhou B, Pache L, Chang M, Khodabakhshi AH, Tanaseichuk O, Benner C, Chanda SK. Metascape provides a biologist-oriented resource for the analysis of systems-level datasets. Nat Commun. 2019;10(1):1523.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  33. Jeong H, Mason SP, Barabasi AL, Oltvai ZN. Lethality and centrality in protein networks. Nature. 2001;411(6833):41–2.

    Article  PubMed  CAS  Google Scholar 

  34. Schmidt E, Raposo P, Vavrek R, Fouad K. Inducing inflammation following subacute spinal cord injury in female rats: a double-edged sword to promote motor recovery. Brain Behav Immun. 2021;93:55–65.

    Article  PubMed  CAS  Google Scholar 

  35. Li D, Dossa K, Zhang Y, Wei X, Wang L, Zhang Y, Liu A, Zhou R, Zhang X. GWAS uncovers differential genetic bases for drought and salt tolerances in sesame at the germination stage. Genes. 2018;9(2):87.

    Article  PubMed Central  CAS  Google Scholar 

  36. Shi K, Bing ZT, Cao GQ, Guo L, Cao YN, Jiang HO, Zhang MX. Identify the signature genes for diagnose of uveal melanoma by weight gene co-expression network analysis. Int J Ophthalmol. 2015;8(2):269–74.

    PubMed  PubMed Central  Google Scholar 

  37. Pineau I, Lacroix S. Proinflammatory cytokine synthesis in the injured mouse spinal cord: multiphasic expression pattern and identification of the cell types involved. J Comp Neurol. 2007;500(2):267–85.

    Article  PubMed  CAS  Google Scholar 

  38. de Rivero Vaccari JP, Lotocki G, Marcillo AE, Dietrich WD, Keane RW. A molecular platform in neurons regulates inflammation after spinal cord injury. J Neurosci. 2008;28(13):3404–14.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  39. Dusaban SS, Purcell NH, Rockenstein E, Masliah E, Cho MK, Smrcka AV, Brown JH. Phospholipase C epsilon links G protein-coupled receptor activation to inflammatory astrocytic responses. Proc Natl Acad Sci U S A. 2013;110(9):3609–14.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  40. Trias E, Kovacs M, King PH, Si Y, Kwon Y, Varela V, Ibarburu S, Moura IC, Hermine O, Beckman JS, et al. Schwann cells orchestrate peripheral nerve inflammation through the expression of CSF1, IL-34, and SCF in amyotrophic lateral sclerosis. Glia. 2020;68(6):1165–81.

    Article  PubMed  Google Scholar 

  41. Brockie S, Hong J, Fehlings MG. The role of microglia in modulating neuroinflammation after spinal cord injury. Int J Mol Sci. 2021;22(18):9706.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  42. Devanney NA, Stewart AN, Gensel JC. Microglia and macrophage metabolism in CNS injury and disease: The role of immunometabolism in neurodegeneration and neurotrauma. Exp Neurol. 2020;329:113310.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  43. Mesquida-Veny F, Del Rio JA, Hervera A. Macrophagic and microglial complexity after neuronal injury. Prog Neurobiol. 2021;200:101970.

    Article  PubMed  CAS  Google Scholar 

  44. Askari JA, Tynan CJ, Webb SE, Martin-Fernandez ML, Ballestrem C, Humphries MJ. Focal adhesions are sites of integrin extension. J Cell Biol. 2010;188(6):891–903.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  45. Barros CS, Nguyen T, Spencer KS, Nishiyama A, Colognato H, Muller U. Beta1 integrins are required for normal CNS myelination and promote AKT-dependent myelin outgrowth. Development. 2009;136(16):2717–24.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  46. Previtali SC, Feltri ML, Archelos JJ, Quattrini A, Wrabetz L, Hartung H. Role of integrins in the peripheral nervous system. Prog Neurobiol. 2001;64(1):35–49.

    Article  PubMed  CAS  Google Scholar 

  47. Hara M, Kobayakawa K, Ohkawa Y, Kumamaru H, Yokota K, Saito T, Kijima K, Yoshizaki S, Harimaya K, Nakashima Y, et al. Interaction of reactive astrocytes with type I collagen induces astrocytic scar formation through the integrin-N-cadherin pathway after spinal cord injury. Nat Med. 2017;23(7):818–28.

    Article  PubMed  CAS  Google Scholar 

  48. Li Y, Zhou D, Ren Y, Zhang Z, Guo X, Ma M, Xue Z, Lv J, Liu H, Xi Q, et al. Mir223 restrains autophagy and promotes CNS inflammation by targeting ATG16L1. Autophagy. 2019;15(3):478–92.

    Article  PubMed  CAS  Google Scholar 

  49. Estebanez B, Jimenez-Pavon D, Huang CJ, Cuevas MJ, Gonzalez-Gallego J. Effects of exercise on exosome release and cargo in in vivo and ex vivo models: a systematic review. J Cell Physiol. 2021;236(5):3336–53.

    Article  PubMed  CAS  Google Scholar 

  50. Khan NZ, Cao T, He J, Ritzel RM, Li Y, Henry RJ, Colson C, Stoica BA, Faden AI, Wu J. Spinal cord injury alters microRNA and CD81+ exosome levels in plasma extracellular nanoparticles with neuroinflammatory potential. Brain Behav Immun. 2021;92:165–83.

    Article  PubMed  CAS  Google Scholar 

  51. Soares LC, Al-Dalahmah O, Hillis J, Young CC, Asbed I, Sakaguchi M, O’Neill E, Szele FG. Novel galectin-3 roles in neurogenesis inflammation and neurological diseases. Cells. 2021;10(11):3047.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  52. Tan Y, Zheng Y, Xu D, Sun Z, Yang H, Yin Q. Galectin-3: a key player in microglia-mediated neuroinflammation and Alzheimer’s disease. Cell Biosci. 2021;11(1):78.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  53. Montresor A, Bolomini-Vittori M, Toffali L, Rossi B, Constantin G, Laudanna C. JAK tyrosine kinases promote hierarchical activation of Rho and Rap modules of integrin activation. J Cell Biol. 2013;203(6):1003–19.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  54. Wang Q, Gong L, Mao S, Yao C, Liu M, Wang Y, Yang J, Yu B, Chen G, Gu X. Klf2-Vav1-Rac1 axis promotes axon regeneration after peripheral nerve injury. Exp Neurol. 2021;343:113788.

    Article  PubMed  CAS  Google Scholar 

  55. Peng J, Pang J, Huang L, Enkhjargal B, Zhang T, Mo J, Wu P, Xu W, Zuo Y, Peng J, et al. LRP1 activation attenuates white matter injury by modulating microglial polarization through Shc1/PI3K/Akt pathway after subarachnoid hemorrhage in rats. Redox Biol. 2019;21:101121.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  56. Kajiwara Y, McKenzie A, Dorr N, Gama Sosa MA, Elder G, Schmeidler J, Dickstein DL, Bozdagi O, Zhang B, Buxbaum JD. The human-specific CASP4 gene product contributes to Alzheimer-related synaptic and behavioural deficits. Hum Mol Genet. 2016;25(19):4315–27.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  57. Mota M, Porrini V, Parrella E, Benarese M, Bellucci A, Rhein S, Schwaninger M, Pizzi M. Neuroprotective epi-drugs quench the inflammatory response and microglial/macrophage activation in a mouse model of permanent brain ischemia. J Neuroinflammation. 2020;17(1):361.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  58. Li ZH, Bresnick AR. The S100A4 metastasis factor regulates cellular motility via a direct interaction with myosin-IIA. Cancer Res. 2006;66(10):5173–80.

    Article  PubMed  CAS  Google Scholar 

  59. Ustun Y, Reibetanz M, Brachvogel B, Nischt R, Eckes B, Zigrino P, Krieg T. Dual role of laminin511 in regulating melanocyte migration and differentiation. Matrix Biol. 2019;80:59–71.

    Article  PubMed  CAS  Google Scholar 

  60. Xu W, Liu R, Dai Y, Hong S, Dong H, Wang H. The role of p38gamma in cancer: from review to outlook. Int J Biol Sci. 2021;17(14):4036–46.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  61. Cattin AL, Burden JJ, Van Emmenis L, Mackenzie FE, Hoving JJ, Garcia Calavia N, Guo Y, McLaughlin M, Rosenberg LH, Quereda V, et al. Macrophage-induced blood vessels guide schwann cell-mediated regeneration of peripheral nerves. Cell. 2015;162(5):1127–39.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

Download references

Acknowledgements

We thank all the authors for their time and effort in this study, and all the institutions for their financial support.

Funding

This work was supported by the Science and Technology Bureau of Nantong Municipality under Grant No. JC2020013 and the Health Committee of Jiangsu Province under Grant No. ZDB2020004. The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.

Author information

Authors and Affiliations

  1. The Second Affiliated Hospital of Nantong University, No.6, North Road, 226000, Haierxiang, Nantong, Jiangsu, People’s Republic of China

    Lei Yan, Jiawei Fu, Xiong Dong, Baishen Chen, Hongxiang Hong & Zhiming Cui

Authors
  1. Lei Yan
    View author publications

    You can also search for this author in PubMed Google Scholar

  2. Jiawei Fu
    View author publications

    You can also search for this author in PubMed Google Scholar

  3. Xiong Dong
    View author publications

    You can also search for this author in PubMed Google Scholar

  4. Baishen Chen
    View author publications

    You can also search for this author in PubMed Google Scholar

  5. Hongxiang Hong
    View author publications

    You can also search for this author in PubMed Google Scholar

  6. Zhiming Cui
    View author publications

    You can also search for this author in PubMed Google Scholar

Contributions

LY conceived and designed the experiments, analyzed the data, prepared figures and/or tables. Data collection were performed by JF, XD, BC and HH. ZC contributed to the study conception and design. The first draft of the manuscript was written by LY and all authors commented on previous versions of the manuscript. All authors read and approved the final manuscript.

Corresponding author

Correspondence to Zhiming Cui.

Ethics declarations

Ethics approval and consent to participate

We declare that this study was conducted in accordance with ARRIVE guidelines, that all methods were conducted in accordance with relevant guidelines and regulations. All study procedures were approved by the institutional animal care and utilization committee of Laboratory Animal Center of Nantong University (Protocol number: S20211229-408).

Consent for publication

Not applicable.

Competing interests

The authors declare that they have 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

Yan, L., Fu, J., Dong, X. et al. Identification of hub genes in the subacute spinal cord injury in rats. BMC Neurosci 23, 51 (2022). https://doi.org/10.1186/s12868-022-00737-5

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doi.org/10.1186/s12868-022-00737-5

Keywords

  • Spinal cord injury
  • WGCNA
  • PPI
  • Bioinformatics analysis

BMC Neuroscience

ISSN: 1471-2202

Contact us