Proteomics Analysis of the Spinal Dorsal Horn in Diabetic Painful Neuropathy Rats With Electroacupuncture Treatment

Background Clinical evidence demonstrates that electro-acupuncture (EA) of the Zu sanli (ST36) and Shen shu (BL23) acupoints is effective in relieving diabetic painful neuropathy (DPN); however, the underlying molecular mechanism requires further investigation, including the protein molecules associated with EA’s effects on DPN. Methods Sprague-Dawley adult male rats (n =36) were randomly assigned into control, DPN, and EA groups (n=12 each). After four weeks of EA treatment, response to mechanical pain and fasting blood glucose were analyzed. A tandem mass tag (TMT) labeling approach coupled with liquid chromatography with tandem mass spectrometry was used to identify potential biomarkers in the spinal dorsal horn. Further, proteomics analysis was used to quantify differentially expressed proteins (DEPs), and gene ontology, KEGG pathways, cluster, and string protein network interaction analyses conducted to explore the main protein targets of EA. Results Compared with the DPN model group, the mechanical pain threshold was significantly increased, while the fasting blood glucose levels were clearly decreased in EA group rats. Proteomics analysis was used to quantify 5393 proteins, and DEPs were chosen for further analyses, based on a threshold of 1.2-fold difference in expression level (P < 0.05) compared with control groups. Relative to the control group, 169 down-regulated and 474 up-regulated proteins were identified in the DPN group, while 107 and 328 proteins were up- and down-regulated in the EA treatment group compared with the DPN group. Bioinformatics analysis suggested that levels of proteins involved in oxidative stress injury regulation were dramatically altered during the EA effects on DPN. Conclusions Our results provide the valuable protein biomarkers, which facilitates unique mechanistic insights into the DPN pathogenesis and EA analgesic, antioxidant stress and hypoglycemic effect.


INTRODUCTION
There are various types of chronic pain that can disturb the physical and mental health of individuals. Diabetic painful neuropathy (DPN) is a painful complication of diabetes, caused by constant high glucose levels. Approximately 20% of all diabetic patients present with neuropathic pain, manifesting as spontaneous pain, hyperalgesia, and allodynia (1). The streptozotocin (STZ)-diabetic rat model can be used to mimic the symptoms and pathological changes that occur during DPN. To date, the majority of studies have focused on peripheral sensory nerve function; however, the spinal cord is also important in the development of DPN (2). A number of morphological studies have demonstrated spinal cord damage and atrophy in patients with diabetes, and peripheral nerve stimulation can reduce blood-oxygen level-dependent activity in the dorsal horn of diabetic rats (3,4).
Proteomics analysis provides a valuable strategy for exploring the pathogenesis of diabetes mellitus, as well as therapeutic targets in this condition. Some proteomic studies have been performed in animal models of DPN or diabetic neuropathy, to search for differentially expressed proteins (DEPs) induced by STZ. Mitochondrial oxidative phosphorylation, metabolic dysregulation of mitochondrial complexes, and antioxidases are involved in DPN progression and represent potentially important therapeutic targets. However, previous studies have focused only on the proteomes of the sciatic nerve, lumbar 4/5 dorsal root ganglia, and trigeminal ganglia, while the exact spinal mechanism underlying DPN is far from clear and additional potential spinal biomarkers are needed (11)(12)(13)(14).

Animals and Ethics Statement
Male Sprague-Dawley rats weighing 160 ± 20g were obtained from the Experimental Animal Center of Fujian University of Chinese Medicine. Rats were housed under standard living conditions at 22°C and a 12 h light/dark cycle with sufficient food and water and acclimatized to the pain-testing environment before the experiments. All animal procedures were in strict accordance with relevant international laboratory animal use guidelines. After adaptive feeding, rats were randomly allocated into three groups as follows: Control group (n=12), DPN group (n=12), and EA group (n = 12).
According to the previous literature reports, diabetes was induced by a single intraperitoneal injection of STZ (65 mg/kg, Sigma Chemicals, USA) dissolved in citrate buffer (10 mmol/L, Na citrate, pH=4.3) (15). The control group received an equivalent volume of citrate buffer only. Two weeks later, tail venous fasting blood glucose (FBG) were measured using a glucometer (OEM, USA). Only rats with FBG >16.7 mmol/L and mechanical paw withdrawal threshold (PWT) <5 g were included in further experiments.
All experimental procedures and protocols were approved by the Animal Research and Ethics Committee of Fujian University of Traditional Chinese Medicine. At the endpoint of the experiment, all rats were euthanized, according to the Care Guidelines.

Electroacupuncture (EA) Procedure
In the EA group, the rats were loosely fixed on a wooden stand, where their head and limbs could move freely. After the rats were stabilized, the acupoints "Zu sanli" (ST36) and "Shen shu" (BL23) were selected for acupuncture and electrical stimulation ( Figure 1A). EA intervention began at two weeks after STZ injection. ST36 is located at the posterolateral knee about 5mm below the fibula head (16), while BL23 is located in the depressions lateral to the lower border of the spinous processes of the second lumbar vertebrae, approximately 8mm from the midline of the adult rat (17). Needles inserted at ST36 and BL23 were connected to a G6805-1A multifunctional EA apparatus (Shanghai Medical Electronic Apparatus Company, Shanghai, China), with a stimulation intensity of 1mA, frequency 10 HZ. According to previous literature reports, EA stimulation was maintained for 30 min once every other day for a consecutive four weeks (15,18).

Measurement of Paw Withdrawal Threshold
The mechanical PWT was examined before and after injection of STZ or citrate buffer (0, 2, 6 weeks after injection). The test environment was a 10 cm × 20 cm × 20 cm plastic cage with a plexiglass floor. After the rats were placed in the cage and allowed to stand for 10 minutes, the center of the feet of the rats were stimulated with different tense by von Frey Hairs (0.4, 0.6, 1.4, 2, 4, 6, 8, 10, and 15g) according to the up-down method described by Dixon (19). Each von Frey hair was held for about 1-2 s, with 10min intervals between each stimulation. Testing was initiated with the 2.0 g hair.
The positive response was defined as a withdrawal of hind paw upon the stimulus. Whenever a positive response appeared, the next lower gram was applied, and whenever a negative response occurred, the next higher gram was applied. The testing consisted of five more stimuli after the first change in response occurred, and the pattern of response was converted to a 50% von Frey threshold using the method described by Dixon and based on our previous literature (20,21). All data are presented as mean ± SE.

Sample Collection
To investigate DEPs in three groups, all rats were sacrificed under low dose isoflurane (0.5-1%) anesthesia, six weeks after STZ injection. The L 4-5 spinal cords were rapidly dissected and immediately placed in -80°C refrigerator for the further proteomics detection.

Protein Extraction and Trypsin Digestion
The L 4-5 spinal cord tissue samples were removed from -80°C freezer, and an appropriate amount of sample was added to a pre-chilled mortar with liquid nitrogen and ground to powder. Proteins were then extracted from samples in lysis buffer (8 M urea, 1% protease inhibitor) using a high-intensity ultrasonic processor. Remaining debris was removed by centrifugation at 12,000g for 10min at 4°C. The protein concentration of the supernatant was detected with BAC kits.  For trypsin digestion, protein solution was reduced with dithiothreitol (final concentration, 5mM) for 30 min at 56°C, then incubated with iodoacetamide (final concentration, 11mM) for 15min at 37°C in darkness. Then, the protein samples were diluted with mM triethylamonium bicarbonate (TEAB), until the urea concentration was < 2 M. Finally, trypsin was added at a mass ratio (trypsin: protein) of 1:50 for the first digestion at 37°C overnight, then of 1:100 for the second digestion (4h, 37°C).

TMT Labeling
After trypsin digestion, peptides were desalted and vacuum-dried with the strata X C18 SPE column. Peptides were reconstituted in 0.5M TEAB and marked according to the TMT kit protocols.

LC-MS/MS Analysis
Separated peptides were subjected to sodium/iodide symporter sources. Tandem mass spectrometry (MS/MS) was performed using a Q ExactiveTM Plus (Thermo) instrument, coupled online to an ultra-performance liquid chromatography system.
Electrospray voltage was 2.0 kV, and intact peptides were detected in the Orbitrap at 70,000 mass resolution. The primary MS scan range was 350-1600 m/z. Collected data were processed using a data-dependent scanning program (DDA). Automatic gain control was set at 5E4, signal threshold as 5000 ions/s, maximum time 200 s, and dynamic exclusion time of the tandem mass scan 15 s, to avoid repeated scanning of precursor ions.

Database Search
MS/MS data were analyzed using the max quant search engine (v.1.5.2.8), with parameters as follows: Rat_Uniprot was first screened, then a reverse library added to calculate the false positive ratio (FPR); Trypsin/P was specified as the cleavage enzyme, allowing up to two missing cleavages; minimum peptide length was seven amino acid residues; maximum number of peptide modifications was five; mass tolerance for the primary precursor ions search was 20 ppm and 5 ppm for the main search; mass tolerance for fragment ions was 0.02 Da; the quantitative method was set to TMT-10plex; and the FPR for peptide-spectrum match identification was 1%.

Bioinformatics Annotation
Gene Ontology (GO) annotations were derived from the UniProt-GOA database (http://www.ebi.ac.uk/GOA/). InterProScansoft was used to elucidate the GO categories: molecular function, cellular component, and biological process. KEGG online service tools were used to annotate protein descriptions, which were matched into corresponding pathways using KEGG mapper. Wolfpsort (https://wolfpsort.hgc.jp/) was used to investigate subcellular localization. KEGG database pathway enrichment analysis was conducted using a two-tailed Fisher's exact test. Pathways were classified according to the KEGG website.
Further cluster analysis of function enrichment was conducted to explore potential connections and differences in specific functions. First, all categories were collected after enrichment, according to their P values, then filtered for those categories with P value < 0.05. This filtered P value matrix was transformed using the function x = −log10 (P value). Finally, these x values were z-transformed for each functional category. Cluster membership was visualized using a heat map.

Protein-Protein Interaction Network
Numbers or sequences of DEPs were compared using the STRING (v.10.5) protein network interaction database. An interaction relationship was extracted according to a confidence score > 0.7 (high confidence).

Statistical Analysis
All data are presented as the mean ± standard deviation (SD). The three groups were analyzed by the Analysis Of Variance (ANOVA) and post hoc Tukey test in the SPSS 21.0 software (SPSS, Armonk, NY, USA). Graphs were generated by using GraphPad Prism 7.0 software.

Electroacupuncture Significantly Reduced Mechanical Hypersensitivity and Fasting Blood Glucose Levels During the Development of DPN.
As shown in the flow chart of the experiment and acupoint diagram ( Figure 1A), DPN model of rats and EA intervention were established and fasting blood glucose and behavioral tests conducted before and after STZ injection and EA intervention. Compared with the control group, body weight and mechanical pain threshold were significantly reduced while fasting blood glucose was dramatically increased at 2 weeks after STZ injection(**P<0.01).
However, compared with the DPN model group, the body weight and mechanical pain threshold were obviously upgraded, while the fasting blood glucose was soothingly downgraded after four weeks of EA treatment(##P<0.01) (Figures 1B-D).

Quantitative of Differentially Expressed Proteins
A total 5744 proteins were identified from spinal dorsal horn samples from the three experimental groups, and 5393 proteins were quantified after proteomic analysis. Molecules with expression ratios showing a > 1.2-fold change and p-values < 0.05 were considered DEPs.
Compared with the control group, levels of 435 proteins were significantly changed (328 down-and 107 up-regulated) in DPN rats and 643 proteins changed in EA rats (p < 0.05 and foldchange > 1.5) (Figure 2A). Further analysis showed that 118 proteins were significantly changed (25 down-and 93 upregulated) in DPN rats while 14 proteins were up-regulated and 83 proteins down-regulated after EA treatment (p < 0.05 and fold-change > 1.5) ( Figure 2B). There were 51 common positive proteins when DPN up-regulated and EA downregulated (p < 0.05 and fold-change > 1.5), while 8 common positive proteins were found which DPN down-regulated and EA-upregulated (p < 0.05 and fold-change > 1.5). There were total 59 common changed proteins shown in the  In terms of biological processes, the top five enriched for DPNassociated proteins were cellular process, single-organism process, biological regulation, metabolic process, and response to stimulus, while the top five enriched for EA-regulated proteins were single-organism process, cellular process, metabolic process, biological regulation, and response to stimulus. Cell component analysis indicated the most enriched for DPN-associated proteins were cell, organelle, membrane, extracellular region, and macromolecular complex, among others, while EA-regulated proteins were mainly involved in the cell, organelle, membrane, macromolecular complex, and extracellular region. Molecular functions analysis suggested that the most common molecular functions of DPNassociated and EA-regulated proteins were binding and catalytic activity ( Figures 3A, B).
According to subcellular structural localization analysis, the main subcategories (> 10%) in DPN-associated proteins were   Figures 3C, D). Last, GO analysis demonstrated that the cellular component, biological process molecular function of DPN and EA treatment were closely related to the mitochondrial respiratory chain and NADH dehydrogenase activity (Figures 3E-J).

Functional Enrichment Analysis With KEGG Pathway and Cluster Analysis of KEGG Pathway Functional Correlation of the DEPs
To further explore pathways significantly enriched for DEPs (DPN vs. Control, n = 643; EA vs. DPN, n = 435), DEPs were mapped to the KEGG database and enrichment levels calculated (Fisher's exact test; -log 10[p-value]), where a higher p-value indicates more DEPs enriched in a category. As shown in Figures 4A, B, KEGG pathway enrichment analysis demonstrated that DPN-associated and EA-regulated DEPs were mainly enriched in several pathways, particularly oxidative phosphorylation ( Figure 4C). Further, remarkable changes were observed in ten proteins closely related to oxidative phosphorylation during EA treatment of DPN rats ( Table 2 and Figure 4D).

DISCUSSION
In our research, we found that four weeks EA can improve weight loss, blood glucose, and sensitivity to mechanical pain in DPN rats. We generated proteomics data from DPN rat spinal dorsal horn samples treated with repeated EA stimulation. A total of 5393 proteins were quantified using the TMT labeling proteomics method. Of these proteins, 93 were down-regulated,  KEGG pathway enrichment analysis suggested that oxidative phosphorylation was a major factor involved in the effects of EA therapy on DPN. The results of this study provide a potentially valuable theoretical underpinning of the mechanisms underlying DPN pathophysiology and the clinical application of acupuncture. The occurrence of diabetes and its complications are closely related to the mechanism of oxidative phosphorylation via the disturbances of mitochondrial ATP, complex IV and reactive oxygen species (ROS) metabolism (22)(23)(24). Upregulated expression of complexes I, III, IV, and V was found in the sciatic nerve of diabetic neuropathy rats (11). In addition, mitochondrial complex III of dorsal root ganglia was significantly decreased, while glutathione peroxidase and peroxidase were increased during STZ induced diabetic peripheral neuropathy (13).
Our results indicate that cytochrome c oxidase (Cox, Complex IV) factors, including Cox4i1, Cox5b, Cox6c2, Cox7b, and Cox7c, were significantly upregulated in spinal dorsal horn during DPN, as were NADH-ubiquinone oxidoreductase (Complex I) factors, including Ndufb5, Ndufb7, andNdufs5. Further, treatment with EA could decrease the up-regulation of Cox and NADH oxidoreductase. As we know, abnormal function of mitochondrial Cox and NADH oxidoreductase could induce ROS overproduction, caspase 3 activation cell apoptosis (25,26). Spinal ROS overproduction and activated caspase-3 induced by Cox play a critical role in microglial activation and neuronal apoptosis in several neuropathic pain models included DPN (27,28). Blocking NADPH oxidase, scavenging ROS, and reducing oxidative phosphorylation by various pharmacotherapy could prevent the induction of spinal long-term potentiation and relieve the persistent pain (29)(30)(31).
In addition, EA can also inhibit oxidative phosphorylation and ROS production through different pathways in neuropathic pain and other chronic diseases. For example, Neuronal ATP receptor P2X7 activation leads to ROS production and subsequent nociceptive pain (32,33). EA can improve neuropathic pain by inhibiting oxidative phosphorylation, reducing ATP production, and suppressing central sensitization in the spinal dorsal horn (34). In middle cerebral artery occlusion rats, EA treatment could significantly improve the respiratory control ratio, and increase the activities of NADH dehydrogenase and Cox in the perifocal penumbra zone (35). In the state of ischemia reperfusion injury, EA can inhibit cardiac Cox and indirectly inhibit the activation of caspase 3 (36). However, EA can alleviate bupivacaine-induced myocardial injury by upregulating enzymatic activity of Cox (37). In a word, EA can biaxially regulate oxidative phosphorylation of Cox and NADH dehydrogenase with bidirectional properties in different diseases and tissues.
According to the above research reports, combined with our results, we speculated that EA can reduce diabetic neuropathic pain by downregulating Cox and NADH dehydrogenase in the spinal cord, which may mainly be through the following ways: (1) inhibiting the spinal ATP release and P2X7 activation, reducing the spinal neuronal abnormal excitability and microglial activation; (2) reducing spinal ROS production and oxidative phosphorylation injury, inhibiting cell injury and even apoptosis. In the next experiment, we will further explore the mechanism of EA regulating ROS and oxidative phosphorylation in DPN.

CONCLUSION
In conclusion, we used a TMT labeling approach coupled with liquid chromatography with tandem mass spectrometry to show that DPN and EA stimulation drive significant changes in spinal proteins, particularly those involved in oxidative phosphorylation. Our results highlight potential candidate protein biomarkers for DPN diagnosis and treatment.

DATA AVAILABILITY STATEMENT
The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found in the article/ supplementary material.

ETHICS STATEMENT
The animal study was reviewed and approved by Fujian University of Traditional Chinese Medicine Animal Ethics Committee.

AUTHOR CONTRIBUTIONS
XY and XC made an equal contribution to this research. ZW and JT are the co-corresponding authors and they completed the project design, data analysis, and drafted part of the manuscript. All authors contributed to the article and approved the submitted version.