Integrated bioinformatics analysis identifies the effects of Sema3A/NRP1 signaling in oligodendrocytes after spinal cord injury in rats

Objective To investigate the effect of Sema3A/NRP1 signaling in oligodendrocytes (OLs) after spinal cord injury. Methods Three analysis strategies, namely differential expression gene analysis, Gene Ontology (GO) enrichment, and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analysis, were applied. The protein-protein interaction (PPI) network was constructed using the STRING website to explore the correlation between Sema3A/NRP1 and oligodendrocytes. Then, the T10 spinal cord segment of rats was injured by the Allen method to establish a spinal cord injury (SCI) model. Real-time quantitative PCR, Western blotting, Nissl staining and immunofluorescence staining were used to detect the effect of Sema3A/NRP1 signaling on oligodendrocytes in vivo. Results After the SCI model was established, significantly fewer oligodendrocytes were observed. At the same time, R software was used to analyze the expression of related genes, and NRP1 expression was increased. PCR also demonstrated similar results, and NRP1 ligand Sema3A was also upregulated. KEGG and GO functional enrichment analysis indicated that the SCI model was mainly related to cytokine interaction, cell proliferation, differentiation and maturation. Interestingly, we found that NRP1 was involved in semaphorin-plexin signaling pathway neuronal projection guidance and axon guidance, mediating cell growth and migration. Moreover, Sema3A/NRP1 signaling was closely associated with platelet-derived growth factor receptor α (PDGFRα) in the PPI network. When Sema3A/NRP1 signaling was specifically blocked at early stages, PDGFRα expression was effectively inhibited, and the expression of OLs was promoted. Furthermore, inhibition of Sema3A/NRP1 signaling increased the Basso-Beattie-Bresnahan (BBB) score of lower limb motor function in SCI rats and promoted the survival of motor neurons in the ventral horn of the injured spinal cord. Conclusion Our data suggest that Sema3A/NRP1 signaling may regulate the development of OPCs and OLs after SCI, thereby affecting functional recovery.


INTRODUCTION
Spinal cord injury (SCI) has become a major public health concern with no effective treatment. Clinically, SCI is mainly manifested by irreversible neurological impairment, resulting in a series of dysfunctions (Almeida et al., 2018). Current evidence suggests that there are approximately 133,000 to 226,000 new cases of SCI per year globally. In the United States, the hospital charges related to SCI are close to US$9.7 million per year, which burdens the medical care system, the patient's family and society (Wyndaele & Wyndaele, 2006;Kumar et al., 2018). To date, little is currently known about the underlying mechanisms, warranting further studies.
It is widely acknowledged that oligodendrocyte precursor cells (OPCs) and oligodendrocytes (OLs) are widely distributed in the central nervous system (CNS). OLs are the glial cells of myelin axons derived from OPCs, playing an important role in the integrity of axons, maintaining the normal structure of axons and optimizing axon conduction (Gritsch et al., 2014;Welsh & Kuenas, 2018). SCI and its secondary impairment may lead to loss of OLs and neurological dysfunction (Baumann & Pham-Dinh, 2001). Therefore, the development and recovery of OL are essential during the recovery process after nerve damage. Neuropilin1 (NRP1) is a transmembrane glycoprotein mainly expressed in dendrites and axons of neurons that plays a key role in the functional recovery of SCI after extracellular binding with Semaphorin3A (Sema3A) (Parker et al., 2012). In recent years (Shim et al., 2013;Nakanishi, Fujita & Yamashita, 2019), the role of Sema3A/NRP1 signaling in directional migration of OLs and axonal regeneration after nerve injury has been reported. In addition, Sema3A/NRP1 signaling is a directional guide for oligodendrocyte migration in the hallmark of multiple sclerosis (Jiang et al., 2010). However, the role of Sema3A/NRP1 signaling in regulating OLs development after SCI remains largely unknown. This study observed the effects of Sema3A/NRP1 signaling in OLs associated with SCI in rats to explore the molecular mechanism underlying functional recovery.

Data download and processing
The gene expression data for the present study was downloaded from Gene Expression Omnibus (GEO) database (https://www.ncbi.nlm.nih.gov/geo/) using the accession number (Chamankhah et al., 2013). The dataset GSE450006 was based on the GPL1355 platform and contained 4 spinal cord samples from sham rats and 20 from SCI model female rats. 24 samples were selected for further data processing. Specifically, the raw expression datasets were log2-transformed for further analysis. Subsequently, according to platform annotation files, the probes were matched with corresponding genes, and non-matched probes to gene symbol were removed. When multiple probe names were assigned to a single gene symbol, the values were averaged. Finally, the data was assembled for subsequent analysis.

Identification of DEGs
Data were collated as previously described in Ren et al. (2020). Specifically, R software and related R packages were used to normalize and analyze differentially expressed genes (DEGs). Log2 fold change (FC) greater than 2 and adjusted p-value < 0.05 were selected as the cutoff criterion for DEGs screening. The analysis results were presented by heat maps and volcano plots drawn in R Studio software (version:1.2.1335;RStudio Team, 2015).

Functional enrichment analysis
To further analyze the biological processes of DEGs in the sham group compared with the SCI group, Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analyses were carried out using the cluster Profiler in R package and R studio software (Yu et al., 2012). Significantly enriched GO terms and KEGG pathways were screened using the criteria: p-value <0.05 and an enriched gene count >2.

Protein-protein interaction (PPI) network analysis
To explore the mutual relationship between proteins encoded by different genes, DEGs were imported into the STRING website (version:11.0) for further analysis (Szklarczyk et al., 2019;Snel et al., 2000;Athanasios et al., 2017). The lowest interaction score was greater than 0.4, and isolated nodes in the network were removed. Then the analysis results were output to the TSV format file, and the detailed processing and module analysis were carried out by Cytoscape software (version:3.7.1). Molecular Complex Detection (MCODE) (Fischer et al., 2015) is a plug-in downloaded from Cytoscape App Store used to find closely connected nodes in a complex network based on topology. In the present study, we applied plug-in MCODE to select critical modules of PPI network.

Animals
All healthy male SD rats (eight-week old, 220 ± 10 g) were provided by Shanghai Xipu Bikai Shanghai Xipu Laboratory Animal Co. (animal license No: SCXk (Shanghai) 2018-0006) and were housed in the Experimental Animal Center of Zhejiang Chinese Medical University (AAALAC, animal license No: SYXK (Zhejiang) 2018-0012). The animals were housed in environmentally-controlled standard conditions (illumination, 12/12 h light/dark cycle; humidity, 50-60%; room temperature, 21-23 • C) with ad libitum access to food and water. 3% pentobarbital sodium (0.15 ml/kg, i.p.) was intraperitoneally injected before the end of the experiment. Rats were divided into a sham group (sham, only underwent laminectomy does not damage the spinal cord), the SCI model group (SCI, underwent laminectomy and damage the T10 spinal cord, details in 'SCI Rat Model Establishment'), the SCI+AAV Sema3A group (spinal cord stereotactic injection of AAV virus, details in 2.4), and the SCI+AAV NC group (spinal cord stereotactic injection of AAV control virus, details in 'Virus Injection'). All experiments were performed in compliance with relevant ethical regulations and approved by the animal ethics committee of Zhejiang Chinese Medical University (IACUC-20190128-02). All the experimental protocols strictly followed the National Institutes of Health guidelines for the care and use of laboratory animals (NIH Publication No. 8023).

SCI rat model establishment
To establish a contusive SCI model at T 10 , contusive injury was induced on the T 10 spinal cord using the MASCIS weight-drop device with a 5 × 10 g/cm gravitational potential energy after T 10 laminectomy. The severity and consistency of injury were verified by checking the bruise on the spinal cord or tail-flick of the rats after weight drop. All animals received intraperitoneal injections of penicillin (100 U/d) for three days to prevent infection.

Nissl staining
Procedures for Nissl staining were as follows. Briefly, after dehydration in 30% sugar solution, the tissues were sliced with a freezing microtome. The slices were incubated with nissl staining solution in a 37 • C constant temperature water bath for 10 min, followed by alcohol gradient decolorization for 1 min each, and washed with distilled water for 30 s immediately. Finally, the slices were sealed with neutral resin. Images were taken with a light microscope.

Real-time quantitative PCR (qPCR) analysis
First, total RNA was extracted from the T10 spinal cord tissue and the quality of RNA tested (Invitrogen, Carlsbad, CA, USA), cDNA was synthesized using PrimeScript RT reagent Kit (TaKaRa). qRT-PCR was conducted by using the Fast Start Universal SYBR Green Master kit (TaKaRa Bio Inc, Beijing, China) with a 20 µl reaction system. The data for qRT-PCR were collected with CFX96 Real-Time System (BioRad, USA). All qRT-PCR results were analyzed using comparative Ct methods (2 − Ct). Primer sequences are listed in Table 1.

Western blotting
Rat spinal cord tissue was homogenized with radioimmunoprecipitation assay (RIPA) buffer. Protein concentrations were determined using the bicinchoninic acid (BCA) assay. Furthermore, equal quantities of protein were separated via 8% SDS-PAGE gel and transferred by electroblotting polyvinyl difluoride (PVDF) membranes. The membranes were blocked with 5% nonfat milk in TBST containing 0.1% Tween 20 at room temperature for 1.5 h and then incubated overnight at 4 • C with primary antibody: semaphorin3A (GTX37550, 1:500), neuropilin1 (ab81321, 1:1,000), Plateletderived growth factor receptor α (PDGFR α) (ab203491, 1:1,000), and β-actin (ab20271, 1:5,000). The following day, the membrane was incubated at room temperature for 2 h with secondary antibodies (CST7074s, 1:2000). Finally, immunoreactive bands were detected by enhanced chemiluminescence and visualized with an Image Quant LAS 4000. For statistical comparisons, primary antibody values were normalized against β-actin for each run in the normal Western blots. The density of each band was measured by ImageJ analysis software.

Behavioral testing
Hindlimb locomotor function was evaluated using Basso-Beattie-Bresnahan (BBB) score. As previously described (Wei et al., 2017), the BBB scale ranges from 0-21 and reflects the extent of hind limb movements , weight support, stepping ability, coordination, and trunk stability. Evaluation and analysis were conducted by a specialist with an independent assistant blinded to the group allocation. Then the motor performance of the rats in an open field was scored, including hind limb joint movement, weight support, plantar stepping, coordination, paw position, and trunk and tail control. Scoring was performed as illustrated in Table 2.

Statistical analyses
The statistical analyses between control and experimental groups were conducted by an unpaired t -test or one-way ANOVA followed by Tukey Kramer tests in GraphPad Prism7. Results are expressed as mean ± SEM. A p-value < 0.05 was statistically significant. Continuous paw and palm bearing movement and consistent FL-HL coordination, no or occasional ground grasping movement in the forward motion of the forelimbs, and the position of the main claw parallel to the body at the initial contact 16 In the gait, the continuous paw landing and the coordinated movement of the front and rear limbs are common in the process of grasping the ground; The main claw position is parallel to the body at initial contact, and rotates after load transfer 17 In the gait, the continuous paw landing and the coordinated movement of the front and rear limbs are common in the process of grasping the ground; The main claw position is parallel to the body at initial contact and load transfer 18 In the gait, the continuous paw touches the ground in a coordinated manner with the front and rear limbs. In the process of progress, the continuous paw grasps the ground. The position of the main paw is parallel to the body at the initial contact 19 In the gait, the continuous paw touches the ground in a coordinated manner with the front and rear limbs. The continuous paw grasps the ground in the process of advancing. The position of the main paw is parallel to the body at the initial contact and load transfer 20 The position of the main claw is parallel to the body during initial contact and weight transfer. The trunk is unstable and the tail kept cocking up 21 The position of the main claw is parallel to the body at the initial contact and load transfer, and the trunk is stable and the tail kept cocking up

The expression of OLs and PDGFR α after SCI
The injured spinal cord tissues were harvested on days 1, 7, and 14 after SCI, and the expression of OLs was detected by qPCR. The results showed that compared with the sham group, olig2, the early specific marker of OLs was significantly reduced after SCI in the injured spinal cord tissue of rats (Fig. 1A). Sox10, the mature OLs specific marker, exhibited low levels after SCI (Fig. 1B). Furthermore, the specific marker of OPCs, Platelet-derived growth factor receptor α (PDGFR α), was significantly increased after SCI (Fig. 1C).
We further detected the expression of PDGFR α in spinal cord tissue of SCI by western blotting, and the results showed that the expression was increased at 7 and 14 days ( 1D). Consistent with previous studies (Assinck et al., 2017) that reported a large number of OL deaths after SCI, and the spontaneous differentiation of OPCs could be stimulated.

Increased expression of Sema3A and NRP1 after SCI
The dataset GSE450006 based on the GPL1355 platform contained four spinal cord samples from sham rats and 20 from SCI model rats. A total of 299 DEGs were obtained, including 176 upregulated and 123 downregulated genes ( Fig. 2A), and the volcano plot  showed upregulated and downregulated gene expression in each sample (Fig. 2B). We observed that OLs were significantly decreased ( Fig. 1) while NRP1 and Sema3A were upregulated (Figs. 3A and 3B) in SCI rats. Based on the fact that NRP1 can only exert its signal transduction function after forming a complex with Sema3A extracellularly (Parker et al., 2012;Hu et al., 2020), we performed immunofluorescence staining on the spinal cord tissues of sham and SCI seven-day rats. We observed that Sema3A and NRP1 were co-expressed in the SCI model (Fig. 3C).

Analysis of the function and pathway of DEGs in SCI model rats
To further explore the mechanism of SCI, we performed GO and KEGG pathway enrichment analyses of DEGs with the cluster profiler package in R software. GO is one of the most widely used ontologies, which annotates genes in terms of their molecular functions, cellular constituent, and biological processes (Feunang et al., 2016). Moreover, the KEGG database is an important tool for identifying functional and metabolic pathways (Sun et al., 2021). In this study, we used the R package ''clusterProfiler'' that automates the process of biological-term classification and the enrichment analysis of gene clusters. The screening criteria for statistically significant GO terms or pathways was a p-value less than 0.05. We found that the biological functions of DEGs were mainly related to wound healing, cell growth, chemotaxis and metabolism regulation, extracellulaar matrix organization in the SCI group (Fig. 4A). KEGG pathway analysis of DEGs between the SCI and the sham groups revealed significant enrichment in the interaction between cells and cytokines, cell proliferation, differentiation and maturation, ECM-receptor interaction (Fig. 4B). We further condensed all enriched pathways processed by KEGG and GO functional enrichment analysis, summarized all NRP1 involved in the results (Table S1), among which unexpected findings include the semaphorin-plexin signaling pathway involved in neuron projection guidance and axon guidance, mediating cell growth and migration (Table 3).

Sema3A regulates the expression of NRP1 after SCI
To further define the relationship between Sema3A and NRP1, we designed SNCA-shRNA sequence plasmids with rad fluorescent tags and integrated them into an adeno-associated virus (AAV). Subsequently, we injected AAV2/9-U6-shRNA (Sema3A)-CAG-tdtomato into the junction between T 9 and T 10 spinal cord rats 21 days before modeling to knock down the expression of Sema3A. Four weeks after injection, we took the thoracic spinal cord for sectioning and observed it under the microscope. It was found that the T10 segment had the best effect of AAV virus transfection (Fig. S1). Next, we observed the changes in NRP1 and Sema3A after SCI. The immunofluorescence analysis results showed that compared with the sham group, the expression of Sema3A and NRP1 in spinal cord tissue of SCI rats increased, and knockdown of Sema3A could significantly downregulate the expression of NRP1 (Figs. 5A-5C). The qPCR results were consistent with the results obtained by immunofluorescence analysis results(Figs. 5E-5F).

Sema3A/NRP1 signaling regulates the expression of PDGFR α and OLs
To investigate the role of Sema3A/NRP1 signaling in SCI, we constructed the PPI network using STRING database and visualized it with Cytoscape. The upregulated genes were highlighted in red, downregulated genes in green, key genes in blue, and gray lines represented the protein interaction relationships. The results imply a possible regulatory relationship between the Sema3A/NRP1 and PDGFR α (Fig. S2). We measured PDGFR α mRNA and PDGFR α protein levels by real-time PCR, immunofluorescence, and western blotting. The results revealed that the overexpression of PDGFR α in SCI model rats was significantly reduced by inhibition of Sema3A/NRP1 signaling (Figs. 5A, 5D, 6C, 6D, and 6G). To better understand the cellular mechanism of Sema3A/NRP1signaling in mediating motor function recovery in SCI, we investigated whether Sema3A/NRP1 signaling was involved in the maturation of OLs after SCI, which is also a key link in neurological function recovery (Tetzlaff et al., 2011). The qPCR results indicated that the expression of Olig2 and Sox10 of SCI group rats (after Sema3A/NRP1 signaling inhibition) significantly increased on day seven (Figs. 6A-6B) but showed no significant change on day 14 (Figs. 6E-6F). The above data indicate that Sema3A/NRP1 signaling may regulate the expression of PDGFR α and OLs.

Inhibition of Sema3A/NRP1 signaling promotes function recovery after SCI
The recovery of hind limb function was assessed using the BBB rating scale to examine the effect of Sema3A/NRP1 signaling on locomotor recovery after SCI. The line graph shows the BBB score for each group within 1 week after SCI. We found that all rats had a BBB score of 21 before SCI and lost locomotor function after SCI. When Sema3A/NRP1 signaling was inhibited, the BBB score was higher than the SCI group after 7days (Fig. 7C).
Nissl staining showed that the color intensity and the number of nissl bodies in the SCI-7d group and the SCI-7d+AAV NC decreased significantly compared with the sham group. Compared with the SCI group, the SCI-7d+AAV Sema3A group showed deeper nissl bodies staining with significantly increased expression of tiger-spotted motor neurons that were smaller than in the sham group (Figs. 7A-7B). These results suggest that suppressing Sema3A/NRP1 signal transduction is beneficial to functional recovery after injury.

DISCUSSION
In this study, we revealed the related mechanism of Sema3A/NRP1 signaling on OPCs and OLs through bioinformatics analysis. Previous studies (Mctigue, Wei & Stokes, 2001) have shown that OLs are important cells involved in nerve repair after SCI, mainly formed by the proliferation, migration and differentiation of OPCs. Massive OLs death after SCI is accompanied by demyelination which restricts axon electrical signal conduction and leads to dysfunction post-SCI (Totoiu & Keirstead, 2005). Interestingly, cell transplantation can enhance the differentiation of OPCs and OL post-SCI, thus improving myelination and promoting spinal cord nerve function recovery. Our study found that the expression of OL markers olig2 and sox10 decreased significantly post-SCI, indicating that significant OL death post-SCI may seriously hinder the repair of nerve function after injury. NRP1 was first documented as a key receptor that inhibits the extension of axons when bound to Sema3A (Goshima et al., 2000). Sema3A/NRP1 signaling is one of the key molecular signals for the formation and growth of myelin sheath (Hamashima et al., 2020). Orr, Fetter & Davis (2017) extracted SM-216289, a potent and selective inhibitor of Sema3A, and injected it into the lesion site of rats with spinal cord transection for 4 weeks. The results showed that the rats treated with SM-216289 exhibited axon regeneration, myelination, and angiogenesis, and damaged nerve function was significantly restored. In addition, inhibition of Sema3A/NRP1 signaling can promote motor function recovery after SCI (Zhang et al., 2014). Our previous studies (Hu et al., 2020) confirmed that Sema3A /NRP1 signaling plays an important role in functional recovery after SCI in rats. Herein, we found that the expression of Sema3A and NRP1 was significantly increased after SCI. However, the mechanism of action between Sema3A/NRP1 signal and OLs after SCI is still unclear.
GO and KEGG analyses showed that the most significant upregulated genes were related to wound healing, cell-cytokine interaction, cell proliferation, differentiation and maturation in the SCI group. Meanwhile, NRP1 is well-established to be involved in the semaphorin-plexin signaling pathway neuronal projection guidance and axon guidance, mediating cell growth and migration. Of note, PPI network analysis found that Sema3A/NRP1 signaling may regulate the expression of PDGFR α. It is well known that PDGFR α, one of the specific markers of OPCs, is produced by both astrocytes and neurons (Hamashima et al., 2020;Duncan et al., 2020). OPCs are tightly correlated with the maturation of OLs and can differentiate into mature OLs . Hence, we speculated that Sema3A/NRP1 plays an important role in post-injury repair and may be involved in the development of OPCs and OL after SCI.
Evidence has verified that Sema3A is a selective and reversible inhibitor of OPCs differentiation in vitro (Syed et al., 2011). Drug inhibition or knockout mice blocking  Sema3A/NRP1 signal transduction can increase OPCs recruitment, reduce demyelination lesions (Wang & Li, 2020) and axon growth-inhibitory factor (Boyd, Zhang & Williams, 2013) expression. To clarify the relationship between 3A/NRP1 and OPCs and OLs, we selectively downregulated the expression of NRP1 ligand Sema3A in the spinal cord at the injury site and observed that PDGFR α expression was downregulated when the expression of Sema3A/NRP1 signaling was inhibited, while olig2 and sox10 were significantly increased on day 7. These results suggest that inhibition of Sema3A/NRP1 signaling may promote the growth and development of OLs after injury. Moreover, we found that motor function improved with downregulation of Sema3A/NRP1 signaling with increased motor neuron survival in the ventral horn of the spinal cord. According to the above findings, decreasing Sema3A/NRP1 signaling can promote the repair of spinal motor neuron function and lower limb motor function post-SCI.