Identification of lncRNA expression profile in the spinal cord of mice following spinal nerve ligation-induced neuropathic pain

Background Neuropathic pain that caused by lesion or dysfunction of the nervous system is associated with gene expression changes in the sensory pathway. Long noncoding RNAs (lncRNAs) have been reported to be able to regulate gene expression. Identifying lncRNA expression patterns in the spinal cord under normal and neuropathic pain conditions is essential for understanding the genetic mechanisms behind the pathogenesis of neuropathic pain. Results Spinal nerve ligation (SNL) induced rapid and persistent pain hypersensitivity, characterized by mechanical allodynia and heat hyperalgesia. Meanwhile, astrocytes and microglia were dramatically activated in the ipsilateral spinal cord dorsal horn at 10 days after SNL. Further lncRNA microarray and mRNA microarray analysis showed that the expression profiles of lncRNA and mRNA between SNL and sham-operated mice were greatly changed at 10 days. The 511 differentially expressed (>2 fold) lncRNAs (366 up-regulated, 145 down-regulated) and 493 mRNAs (363 up-regulated, 122 down-regulated) were finally identified. The expression patterns of several lncRNAs and mRNAs were further confirmed by qPCR. Functional analysis of differentially expressed (DE) mRNAs showed that the most significant enriched biological processes of up-regulated genes in SNL include immune response, defense response, and inflammation response, which are important pathogenic mechanisms underlying neuropathic pain. 35 DE lncRNAs have neighboring or overlapping DE mRNAs in genome, which is related to Toll-like receptor signaling, cytokine–cytokine receptor interaction, and peroxisome proliferator-activated receptor signaling pathway. Conclusion Our findings uncovered the expression pattern of lncRNAs and mRNAs in the mice spinal cord under neuropathic pain condition. These lncRNAs and mRNAs may represent new therapeutic targets for the treatment of neuropathic pain.


Background
Neuropathic pain is one of the most common chronic pain in humans and characterized by an increase in the responsiveness of nociceptive neurons in the peripheral and central nervous system (CNS) [1]. Peripheral and central sensitization represents the altered functional status of nociceptive neurons and results from changes of a vast amount of functional protein and signaling pathways in the neuron and glial cell [2,3]. Recent pharmaceutical research and discovery activities focus on well-characterized molecular targets, such as ion channels, G-protein-coupled receptors, and kinases in neurons and glial cells localized along the nociceptive pathways, which are regarded as direct contributors to the sensitization of pain signaling systems [4,5]. However, the transcriptional Open Access *Correspondence: gaoyongjing@hotmail.com or translational regulatory mechanisms underlying the expression and functional changes of these molecules are poorly defined.
RNAs that do not code for a protein (noncoding RNAs, ncRNAs) consist of two major classes: the small ncRNAs, which include microRNAs (miRNAs) and other noncoding transcripts of less than 200 nucleotides, and long noncoding RNAs (lncRNAs), which are a novel class of non-protein coding transcripts longer than 200 nucleotides [6]. LncRNAs were initially considered as transcriptional by-products, but recent data suggest that lncRNAs can regulate gene expression via interfering with transcription, post-transcriptional processing, chromatin remodeling, miRNA sequestration, and generating small ncRNAs [7,8]. Also, lncRNAs are involved in various aspects of cell biology and disease etiology, such as development [9], immune [10], cardiovascular disease [11], oncogenesis [12], and neurological disease [13]. LncR-NAs are highly expressed in the CNS, and their expression profiles are associated with specific neuroanatomical regions, cell types, or subcellular compartments suggesting their potential functional roles in the nervous system [14][15][16]. It was reported that sciatic nerve resection induced differential expression of lncRNAs in dorsal root ganglia (DRG) [17]. Moreover, Zhao et al. have recently identified a functional lncRNA Kcna2, which contributed to neuropathic pain by silencing Kcna2 in DRG neurons [18]. These findings indicate the involvement of lncRNAs in neuropathic pain.
The spinal cord is responsible for receiving input from nociceptors and projecting to the brain, and plays an important role in the integration and modulation of pain-related signals. To clarify the molecular mechanisms underlying neuropathic pain and explore novel approaches for analgesic strategies, herein, we investigated the genome-wide expression of lncRNAs in the spinal cord following L5 spinal nerve ligation (SNL)-induced neuropathic pain. We found a large number of differentially expressed (DE) lncRNAs and mRNAs in the spinal cord after SNL. Among them, 39 correlated lncRNA-mRNA pairs, consisting of DE lncRNAs and mRNAs with adjacent or overlapping position relationship, were screened out. Our findings will provide new insights into the roles of lncRNAs in the regulation of neuropathic pain-associated genes.

Model identification of neuropathic pain
The SNL model has been widely used in the investigation of the mechanisms underlying neuropathic pain [19]. Here we also found that SNL induced rapid (1 d) and persistent (>21 d) mechanical allodynia ( Figure 1a) and heat hyperalgesia (Figure 1b) in mice. We then harvested the spinal cord at 10 days (maintenance phase) after SNL and checked the expression of astrocytic marker GFAP and microglial marker IBA-1, which are known to be upregulated in the spinal cord under neuropathic pain condition [20,21]. As shown in Figure 1c, d, GFAP expression and IBA-1 expression were both increased in the ipsilateral dorsal horn in SNL animals but not in sham-treated animals, indicating that glial activation was induced in the spinal cord by SNL.

Overview of lncRNAs and mRNA expression profiles after SNL
We then detected the expression profiles of lncRNAs and mRNAs in the L5 spinal cord at 10 days after SNL by microarray. First, we obtained a graphically overview of the expression signatures of lncRNAs and mRNAs by using scatter plot and hierarchical clustering analyses. The scatter plots showed that a large number of lncRNAs and mRNAs were differentially expressed between SNL and sham-operated mice (Figure 2a, b). Hierarchical cluster analysis of all lncRNAs or mRNA showed that the 3 sham or 3 SNL samples were clustered together respectively, and signal intensity was consistent in sham or SNL group (Figure 2c, d). The heatmap of DE lncRNAs or mRNAs whose expression were up-regulated or down-regulated by twofold were magnified (Figure 2e, f), indicating the high level of concordance in either SNL or sham samples. These data suggest that neuropathic pain is associated with the changes of lncRNAs and mRNAs in the spinal cord.

Real-time quantitative PCR (qPCR) validation of lncRNA and mRNA expression
To validate the reliability of the microarray results and also analyze the temporal changes of lncRNA and mRNA . Data are expressed as mean ± SEM (n = 5 for each group). ***P < 0.001, two-way repeated measures ANOVA. c, d Representative images of GFAP and IBA-1 immunofluorescence in the L5 spinal cord from sham and SNL mice. GFAP and IBA-1 immunoreactive were very low in sham-treated mice, but significantly increased in the ipsilateral superficial dorsal horn at 10 days after SNL. or mRNA (f), whose expression changes were more than twofold. In clustering analysis, up-and down-regulated genes are colored in red and green, respectively. expression after SNL, the up-regulated lncRNAs including Speer7-ps1 and uc007pbc.1, the down-regulated lncRNAs, including ENSMUST00000171761 and ENS-MUST00000097503, the up-regulated mRNA Cyp2d9, and the down-regulated mRNA Mnx1 were randomly selected and analyzed by qPCR. The spinal cord tissues were collected from naïve animals, and SNL animals at 1, 3, 10, and 21 days. Speer7-ps1 and uc007pbc.1, which are intergenic lncRNAs, were both significantly increased at 10 days and peaked at 21 days (Figure 3a, b). ENS-MUST00000171761 and ENSMUST00000097503 are antisense overlap and bidirectional lncRNA with matching gene Tagap (T-cell activation Rho GTPase-activating protein) and Zfp236 (zinc finger protein 236). They were significantly decreased at 10 days and persisted to 21 days (Figure 3c, d). Cyp2d9, a member of cytochrome P450, family 2, subfamily d, was increased more than 12-fold at 10 days (Figure 3e). Mnx1 is a sequence-specific DNA binding transcription factor. It decreased from 1 to 21 days (Figure 3f ). In addition, the fold changes of these lncRNAs and mRNAs detected by qPCR at SNL 10 days were consistent with the results from microarray (Figure 3g), further supporting the reliability of the array data.

Class distribution of changed LncRNAs
lncRNAs were shown to regulate the expression of adjacent or overlapping mRNAs in genome [18,27,28].
Thus, the associations of DE lncRNAs with coding genes were analyzed and classified according to the method described by Li et al. [29]. LncRNAs are classified into four groups: intergenic lncRNAs (lncRNAs are located and transcribed from intergenic regions, and do not overlap with known protein coding genes or other types of genes in genome. It is also called lincRNAs), antisense lncRNAs (LncRNA exon is transcribed from the antisense strand and overlaps with a coding transcript exon), sense lncRNAs (LncRNA exon overlaps with a coding transcript exon on the same genomic strand), and bidirectional lncRNAs (LncRNA is oriented head to head with a coding transcript within 1,000 bp). As shown in

Functional prediction of DE mRNAs in SNL
To explore the molecular mechanism in neuropathic pain, we further did GO and pathway analysis of deregulated genes in SNL versus sham. The GO results showed that the most significant enriched molecular function of up-regulated genes in SNL was chemokine activity, CCR chemokine receptor binding, chemokine receptor binding, and cysteine-type endopeptidase inhibitor activity ( Figure 5a). The most significant enriched biological processes of up-regulated genes in SNL were immune response, immune system process, defense response, and regulation of immune system process ( Figure 5b).
The most noteworthy enriched cellular components of up-regulated genes in SNL were extracellular region, extracellular space, extracellular region part, and external side of plasma membrane ( Figure 5c). The most significant enriched molecular function of down-regulated genes in SNL were binding, receptor binding, calcium ion binding, and tropomyosin binding (Figure 5d). The most significant enriched biological processes of downregulated genes in SNL were regulation of ATPase activity, monovalent inorganic cation transport, glucosamine-containing compound catabolic process, and amino sugar catabolic process (Figure 5e). The most significant enriched cellular components of down-regulated genes in SNL were extracellular region, striated muscle thin filament, extracellular space, and cell part (Figure 5f ). Similarly, different genes were analyzed in KEGG. The results showed that the up-regulated genes in SNL are involved in complement and coagulation cascades, Tolllike receptor signaling pathway, chemokine signaling pathway, cytosolic DNA-sensing pathway, and cytokinecytokine receptor interaction, Changas disease, and , lncRNA ENSMUST00000171761 (c), and lncRNA ENSMUST00000097503 (d) were significantly deregulated at 10 and 21 days after SNL. e The expression of Cyp2d9 mRNA was markedly up-regulated at 10 days after SNL. f The expression of Mnx1 mRNA was significantly down-regulated at 1, 3, 10 and 21 days after SNL. One-way ANOVA followed by Tukey's multiple comparison test. *P < 0.01, **P < 0.01, ***P < 0.001. g Log 10 value of signal intensity detected by microarray.

Comparison of our DE mRNAs with previously published microarrays
Previous studies have shown differential gene expression profile in the spinal cord in rats with neuropathic pain [30,31]. In order to compare neuropathic pain-associated gene expression patterns in mice and rats, we did the overlap analysis between other's microarray data from rat [30] and our current data from mice ( Figure 7a). LaCroix-Fralish et al. reported that 88 genes were upregulated and 83 genes were downregulated in the spinal cord 7 days after L5 nerve root ligation in rats [30]. Surprisingly, compared to 361 up-regulated genes and 119 down-regulated genes in mouse, only 1 gene (Cd74) was upregulated and 2 genes (Nefm, Aco2) were downregulated in both rats and mice (Figure 7b). In addition, we compared our array data with 79 significantly regulated genes which were identified by meta-analysis from 20 independent microarray experiments from rats and mice after tissue inflammation or nerve injury [2]. We observed an overlap of 15 genes with the meta-analysis dataset (Figure 7c). These genes included 14 up-regulated genes (Ctss, C1qb, C1qc, Npy, Cd74, Gal, Aif1, Calca, Cxcl10, Atf3, Ccl2, Ctsh, Fcgr2b and Sprr1a) and 1 down-regulated gene (Nefm) (Figure 7d).

Relational analysis of lncRNAs and mRNAs
As some lncRNAs have been suggested to play key roles in regulating the expression of their neighboring or overlapping genes in genome wide, we further screened out DE mRNAs related to DE lncRNAs based on their location distributions on mouse chromosomes by UCSC Genome Browser. In the spinal cord, there are 39 DE lncRNA-mRNA pairs for 35 DE lncRNAs and 35 DE mRNAs. Among them, 32 pairs exhibited coordinated expression changes, and 7 pairs were non-coordinated, which may suggest a complex and various regulatory mechanisms across different lncRNAs and their target mRNAs. Intriguingly, all the seven non-coordinated lncRNA-mRNA pairs belong to intergenic lncRNA-mRNA pairs (Table 3). Further GO and pathway analysis showed that the high enriched molecular functions include pheromone binding, chemokine activity, highdensity lipoprotein binding, and phosphatidylcholinesterol O-acyltransferase activator activity (Figure 8a). Based on gene-pathway network graph analysis, we found that the DE mRNAs from lncRNA-mRNA pairs, such as Cxcl9 (chemokine (C-X-C motif ) ligand 9), Cxcl10 (chemokine (C-X-C motif ) ligand 10), Cxcl11 (chemokine (C-X-C motif ) ligand 11), Trhr (thyrotropin releasing hormone receptor), and Apoa2 (apolipoprotein A-II), might involve in toll-like receptor signaling pathway, calcium signaling pathway, and PPAR signaling pathway (Figure 8b; Table 3), which have been proven to be involved in neuropathic pain pathogenesis [32][33][34].

Discussion
Chronic neuropathic pain is a somatosensory disorder caused by nerve injury or disease that affects the nervous system [35]. Evidence suggested that the particular patterns of gene expression at different levels of the nociceptive system play important roles in the development and maintenance of neuropathic pain [2,36]. Over the past decades, the molecular mechanisms underlying neuropathic pain have been extensively studied; however, the pathophysiological process of pain is still vague. LncRNAs were recently shown to regulate gene expression [37] and traffic cellular protein complexes, genes, and chromosomes to appropriate locations [8]. Their function in regulating gene expression switching in the maintenance phase of neuropathic pain is poorly understood. In this study, we for the first time identified the global expression changes in lncRNAs and analyzed their characteristics and possible relation with coding genes in the spinal cord under neuropathic pain condition. The 24,833 lncRNAs were detected in the spinal cord of mice. Among them, 366 lncRNAs were up-regulated and 145 lncRNAs were down-regulated at 10 days after SNL. These DE lncRNAs are consistently altered in a high percentage of analyzed spinal cords from SNL and sham mice, suggesting that lncRNAs may be involved in neuropathic pain processing. So far, most DE lncRNAs have not been functionally characterized. Although it was still too early to translate this knowledge into the development of novel analgesic agents for better pain relief, these findings may likely provide novel insight into the molecular basis of pain.
In this study, the expression profiles of mouse genomewide mRNAs were also detected using lncRNA Microarray Chip at the same time. Among DE mRNAs, the up-regulated mRNAs are far more numerous than the down-regulated in SNL samples, which reflects the emergence of new biology processes and pathways in pathological conditions. A number of reported pain-related genes, including Cacna1g, Trpv1, Ccl5, Cx3cr1 and Irf5 were dramatically increased after SNL. Moreover, a lot of other mRNAs, such as Sprr1a, Anxa10, Kng1, and Gpr151 (G-protein-coupled receptor 151), whose functions are unclear in the spinal cord were also screened out. As the expression changes for some genes may be related to nerve damage and homeostatic responses to denervation, further studies are needed to identify whether they are involved in neuropathic pain processing.     Based on the GO term enrichment analyses of DE mRNA, we found that significantly enriched molecular functions and biological processes of up-regulated gene in SNL vs sham were mainly involved in chemokine activity, inflammation, and immunity. These findings are consistent with previous studies showing that neuroinflammation, manifested as infiltration of immune cells [38], activation of glial cells [39] and production of inflammatory mediators [40] in the peripheral and CNS, plays an important role in the induction and maintenance of chronic pain [41]. Additionally, our immunostaining of GFAP and IBA-1 showed dramatic glial activation in the spinal cord at 10 days after SNL. From significant pathway analyses of DE gene, the third most significant enriched pathway of the up-regulated genes in SNL vs sham is the toll-like receptor signaling pathway. Indeed, Tlr2 [42], Tlr4 [43], and Tlr7 [44] have been implicated as potential therapeutic targets in neuropathic and other pain models. The data collectively indicate that anti-neuroinflammation may be an effective strategy for the treatment of neuropathic pain.
Previous studies utilizing cDNA microarrays to analyze gene expression profiles primarily focus on pain models in rats, rarely in mice [2]. The overlap analysis showed little overlap between rat and mice spinal cord gene expression patterns under neuropathic pain states, suggesting the species difference in gene expression. However, we found that there were 15 overlap genes between our current data and meta-analysis results reported by LaCroix-Fralish et al. [2]. These overlap genes including Atf3, Sprr1al and Nefm can be induced by nerve damage, which contribute to chronic pain [45][46][47]. In addition, gene ontology-based functional annotation clustering analyses of the previous gene chip study revealed strong evidence for regulation of immune-related genes in pain states, which was consistent with our data.
Although lncRNAs play important roles in the regulation of gene expression [48], there is a large gap between the number of existing lncRNAs and their known association with a particular molecular or cellular function [49]. Regulatory mechanisms and major functional principles of lncRNAs are complex and quite obscure. Unlike microRNA, there are no common languages that can be used to predict lncRNAs' target genes and function by their sequence information or secondary structure. Accumulating evidence suggests that a number of lncRNAs function locally to activate or repress their neighboring or overlapping genes' expression [18,27,50]. In this study, we found that intergenic lncRNAs (lincRNAs) were the largest category in all DE lncRNAs after SNL. In reality, lincRNAs are found to be conserved across multiple vertebrate species [51] and perform important functions  Table 3. The DE lncRNAs-related genes and the corresponding pathways were shown in the circles and boxes, respectively. The color of pathway terms is defined by the enrichment P value.
in many cellular processes, from cell proliferation to cancer progression [52]. Furthermore, lincRNAs can function through different types of mechanisms, including cis or trans transcriptional regulation, translational control, splicing regulation, and other post-transcriptional regulation [33]. We examined whether their neighboring or overlapping protein-coding genes in the genome are simultaneously DE in the spinal cord after SNL, and found that there are 39 DE lncRNA-mRNA pairs. Our further analysis showed that an up-regulated lincRNA, MM9LINCRNAEXON10576− in the spinal cord after SNL was found to be located near Cxcl10, Cxcl9 and Cxcl11 gene cluster in mice chromosome 5. All the four RNAs have the same expression trends and increased more than twofold after SNL. Recently, studies using animal models have shown that upregulation of chemokines in the spinal cord play a vital role in the development and maintenance of chronic pain [41,53,54]. Indeed, recent research found that Cxcl10 and its receptor Cxcr3 were involved in inflammatory pain and cancer pain [55][56][57]. Therefore, lncRNA MM9LINCRNAEXON10576− may contribute to neuropathic pain through regulation of chemokines Cxcl10, Cxcl9 and Cxcl11. In our microarray results, 12 DE mRNA have their corresponding DE sense-overlap lncRNAs, and the change patterns of these lncRNA were same as that of their accompanying protein-coding genes. Di et al. found that a sense-overlap lncRNA arising from the CCAAT/ enhancer-binding protein alpha (Cebpa) gene locus can bind to DNA methyltransferase 1 (DNMT1) and prevent Cebpa gene locus methylation, then to increase the expression of Cebpa gene. Their deep sequencing of transcripts associated with DNMT1 combined with genomescale methylation and expression profiling extend the generality of this finding to numerous gene loci. [27]. Given that the 12 DE mRNA and their DE sense-overlap lncRNAs were both increased after SNL, it's possible that the DE sense-overlap lncRNAs regulate the expression of their sense-overlapping mRNAs via demethylation after SNL.

Conclusion
Our results demonstrated that lncRNA transcripts were highly enriched and hundreds of lncRNAs were differentially expressed in the spinal cord after SNL. Dozens of DE lncRNAs were observed to have neighboring or overlapping DE mRNAs in genome. These lncRNAs may locally regulate their related protein-genes expression and play key roles in the pathogenesis of neuropathic pain. Further studies are required to clarify the molecular and cellular functions of DE lncRNAs and determine whether they can serve as novel analgesic targets in neuropathic pain.

Animals and surgery
Adult male ICR mice (male, 8 weeks) were maintained on a 12:12 light-dark cycle at a room temperature of 22 ± 1°C with free access to food and water. The experimental procedures were approved by the Animal Care and Use Committee of Nantong University and performed in accordance with the guidelines of the International Association for the Study of Pain. To produce a SNL, animals were anesthetized with isoflurane and the L6 transverse process was removed to expose the L4 and L5 spinal nerves. The L5 spinal nerve was then isolated and tightly ligated with 6-0 silk threads [58]. For sham operations, the L5 spinal nerve was exposed but not ligated.

Behavioral test
Animals were habituated to the testing environment daily for at least 2 days before baseline testing. The room temperature remained stable for all experiments. For testing mechanical sensitivity, animals were put in boxes on an elevated metal mesh floor and allowed 30 min for habituation before examination. The plantar surface of each hindpaw was stimulated with a series of von Frey hairs with logarithmically incrementing stiffness (0.02-2.56 g, Stoelting, Wood Dale, IL, USA), presented perpendicular to the plantar surface (2-3 s for each hair). The 50% paw withdrawal threshold was determined using Dixon's updown method [59]. For testing heat sensitivity, animals were put in plastic boxes and allowed 30 min for habituation. Heat sensitivity was tested by radiant heat using Hargreaves apparatus (IITC Life Science Inc., Woodland Hills, CA, USA) and expressed as paw withdrawal latency (PWL). The radiant heat intensity was adjusted so that basal PWL is between 10 and 14 s, with a cutoff of 18 s to prevent tissue damage.

Immunohistochemistry
At 10 days after SNL or sham-operation, animals were deeply anesthetized with isoflurane and perfused through the ascending aorta with PBS followed by 4% paraformaldehyde with 1.5% picric acid in 0.16 M PB. After the perfusion, the L4-L5 spinal cord segments were removed and postfixed in the same fixative overnight. Spinal cord sections (30 μm, free-floating) were cut in a cryostat. The sections were first blocked with 5% goat serum for 2 h at room temperature. The sections were then incubated overnight at 4°C with the following primary antibodies: GFAP antibody (mouse, 1:6,000; Millipore, Billerica, MA, USA), IBA-1 antibody (Mouse, 1:3,000, Serotec, Kidlington, UK). The sections were then incubated for 2 h at room temperature with FITC-conjugated secondary antibodies (1:1,000, Jackson ImmunoResearch). The stained sections were examined with a Leica fluorescence microscope, and images were captured with a CCD Spot camera.

Tissue collection and RNA isolation
We prepared nine mice for SNL and nine mice for shamoperation. At 10 days after operation, the animals were deeply anesthetized with isoflurane and perfused through the ascending aorta with saline. After the perfusion, the L4-L5 spinal cord segments were collected. Total RNA was extracted from the spinal cord dorsal horn tissue using Trizol reagent (Invitrogen, Carlsbad) according to the manufacturer's protocol. The RNA concentration and purity were assayed by the absorbance values at 260 and 280 nm using the NanoDrop 1000 Spectrophotometer (Thermo). RNA integrity was checked by electrophoresis on 2% (m/v) agarose gels. After these testing, equal mRNA from three mice under the same treatment was mixed as one sample. Therefore, six samples (3 for SNL and 3 for sham) were sent for microarray analysis.

Microarray assay
The gene chip of the mouse lncRNA microarray V2.0 (8 × 60K, Arraystar), which includes 25,376 lncRNA probes and 31,423 coding gene probes, was used in the experiments. The total RNAs of sham and SNL groups were individually hybridized with gene chips. Briefly, RNA was purified from 1 μg total RNA after removing rRNA. The RNA sample was then transcribed into fluorescent cRNA along the entire length of the transcripts without 3′ bias utilizing random primers. The labeled cRNAs were hybridized to mouse lncRNA microarray. Finally, arrays were scanned by Agilent Scanner G2505B. The array images were analyzed by Agilent Feature Extraction software (version 10.7.3.1). The GeneSpring GX v11.5.1 software package (Agilent Technologies) was utilized to analyze quintile normalization and subsequent data processing. The microarray hybridization was carried out by Kangchen Bio-tech, Shanghai, China.

Bioinformatics analysis
Differentially expressed lncRNAs and mRNAs with statistical significance were identified through Volcano Plot filtering. The threshold used to screen up-or down-regulated RNAs was fold-change >2.0 (P < 0.05). Hierarchical clustering was carried out by Cluster 3.0, and the heat maps were generated in Java Treeview. The DE mRNAs which were adjacent to or overlap with the DE lncRNAs were recognized as DE lncRNAs related mRNAs using UCSC Genome Browser. The differentially expressed mRNAs or DE lncRNAs related mRNAs were analyzed by pathway annotation and gene ontology (GO) functional enrichment using CapitalBio ® Molecule Annotation System V3.0 (MAS3.0). The −log 10 (P-value) of the GO and pathway results were shown in the histogram.

Real-time reverse transcription-polymerase chain reaction (RT-PCR)
The microarray results were confirmed by RT-PCR. Total RNA was extracted from the spinal cord tissue as described above and total RNA was reverse transcribed using random hexamers primer (TaKaRa Bio Inc) according to the manufacturer's description. The expression level of six genes was checked, including Speer7-ps1, uc007pbc.1, ENSMUST00000171761, ENS-MUST00000097503, Cyp2d9, and Mnx1. The Gapdh was used as house-keeping gene. The sequences of all primers were shown in Table 4. RT-PCR was performed using the Fast Start Universal SYBR Green Master (TaKaRa Bio