Network and pathway‐based analysis of microRNA role in neuropathic pain in rat models

Abstract The molecular mechanisms underlying neuropathic pain (NP) remain poorly understood. Emerging evidence has suggested the role of microRNAs (miRNAs) in the initiation and development of NP, but the specific effects of miRNAs in NP are largely unknown. Here, we use network‐ and pathway‐based methods to investigate NP‐induced miRNA changes and their biological functions by conducting a systematic search through multiple electronic databases. Thirty‐seven articles meet the inclusion criteria. Venn analysis and target gene forecasting are performed and the results indicate that 167 overlapping target genes are co‐regulated by five down‐regulated miRNAs (rno‐miR‐183, rno‐miR‐96, rno‐miR‐30b, rno‐miR‐150 and rno‐miR‐206). Protein‐protein interaction network analysis shows that 77 genes exhibit interactions, with cyclic adenosine monophosphate (cAMP)‐dependent protein kinase catalytic subunit beta (degree = 11) and cAMP‐response element binding protein 1 (degree = 10) having the highest connectivity degree. Gene ontology analysis shows that these target genes are enriched in neuron part, neuron projection, somatodendritic compartment and nervous system development. Moreover, analysis of Kyoto Encyclopedia of Genes and Genomes reveals that three pathways, namely, axon guidance, circadian entrainment and insulin secretion, are significantly enriched. In addition, rno‐miR‐183, rno‐miR‐96, rno‐miR‐30b, rno‐miR‐150 and rno‐miR‐206 are consistently down‐regulated in the NP models, thus constituting the potential biomarkers of this disease. Characterizing these miRNAs and their target genes paves way for their future use in clinical practice.


| INTRODUC TI ON
Neuropathic pain (NP) results from damage to or diseases of the somatosensory system and is often accompanied by maladaptive changes in the nervous system. 1 According to a systematic review of epidemiological studies released in 2014, 7%-10% of the general population had suffered from NP. 2 Neuropathic pain has multiple causes, including trauma, metabolic diseases, infection, tumour invasion or neurotoxic chemicals. 3 The characteristics of chronic NP are complex symptoms, difficult treatments and poor outcomes, which increase the patients' burden and cause anxiety and depression. 4 To date, the molecular mechanisms on NP have remained poorly understood. Mechanism-based treatments having an advantage over disease-or cause-based treatments might be the main reason why NP is often treated inadequately or ineffectively.
MicroRNAs (miRNAs) are endogenous small non-coding RNAs that regulate gene expression by inhibiting protein synthesis. 5 Evidence suggests that the mechanism for NP involves miRNAs. 6 Recent studies revealed that miRNAs are highly expressed in the sensory organs of the nervous system, such as dorsal root ganglion (DRG) and spinal dorsal horn (SDH). [7][8][9] For instance, miR-183 is significantly down-regulated in the DRG of a spinal nerve ligation (SNL) rat model. The overexpression of miR-183 attenuates SNL-induced mechanical allodynia and this effect is closely connected with the substantial down-regulation of voltage-gated sodium channel Nav1.3 and brain-derived neurotrophic factor (BDNF). 8 These findings suggest that miRNAs play important roles in the developmental and biological functions of the nervous system for NP.
In 2018, a bioinformatic analysis of miRNAs related to NP used a microarray profile from Gene Expression Omnibus database and included peripheral blood samples from patients with NP after spinal cord injury. 10 Our study focused on frequently used rat models, such as SNL, spared nerve injury (SNI) and sciatic chronic constriction injury (CCI), to study the miRNA mechanism of NP. We performed a wide literature search on the subject and a subsequent bioinformatics analysis. To the best of our knowledge, this comprehensive bioinformatics analysis is the first to explore the biological functions of miRNAs and identify the potential therapeutic targets for NP.

| Search strategy
A systematic search was conducted from inception to May 2018.
We searched PubMed, EMBASE, Web of Science and Ensco.
Studies were identified using the following criteria: (a) type of studies: such as only original articles investigating the role of miR-NAs in NP by comparing the animal models of NP to those without pain; (b) type of animal models: rat models of NP including SNL, SNI, CCI, partial sciatic nerve injury and chronic compression of DRG (CCD); (c) type of samples such as nervous tissues (eg, sciatic nerve, DRG, spinal cord and brain) and nervous cells (eg, DRG neurons, microglia and astrocytes); and (d) type of measurements such as miRNA expression assessed by polymerase chain reaction, TaqMan low density array (TLDA) or microarray analysis.

| Data extraction
Two researchers (J-BG and YZ) independently reviewed and extracted data from studies to evaluate their eligibility for inclusion.
We extracted data, including the first author, publication year, country, experimental design (eg, experimental models, region used) and information on miRNAs (eg, expression change, target genes and functions). Any discrepancy for selection and extraction was resolved by a third researcher (B-LC).

| Bioinformatics analysis
To examine the functional roles of miRNAs, we predicted the target genes of miRNAs by using the TargetScan software (http://www.targe tscan.org/). Venn diagram analysis showing the number of overlapping miRNAs and targets was based on Functional Enrichment analysis tool (FunRich; http://www.funri ch.org/). The minimum fold value for up-regulated and down-regulated miRNAs in the matrix table is 2.
A protein-protein interaction (PPI) network was used to further understand the correlations between the overlapping targets of differentially expressed miRNAs. Protein-protein interaction data sources were obtained from the String database (http://string-db.org/), 11 and maps were drawn with Cytoscape software v.3.6.0. 12

| RE SULTS
We identified 1438 articles through electronic search. After the exclusion of duplicate records, 1160 articles remained. The researchers then examined the title and abstract and 91 articles were identified for further full-text review. Finally, 37 articles fulfilled the eligibility criteria and were included. [7][8][9] Among these studies, six [11][12][13][14][15][16] were about miRNA profiles and 32 7-9,16-44 were on miRNA experimental verification, of which one 16 was conducted miRNA profile and experimental verification. The flowchart of the study selection procedure is detailed in Figure 1.

| Study characteristics
The main characteristics of the included articles are shown in Tables 1 and 2. As shown in Table 1, six articles 11-16 examined the miRNA changes in NP rat models through TLDA card or microarray analysis. The numbers of the significantly dysregulated miRNAs in these seven studies varied from 1 to 22. In the earliest study, 11 the authors investigated and compared the miRNA expression profile in the spinal cord of rats with CCI with that of sham-operated rats. The results indicated six down-regulated miRNAs. In the latest study, 16 miRNA changes in anterior cingulate cortex (ACC) after CCI were examined by using microarray. Nine miRNAs were significantly up-regulated and 11 were significantly down-regulated.

| Target prediction and Venn diagram analysis
Using the array data from Table 1, we generated a matrix table with a FunRich tool ( Figure 2) and it shows the number and percentage of co-regulated miRNAs through pair-wise comparison. In the matrix table, we found that rno-miR-221, rno-miR-21 and rno-341 were up-regulated in two or more studies, whereas mmu-miR-151-3p and rno-miR-214 were down-regulated in two or more studies. Target Scan software was subsequently used to forecast the target genes of miRNAs and Venn diagrams were drawn. We found 20 overlapping target genes in the three up-regulated miRNAs (rno-miR-221, rno-miR-21 and rno-341), but none in the two down-regulated miRNAs (mmu-miR-151-3p and rno-miR-214) (Figures S1 and S2).
The miR-183 cluster comprises miR-183, miR-96 and miR-182 and shares the same sequence homology. Therefore, we combined the target genes of rno-miR-183 and rno-miR-96 for Venn diagrams analysis.

| PPI network analysis
The PPI data of 167 overlapping target genes were obtained from String database and the network was displayed by using Cytoscape software. The results of PPI analysis are shown in Figure 4. A total of 77 genes exhibited interactions. The sizes and colours of each node represent the degree of functional connection with the 167 genes.
The colours of each edge indicate the strength of data support, as evaluated by combined scores. A low value is represented by small sizes and bright colours in the map and a high value is represented by large sizes and dark colours. The cyclic adenosine monophosphate (cAMP)-response element binding protein1 (CREB1) and cAMP-dependent protein kinase catalytic subunit beta (PRKACB) were the two largest and darkest nodes in the network. Hence, PRKACB (degree = 11) and CREB1 (degree = 10) exhibited the highest connectivity degree. Furthermore, we mapped the miRNAs of Table 2 to their target genes ( Figure 5). Three down-regulated miRNAs (rno-miR-183, rno-miR-96 and rno-miR-30b) directly targeted Nav1.3 and BDNF.
According to the analysis of the up-regulated miRNA and their target gene network, the suppressor of cytokine signalling-1 (SOCS1) is a direct target of rno-miR-155, rno-miR-19a and rno-miR-221.
F I G U R E 1 Flow chart of the study selection procedure (for details of study identification)

| Functional enrichment analysis
We analysed 167 overlapping target genes from five down-regulated miRNAs (rno-miR-183, rno-miR-96, rno-miR-30b, rno-miR-150 and rno-miR-206), as shown in Figure 3. A total of 128 GO terms were significantly enriched (P < 0.05) and the top 10 high enrichment score pathways are shown in Figure 6. The results of the GO analysis comprised biological processes (eg, single-organism developmental  Table S1).

| D ISCUSS I ON
In this study, we focused on miRNA expression in DRG, SDH and ACC. These regions play important roles in the somatosensory pathway from primary sensory neurons to the central nervous system in NP. DRG neurons receive nociceptive afferents carrying peripheral inputs (eg, heat, cold, pressure and chemicals) and transmit information to second-order neurons mostly in SDH. After the integration and processing of the sensory inputs in SDH, the outputs from the spinal networks are carried to the higher cortical centres. 45 Anterior cingulate cortex is a crucial brain region of the limbic system and is associated with the anticipation of pain and attention to pain. 46,47 MicroRNAs in NP are involved in neuroinflammation, neuronal excitability, neuronal plasticity and DNA methylation (Figure 7). Several studies focused on neuronal excitability mechanisms. Peripheral nerve injury induces the hyperexcitability of injured afferent neurons, thereby contributing to ectopic discharge. 48 Table 1 and shows the number and percentage of coregulated miRNAs. The minimum fold value for up-regulated and down-regulated miRNAs is 2. miRNAs, microRNAs F I G U R E 3 Venn diagram analysis. Overlapping target genes of rno-miR-183/96, rno-miR-30b, rno-miR-150 and rno-miR-206. These five down-regulated miRNAs were identified in Table 2 and have been observed in two or more studies. miR-183 and miR-96 belong to the miR-183 cluster, so we combined their target genes to easily create the Venn diagram. miRNAs, microRNAs F I G U R E 4 Protein-protein interaction network analysis. A total of 167 overlapping target genes were obtained from Figure 3. The sizes and colours of each node represent the degree of functional connection with these genes. The colours of each edge indicate the strength of data support as evaluated by combined scores. A low value is represented by small sizes and bright colours in the map and a high value is represented by large sizes and dark colours F I G U R E 5 Dysregulated miRNA-target gene network. The network is based on the dysregulated miRNAs and their target genes identified in Table 2. Green colours represent down-regulated miRNAs and red colours represent up-regulated miRNAs. miRNAs, microRNAs network analysis with pathway analysis is able to be more robust for possible false positives resulted from various miRNAs in different studies.
In summary, our study provides an important step for describing the NP sensitivity of miRNA expression and the target gene regulatory consequences of related expression changes.
F I G U R E 6 GO annotation enrichment analysis. The vertical axis is the description of GO terms and the horizontal axis is the enrichment score (-log10[P-value]) of the pathways; log10[Pvalue] is the logarithm of the P-value and P < 0.05 was considered significant. GO, gene ontology F I G U R E 7 Ascending NP pathway and miRNA regulation of pain genes. Noxious stimulation reaches SDH though afferent nerve fibre. SDH is the site responsible for integrating and processing information of sensory inputs and carries the output to the brain by several pathways. The thalamus and limbic system are important sites in the brain for ascending the NP pathway. miRNA modulation is reflected in DRG, SDH and several brain areas. A, miRNAs, their target genes and functions in DRG. B, miRNAs, their target genes and functions in spinal cord. C, miRNAs, their target genes and functions in brain. miRNAs, microRNAs; NP, neuropathic pain; DRG, dorsal root ganglion; SDH, spinal dorsal horn BDNF, brain-derived neurotrophic factor; RASA1, RAS p21 protein activator 1; DNM3a, DNA methyltransferase 3a; AMPA receptor subunit, GluA1 and GluA2; HMGB1, High mobility group box 1; TREK-1, TWIK-Related K+ Channel 1; NR2B, N-methyl-D-Aspartate receptors 2B; ATG14, Autophagy Related Gene 14; Dusp5, Dual-specificity phosphatase 5; Rap1a, Ras-related protein Rap-1A; SOCS, Suppressor of cytokine signaling; RREB1, Ras responsive element binding protein 1; p-AKT, phosphorylated-protein kinase B; mTOR, mammalian target of rapamycin; ZEB1, Zinc finger E box binding protein 1; STAT3, signal transducer and activator of transcription 3; TNFAIP1, tumor necrosis factor alpha-induced protein 1; TLR5, toll-like receptor 5 Bioinformatics analysis elucidates the functional roles of five miRNAs (rno-miR-183, rno-miR-96, rno-miR-30b, rno-miR-150 and rno-miR-206) and their targets involved in the known relevant pathways for NP. They may serve as potential biomarkers and novel strategies for prevention and treatment of NP following peripheral nerve injury. Further studies are necessary to confirm the involvement of the five miRNAs in different stages of NP and translate the findings to clinical practice.

ACK N OWLED G EM ENTS
This work was supported by the National Natural Science

CO N FLI C T S O F I NTE R E S T
The authors declare that they have no competing interests.