Electroacupuncture Induces Bilateral S1 and ACC Epigenetic Regulation of Genes in a Mouse Model of Neuropathic Pain

Clinical and animal studies have shown that acupuncture may benefit controlling neuropathic pain. However, the underlying molecular mechanisms are poorly understood. In a well-established mouse unilateral tibial nerve injury (TNI) model, we confirmed the efficacy of electroacupuncture (EA) in reducing mechanical allodynia and measured methylation and hydroxy-methylation levels in the primary somatosensory cortex (S1) and anterior cingulate cortex (ACC), two cortical regions critically involved in pain processing. TNI resulted in increased DNA methylation of both the contra- and ipsilateral S1, while EA only reduced contralateral S1 methylation. RNA sequencing of the S1 and ACC identified differentially expressed genes related to energy metabolism, inflammation, synapse function, and neural plasticity and repair. One week of daily EA decreased or increased the majority of up- or downregulated genes, respectively, in both cortical regions. Validations of two greatly regulated genes with immunofluorescent staining revealed an increased expression of gephyrin in the ipsilateral S1 after TNI was decreased by EA; while TNI-induced increases in Tomm20, a biomarker of mitochondria, in the contralateral ACC were further enhanced after EA. We concluded that neuropathic pain is associated with differential epigenetic regulations of gene expression in the ACC and S1 and that the analgesic effect of EA may involve regulating cortical gene expression.


Introduction
Neuropathic pain (NP) is a chronic neurological condition that is extremely challenging to manage and often affects the quality of life in patients [1]. Current pharmacotherapies often utilize anticonvulsants, antidepressants, topical lidocaine, and opioids, but these are not effective for all NP patients and usually come with undesirable side effects [2]. Extensive research is engaged in understanding NP's mechanism and discovering novel targets to inform development of therapeutics with better efficacy and fewer side effects [3].
Electroacupuncture (EA) has been increasingly adopted for pain relief, and many clinical trials and animal experiments have demonstrated its analgesic effect [4][5][6][7][8][9]. It has been shown to induce the release of endogenous opioids, which desensitize peripheral nociceptors and inhibit inflammation [10]. Functional and molecular changes induced by EA have been studied in several pain-related regions such as the dorsal root ganglion (DRG), spinal cord, periaqueductal gray (PAG), amygdala, medial prefrontal cortex (mPFC), anterior cingulate cortex (ACC), and primary somatosensory cortex (S1) [7,[11][12][13][14][15][16], suggesting these nociceptive components are involved in the effect of EA. The fact that an initial peripheral injury can cause lasting peripheral and central sensitization and pain suggests a possible chronical embedment of modifications in the transcription machinery within the nociceptive pathways. Epigenetic alterations refer to changes in gene expression without altering the underlying DNA sequence while being affected by environmental cues and genetic changes [17,18]. DNA methylation is considered to be the most stable epigenetic mechanism and has been identified as an experimentally accessible biomarker candidate among several types of epigenetic markers [18,19]. DNA methylation partakes in the regulation of neuronal intrinsic membrane excitability [20], while restoration of epigenetic features in ion channels has been shown to alleviate NP [21], granting a strong evidence base for the underlying epigenetic mechanism. NP was shown to cause or be associated with epigenetic modulation in the DRG, mPFC, and PAG [22,23] accompanied by robust transcriptomic profile changes involving multiple pathways and biological processes. However, few studies have directly investigated epigenetic regulation and gene expression profiles in S1 and ACC, which are critically involved in the perception of "physical" and "emotional" aspects of pain, the descending modulation of pain, and the pathological development of chronic pain.
To explore the potential alterations in the epigenome and transcriptome in the S1 and ACC under NP conditions as well as the regulation of EA, we analyzed the global methylation and hydroxy-methylation and transcriptomic profiles in the S1 and ACC in a unilateral tibial nerve injury (TNI) model and after EA treatment. Since it is still unclear how acupuncture affects cortical circuits, we examined its effects on cortical S1 regions that are both contralateral and ipsilateral to the injured hindlimb. Our results showed that the TNI resulted in region-and hemisphere-specific changes in DNA methylation and hydromethylation and in transcriptomic profiles, while the effects of EA included both normalizing and enhancing TNI-induced changes in gene expression in these brain regions. We further used immunostaining to confirm significant changes in the expression of gephyrin in S1 and Tomm20 (a mitochondrial marker protein) in the ACC.

Animals
A total of 49 male C57BL/6 mice at the age of 2 months were used in this experiment. These animals were randomly assigned to four groups, namely the Sham (13), Sham-EA (5), TNI (15), and TNI-EA (16) groups. A relatively small number was used for the Sham-EA group because our preliminary results showed a similar pain threshold in the Sham and Sham-EA groups, thus we did not include the Sham-EA group in the histological or RNA sequencing experiments. The animals were housed (no more than 5 mice per cage) in a temperature-and humidity-controlled animal facility on a 12 h light/dark cycle, fed with Teklad Global extruded rodent diet, and had water supplied ad libitum. All procedures were approved by the Animal Care and Use Committee of the Institutional Guide for the Care and Use of Laboratory Animals at the Indiana University School of Medicine (Protocol code: 20072, approval date: 28 September 2020).

Tibial Nerve Injury (TBI) Model of NP
The tibial nerve injury model is well established and widely used to study neuropathic pain. Mice with TNI develop significant mechanical hypersensitivity from 1 week after surgery that persists for several weeks [24]. Briefly, after a mouse was anesthetized with ketamine/xylazine (87.7 mg/kg/12.3 mg/kg, i.p.), an incision was made below the femur on the right hindlimb to expose the biceps femoris muscle. A blunt dissection was then made to expose the trifurcation of the sciatic nerve. The right tibial nerve was cut distally, and the end of the transected tibial nerve was sutured to the anterior surface of the biceps femoris muscle to prevent nerve regeneration [25]. For sham animals, the sciatic nerve was exposed without transection.

Mechanical Nociception Measurement
A mechanical pain threshold test was performed at baseline, day 7, and day 14 after TNI or sham surgery using a simplified up-down (SUDO) method [26] via a series of 5 von Frey filaments. Briefly, mice were placed in acrylic chambers suspended above a wire mesh grid for about 1 h to acclimatize to the environment. During a test, a von Frey filament (North Coast Medical, Inc., Morgan Hill, CA, USA) was pressed against the plantar surface of the hind paw until the filament bent and held for a maximum of 3 s. A positive response was recorded if the mouse moved its foot. We started the tests with the filament number 5 (0.16 g) and progressed accordingly based on an up-down paradigm. We stopped a test if a positive response to the lowest filament or a negative response to the highest possible filament was observed. For each animal, two rounds of tests were performed under each experimental condition to ensure reproducibility. The mechanical threshold was expressed as the paw withdrawal threshold (g).

Electroacupuncture Stimulation
Mice assigned to the EA group were treated with 2/100 Hz EA at the bilateral Yanglingquan (GB34) acupoints located between the heads of the tibia and fibula ( Figure 1A) for 15 min daily from day 7 to day 14 after TNI or sham surgery. This treatment period has been shown to induce a stable antinociceptive effect on neuropathic pain [27][28][29]. Constant current square-wave electrical stimulation generated by a HANS acupoint nerve stimulator (HANS-200A, Nanjing Jisheng Medical Technology Co., Ltd., Nanjing, China) was administered via two acupuncture needles. The intensity of the stimulation was set at 0.1 mA. This intensity can cause noticeable muscle twitching, which is a common criterion used in clinical practice. anterior surface of the biceps femoris muscle to prevent nerve regeneration [25]. For sham animals, the sciatic nerve was exposed without transection.

Mechanical Nociception Measurement
A mechanical pain threshold test was performed at baseline, day 7, and day 14 after TNI or sham surgery using a simplified up-down (SUDO) method [26] via a series of 5 von Frey filaments. Briefly, mice were placed in acrylic chambers suspended above a wire mesh grid for about 1 h to acclimatize to the environment. During a test, a von Frey filament (North Coast Medical, Inc., Morgan Hill, CA, USA) was pressed against the plantar surface of the hind paw until the filament bent and held for a maximum of 3 s. A positive response was recorded if the mouse moved its foot. We started the tests with the filament number 5 (0.16 g) and progressed accordingly based on an up-down paradigm. We stopped a test if a positive response to the lowest filament or a negative response to the highest possible filament was observed. For each animal, two rounds of tests were performed under each experimental condition to ensure reproducibility. The mechanical threshold was expressed as the paw withdrawal threshold (g).

Electroacupuncture Stimulation
Mice assigned to the EA group were treated with 2/100 Hz EA at the bilateral Yanglingquan (GB34) acupoints located between the heads of the tibia and fibula ( Figure  1A) for 15 min daily from day 7 to day 14 after TNI or sham surgery. This treatment period has been shown to induce a stable antinociceptive effect on neuropathic pain [27][28][29]. Constant current square-wave electrical stimulation generated by a HANS acupoint nerve stimulator (HANS-200A, Nanjing Jisheng Medical Technology Co., Ltd., Nanjing, China) was administered via two acupuncture needles. The intensity of the stimulation was set at 0.1 mA. This intensity can cause noticeable muscle twitching, which is a common criterion used in clinical practice. in TNI mice. TNI-induced PWT reduction was improved by EA treatment, while EA did not affect the pain threshold in sham mice (error bars are for standard deviations (same for the following figures); *, ***: p < 0.05, 0.001 compared with TNI group, respectively; two-way ANOVA followed by Bonferroni post hoc tests; n = 5−8).

Tissue Harvest
At day 15 after TNI, half of the mice in each group were anesthetized with ketamine/xylazine (87.7 mg/kg/12.3 mg/kg, i.p.) and transcardially perfused with 30 mL of ice-cold artificial cerebrospinal fluid (ACSF) containing (in mM): 124 NaCl, 2 KCl, 2 MgSO4, 1.24 NaH2PO4, 26 NaHCO3, 2 CaCl2, and 25 dextrose. The animals were then decapitated and the brains were removed. The ACC and S1 regions of both hemispheres were collected, frozen in liquid nitrogen, and stored at −80 °C before RNA extraction. The location of the selected brain regions is illustrated in Figure 1A. in TNI mice. TNI-induced PWT reduction was improved by EA treatment, while EA did not affect the pain threshold in sham mice (error bars are for standard deviations (same for the following figures); *, ***: p < 0.05, 0.001 compared with TNI group, respectively; two-way ANOVA followed by Bonferroni post hoc tests; n = 5-8).

Tissue Harvest
At day 15 after TNI, half of the mice in each group were anesthetized with ketamine/xylazine (87.7 mg/kg/12.3 mg/kg, i.p.) and transcardially perfused with 30 mL of ice-cold artificial cerebrospinal fluid (ACSF) containing (in mM): 124 NaCl, 2 KCl, 2 MgSO 4 , 1.24 NaH 2 PO 4 , 26 NaHCO 3 , 2 CaCl 2 , and 25 dextrose. The animals were then decapitated and the brains were removed. The ACC and S1 regions of both hemispheres were collected, frozen in liquid nitrogen, and stored at −80 • C before RNA extraction. The location of the selected brain regions is illustrated in Figure 1A.
The other half of the animals in each group were similarly anesthetized with ketamine/xylazine (87.7 mg/kg/12.3 mg/kg, i.p.) and then perfused transcardially with 0.9% NaCl solution followed by 4% paraformaldehyde. The brains were removed, postfixed overnight, and transferred to a 30% sucrose solution until they sank. Then, the brains were sectioned at 30 µm thickness using a Leica CM1950 cryostat for immunofluorescent staining.

Fluorescence Microscopy
Fluorescence microscopy images were collected using a high-content imaging system (Molecular Device ImagXpress Micro Confocal) with 10×/NA = 0.45 and 60×/NA = 0.94 Nikon objectives. All collected images were analyzed using MetaXpress (Molecular Device, San Jose, CA, USA). Confocal images were collected using a z-increment of 1 µm and projected via maximal projections. For each animal, at least 2 brain slices containing the ACC and S1 regions were prepared and imaged for the quantification of 5mC, 5hmC, gephyrin, and tomm20. The colocalization of tomm20 and Iba1 was quantified using ImageJ (Image Processing and Analysis in Java; National Institutes of Health, Bethesda, MD, USA) and reported as Mander's overlap coefficients [32].

RNA Extraction
Total RNA was extracted directly from the collected tissues using the PureLink RNA Mini Kit (ThermoFisher, 12183018A) following the manufacturer's protocol. Isolated total RNA concentrations were determined via a Nanodrop Spectrophotometer (Invitrogen) and found to be 40-60 ng/µL. Isolated RNA samples were further analyzed via a Bioanalyzer 2000 to reveal their integrity. RNA samples with an RNA integrity number (RIN) larger than 6.0 were proceeded to sequencing at Novogene (Novogene, Sacramento, CA, USA).

RNA Sequencing (RNA-Seq)
RNA sequencing was performed using a NovaSeq sequencer (Illumina, San Diego, CA, USA) with average clean bases of 7G per sample and an error rate ≤ 0.02%. Mouse genome (GRCm38/mm10) was used as reference for mapping with an average mapping rate > 95%. Differentially expressed genes (DEGs) were calculated using EdgeR. All DEGs with FDR < 0.05 were included in the gene enrichment analysis as indicated in the volcano plot in Figure S2 (Supplementary Materials). Ingenuity Pathway Analysis (IPA) and Gene Ontology (GO) were used to perform the pathway analysis.

Data Analysis and Statistics
Data were expressed as mean ± standard deviation (SD) and analyzed using Graph-Pad Prism 8.3.0 (GraphPad Software, Boston, MA, USA) or OriginPro 2021b software (OriginLab, Northampton, MA, USA).
The values of PWT were compared via two-way ANOVA with Bonferroni post hoc tests. For data in other experiments, all mice samples were collected in at least quadruplets with >2 technical replicates. Statistical difference was determined using a one-way ANOVA followed by Fisher's LSD post hoc test for mean comparisons or Tukey's multiple comparisons test. The significance level was set as p < 0.05.

EA at GB34 Attenuated Allodynia in a Mouse TNI Model
An NP model was established via tibial nerve transection on the right hindlimb in this study. Consistent with previous studies [33,34], animals showed significant allodynia at 7 days after TNI with a lower paw withdrawal threshold (PWT) compared with sham mice, which maintained on day 14 ( Figure 1B, F (2, 30) = 11.57, p < 0.001). Daily EA on bilateral GB34 acupoints for one week resulted in a significant increase in PWT in TNI mice (F (2, 39) = 4.627, p < 0.05), indicating successive EA treatment attenuated TNI-caused allodynia. In contrast, EA stimulation in sham mice did not affect the pain threshold ( Figure 1B, F (2, 24) = 0.1839, p > 0.05).

EA Modulated DNA Methylation and Hydroxy-Methylation in S1 and ACC of TNI Mice
Global changes in DNA methylation and hydroxy-methylation in S1 and ACC of TNI mice were determined by examining the levels of 5mC and 5hmC at 15 days after TNI or sham surgery. Typical images of brain slices stained for 5mC and 5hmC are shown in Figure 2A,B and Figure 2C,D, respectively, with more images available in Figure S1 (Supplementary Materials). The relative changes in 5mC and 5hmC in the ACC and S1 regions were quantified by normalizing to the level of the corresponding sham groups. There was a significant increase (p < 0.05 for both ANOVA and post hoc test) of DNA methylation in the S1 area after TNI in both hemispheres that was partially restored by EA treatment in the contralateral S1 but persisted in the ipsilateral S1 ( Figure 2B,E). The levels of 5mC in the ACC were similar among the three groups ( Figure 2A,E). These results suggested the global DNA methylations were more likely to be altered in S1 by TNI and EA rather than that in ACC.
No significant changes in DNA hydroxy-methylation changes were observed in either the S1 or ACC regions after TNI (see also Figure 2C,D). However, after EA treatment, the level of 5hmC was increased in the contralateral S1 and decreased in the ipsilateral ACC ( Figure 2C,D,F), suggesting regulations of the demethylation mechanism triggered by EA in these brain regions.

EA Regulated Inflammatory and Nervous System Development-Related DEGs in Contralateral S1 of TNI Mice
Since epigenetic changes were observed in the cortex of TNI mice with or without EA treatment, we proceeded to determine transcriptomic alterations in the S1 and ACC via RNA-seq. Overall, our results showed that 32 (9 up and 23 down), 35 (30 up and 5 down), and 33 (26 up and 7 down) differentially expressed genes (DEGs) were identified between the TNI and sham groups in the contralateral and ipsilateral S1 and contralateral ACC, respectively (Figures 3 and 4; see also Figure   and 33 (26 up and 7 down) differentially expressed genes (DEGs) were identified between the TNI and sham groups in the contralateral and ipsilateral S1 and contralateral ACC, respectively (Figures 3 and 4; see also Figure S2 and Tables S1, S3 and S5 (Supplementary Materials)). Similarly, a total of 8 (6 up and 2 down), 24 (2 up and 22 down), and 22 (7 up and 15 down) DEGs were found between the TNI-EA and TNI groups in the contralateral and lateral S1 and contralateral ACC, respectively (Figures 3 and 4; see also Figure S2 and Tables S2, S4 and S6 (Supplementary Materials)).  S1 regions in mice cortex. Venn diagrams of DEGs arising from TNI and EA are summarized in (B) and (D) (top) with an expression map for common genes (bottom) for the contralateral and ipsilateral S1, respectively. (E) Venn diagram showing common genes among different treatment groups in the contralateral and ipsilateral S1. (F) Common DEGs shared between contralateral and ipsilateral S1. N = 4. A total of 32 DEGs were found in the contralateral S1 compared to the sham group (Table S1 in the Supplementary Materials) potentially conferring NP. We then sought to characterize the biological pathways and molecular functions that were significantly enriched by these DEGs. Ingenuity Pathway Analysis (IPA) revealed changes in pathways related to molybdenum cofactor biosynthesis, IL-15 production, human embryonic stem cells, and the role of Oct4 in mammalian embryonic stem cell pluripotency ( Figure 3A), indicating the modulation of neuroplasticity. Gene Ontology (GO) suggested impairments in tissue development with wound healing as the top enriched biological process. The analysis of cellular components and molecular functions further revealed that several genes (e.g., Col5a1 and Fgfr2) involved in extracellular matrix and growth factor binding processes were downregulated ( Figure S3).
To understand the molecular mechanism underlying the analgesic effect of EA on TNI mice, we compared the transcriptomic profiles of the animals in TNI+EA group and TNI only group. Only eight genes were altered as summarized in Table S2 (Supplementary A total of 32 DEGs were found in the contralateral S1 compared to the sham group (Table S1 in the Supplementary Materials) potentially conferring NP. We then sought to characterize the biological pathways and molecular functions that were significantly enriched by these DEGs. Ingenuity Pathway Analysis (IPA) revealed changes in pathways related to molybdenum cofactor biosynthesis, IL-15 production, human embryonic stem cells, and the role of Oct4 in mammalian embryonic stem cell pluripotency ( Figure 3A), indicating the modulation of neuroplasticity. Gene Ontology (GO) suggested impairments in tissue development with wound healing as the top enriched biological process. The analysis of cellular components and molecular functions further revealed that several genes (e.g., Col5a1 and Fgfr2) involved in extracellular matrix and growth factor binding processes were downregulated ( Figure S3).
To understand the molecular mechanism underlying the analgesic effect of EA on TNI mice, we compared the transcriptomic profiles of the animals in TNI+EA group and TNI only group. Only eight genes were altered as summarized in Table S2 (Supplementary Materials). These genes were particularly enriched in inflammatory-response-related signaling pathways such as the FXR/RXR activation, LXR/RXR activation, RhoA signaling, and acute phase response signaling pathways ( Figure 3A). GO analysis indicated these DEGs were mainly enriched in nervous system development-related pathways and the ion channel complex ( Figure S4). Two genes were identified as shared between the comparisons of TNI vs. sham and that of TNI-EA vs. TNI ( Figure 3B, Venn diagram), both of which were non-coding RNA. One was 9330111N05Rik, which was downregulated in TNI mice but reversed by EA treatment (Figure 3B, heatmap). The other was RN7sk, which was downregulated by TNI and further reduced by EA treatment. RN7sk is a long non-coding RNA (lncRNA) that is part of the 7SK nuclear ribonucleoprotein; it is upregulated during neural differentiation, suggesting its possible involvement in the acquisition of a neural fate [35].
3.3.2. EA Regulated Biosynthesis, Inflammation, and Glycinergic Synapse-Related Pathways in Ipsilateral S1 of TNI Mice The transcriptomic profile was also characterized in the ipsilateral S1 region, which has also been shown to be reorganized in peripheral nerve injury models along with the contralateral S1 [36]. A total of 35 DEGs were identified in the ipsilateral S1 after TNI (Table S3 in the Supplementary Materials); the enriched canonical pathways included molybdenum cofactor biosynthesis and IL-15 production, which were also enriched by the DEGs in contralateral S1, as well as the kinetochore metaphase signaling pathway, STAT3 pathway, and tRNA charging pathway ( Figure 3C). The GO analysis suggested that these DEGs were associated with muscle proliferation and growth, glycinergic synapse and inhibitory synapse, as well as protein tyrosine kinase activity ( Figure S5 in the Supplementary Materials).
After EA treatment, 24 genes were differentially expressed in the ipsilateral S1 of the TNI mice. They were mostly enriched in the same canonical pathways as the TNI-induced DEGs: molybdenum cofactor biosynthesis, IL-15 production, the STAT3 pathway, and tRNA charging. Moreover, the role of Oct4 in mammalian embryonic stem cell pluripotency was also noted ( Figure 3C). Similarly, the GO analysis indicated DEGs were enriched in muscle proliferation, glycinergic synapse and inhibitory synapse, and protein tyrosine kinase activity ( Figure S6 in the Supplementary Materials).
Interestingly, among these 24 DEGs, 19 of them were upregulated in the ipsilateral S1 of TNI mice and reversed after EA treatment ( Figure 3D), indicating these abnormally expressed genes may be involved in the pain mechanism and associated with the analgesic effect of EA. The pathway analysis suggested that these genes are particularly enriched in pathways related to transcriptional regulation such as stem cell pluripotency as well as cellular growth and development, such as the STAT3 pathway ( Figure 3C). It is worth noting that Gphn, which encodes gephyrin, a neuronal assembly protein that anchors inhibitory neurotransmitter receptors to postsynaptic cytoskeleton, was increased in the ipsilateral S1 of TNI mice then downregulated by EA, implying a trend of enhanced synaptic inhibition in this region after TNI and its reverse by EA treatment.

Different Patterns of Transcriptomic Profiles between Interhemispheric S1 after TNI
We further explored the overlapping of DEGs among groups ( Figure 3E). There was one gene (Rn7sk) that was shared among all groups. Ten genes were shared between both side of S1 after injury, seven of which were downregulated in the contralateral S1 but upregulated in the ipsilateral S1 ( Figure 3F), indicating a possible compensation of the ipsilateral S1 in response to TNI. The three genes that were increased in both S1 regions were all pseudogenes. Among them, two do not have well-established functions and Tma7-ps is known as a pseudogene for translational machinery associated 7 homologs.
When comparing the EA effects between the ipsilateral and contralateral S1, there were no common genes except Rn7sk that were shared among all groups, suggesting different mechanisms were involved between hemispheres.

EA Regulated Biosynthesis, Inflammation, and Glycinergic-Synapse-Related Pathways in Contralateral ACC of TNI Mice
The ACC is an important brain region in the perception of emotional and motivational components of pain [37,38]. In the contralateral ACC, 33 DEGs were observed between TNI and sham, and 22 were observed between TNI+EA and TNI. Figure 4A summarizes the top five canonical pathways that were primarily regulated by TNI and EA. The regulations of these pathways collectively suggested an altered metabolic state with changes enriched inside mitochondria, which was consistent with our GO analysis (see Figures S7 and S8  (Supplementary Materials)). Comparing DEGs identified from the TNI vs. sham and TNI-EA vs. TNI groups showed that 10 DEGs were shared between the two comparisons ( Figure 4B), all of which were oppositely regulated by TNI and EA treatment. Among them, Col6a1 and Col6a4, which encode the alpha 1 and alpha 4 chain of collagen type VI, were reduced in TNI mice and restored after EA, while the other eight DEGs were elevated in TNI mice but reduced after EA. Among them, mt-Co2 and mt-Atp6 are both mitochondrial protein-coding genes that participate in the oxidation-reduction process, ATP synthesis, and electron transport. These mitochondrial protein-coding genes (such as Mmp12, Cd36, and Gpnmb) were found in the oxidative phosphorylation and mitochondria dysfunction pathway as we showed previously in Figure 4A. The opposite effects of TNI and EA treatment on these shared DEGs implied potential analgesic mechanisms of EA on TNI mice via altering mitochondrial functions.
To examine the region specificity in the effects of TNI and EA treatment, we also compared the DEGs in the contralateral S1 and ACC between the TNI and sham groups as well as the TNI-EA and TNI groups ( Figure 4C). In the comparison of the TNI and sham groups, there was no shared DEG in the contralateral ACC and S1, and only one gene (Crocc) was altered in both the contralateral S1 and ACC after EA treatment, suggesting a strong region-specific response to TNI and EA treatment.

EA Reversed the Augmentation of Gephyrin Expression in Ipsilateral S1 of TNI Mice
Based on the transcriptomic changes in the gephyrin gene in the TNI and TNI-EA groups, we further validated gephyrin protein expression via immunofluorescent staining in brain slices of TNI mice. MAP2 was used as a counter-stain to facilitate the identification of neurites. Typical confocal images of gephyrin staining in the contralateral and ipsilateral S1 are shown in Figure 5A. The density of gephyrin puncta was slightly increased in contralateral S1, but this increase became significant after EA treatment ( Figure 5A,B). In contrast, the density of gephyrin puncta was significantly elevated in the ipsilateral S1 of TNI mice, but it was reversed by EA treatment (Figure 5A,C).

The Mitochondrial Intensity in Neurons Wase Elevated by EA Treatment in ACC of TNI Mice
The overall mitochondrial intensity per cell was increased after TNI and further increased by EA treatment on both the contralateral and ipsilateral ACC ( Figure 6A,B). The trend in overall mitochondrial amount was quantified as the overall mitochondrial density normalized to cell density for the contralateral and ipsilateral regions, respectively ( Figure 6C,D). To further determine cell-type-specific mitochondrial alterations, we colocalized Tomm20 with MAP2 or Iba1 to identity the mitochondria in neurons or microglia ( Figure 6A,B and Figure 7A,B). Although there were no significant changes in the staining of the Tomm20 protein in neurons of both the contralateral and ipsilateral ACC in TNI mice, their expressions were significantly higher after EA treatment ( Figure 6C). On the contrary, Tomm20 expression in microglia was significantly elevated in the contralateral ACC and sustained after EA treatment ( Figure 7C), but its expressions in microglia of the ipsilateral ACC were not altered in either the TNI or TNI+EA conditions ( Figure 7D). These findings suggested that increased mitochondria protein expressions in neurons and microglia may underly the effect of EA.

The Mitochondrial Intensity in Neurons Wase Elevated by EA Treatment in ACC of TNI Mice
The overall mitochondrial intensity per cell was increased after TNI and further increased by EA treatment on both the contralateral and ipsilateral ACC ( Figure 6A,B). The trend in overall mitochondrial amount was quantified as the overall mitochondrial density normalized to cell density for the contralateral and ipsilateral regions, respectively ( Figure 6C,D). To further determine cell-type-specific mitochondrial alterations, we colocalized Tomm20 with MAP2 or Iba1 to identity the mitochondria in neurons or microglia ( Figures 6A,B and 7A,B). Although there were no significant changes in the staining of the Tomm20 protein in neurons of both the contralateral and ipsilateral ACC in TNI mice, their expressions were significantly higher after EA treatment ( Figure 6C). On the contrary, Tomm20 expression in microglia was significantly elevated in the contralateral ACC and sustained after EA treatment ( Figure 7C), but its expressions in microglia of the ipsilateral ACC were not altered in either the TNI or TNI+EA conditions ( Figure 7D). These findings suggested that increased mitochondria protein expressions in neurons and microglia may underly the effect of EA.

Discussion
In this study, we evaluated epigenetic and transcriptomic alterations in the S1 and ACC in a mouse TNI model of NP and the effect of EA. We identified several DEGs that were regulated by TNI and EA treatment, analyzed the involved signaling pathways, and validated the expressions of gephyrin and mitochondria proteins. Our results revealed the epigenetic regulation of genes in the bilateral S1 and contralateral ACC by TNI-caused neuropathic pain and after EA treatment. Figure 7. Expression of Tomm20 in microglia in contralateral (A) and ipsilateral (B) ACC in sham, TNI, and TNI-EA mice. Colocalization of Iba1 and Tomm20 was increased in the contralateral ACC of TNI mice, and EA treatment did not alleviate this alteration (C). The expression of Tomm20 in microglia was similar in the ipsilateral ACC of sham, TNI, and TNI-EA mice (D). Scale bars in panel A: 50 µm (left) and 20 µm (right). *, **: p < 0.05, 0.01 compared with Sham group; one-way ANOVA followed by Tukey's multiple comparisons test; N = 7-8.

Discussion
In this study, we evaluated epigenetic and transcriptomic alterations in the S1 and ACC in a mouse TNI model of NP and the effect of EA. We identified several DEGs that were regulated by TNI and EA treatment, analyzed the involved signaling pathways, and validated the expressions of gephyrin and mitochondria proteins. Our results revealed the epigenetic regulation of genes in the bilateral S1 and contralateral ACC by TNI-caused neuropathic pain and after EA treatment.

Brain-Region-Specific Alteraions in DNA Methylation and Hydroxy-Methylation by TNI and EA
DNA methylation has been demonstrated to be involved in various comorbid conditions in neuropathic pain, and alterations in DNA methylation have been reported in both spinal and supraspinal regions following persistent pain states. Spared nerve injury was accompanied by a global DNA hypermethylation status in the DRG and spinal dorsal horn (SDH) in rats [39]. It decreases global methylation in the PFC and amygdala but not in the visual cortex or thalamus in mice [40]. Here, we found that TNI caused increased DNA methylation in the bilateral S1 but not the ACC (Figure 2A,B,E). These findings suggested that changes in global DNA methylation in neuropathic pain models may be region-specific. Although there is a growing body of studies exploring the molecular mechanisms underlying the analgesic effects of EA, the role of DNA methylation is still unclear. Jang et al. found that acupuncture reversed the decrease in global DNA methylation in the PFC in a partial sciatic nerve ligation model and restored the reduction in DNA methylation of Nr4a1, Chkb, and Ras pathway-related genes [41]. Li et al. found that EA ameliorated hypermethylation of hippocampal neural stem cells in mice with peripheral nerve injury by restoring inhibited Tet1, which contributed to an antidepressant effect of EA [42]. In our study, EA treatment partially restored increased DNA methylation in the contralateral S1 in TBI mice while increasing the level of DNA hydroxy-methylation in the contralateral S1 and ipsilateral ACC ( Figure 2). Together, these findings indicated the effect of EA on DNA methylation/hydroxy-methylation in mediating the perceptive, sensory, cognitive, and neurogenic impacts of neuropathic pain in a region-specific manner.

Altered Transcriptional Profiles in S1 and ACC after TNI
The S1 is a brain region critical in the processing of pain signaling that has been reported to show functional or structural changes in chronic pain conditions [43]. The S1 of chronic pain patients was actively displayed by fMRI [44]. Our previous study also showed increased neuronal activity in S1 of NP mice [45]. Imaging studies further suggested that remodeling of the S1 synaptic structure following a peripheral nerve injury contributed to the development of NP [46,47]. Increased axonal bouton turnover was observed in the S1 cortex of mice with partial sciatic nerve injury during the early development phase of neuropathic pain [47]. In addition to the contralateral S1, the ipsilateral S1 also experiences functional changes in NP [36]. Our RNA-seq results indicated that the DEGs of the contralateral S1 were mainly enriched in neuroplasticity, inflammatory responses, wound healing, extracellular, and growth factor binding-related pathways or biological processes. The DEGs determined in the ipsilateral S1 were enriched for biosynthesis, inflammatory response, the kinetochore metaphase signaling pathway, muscle proliferation and growth, and inhibitory synapse and protein tyrosine kinase activity. These findings were consistent with previous studies regarding the contribution of S1 in NP. Ten DEGs were shared in the ipsi-and contralateral S1, seven of which were reduced in the contralateral while elevated in the ipsilateral S1 (Fgfr2, Rn7sk, Jarid2, Zc3h7a, CT010467.1, Gphn, and Rpph1). Among them, Fgfr2 and Gphn have been demonstrated to be downregulated in the DRG of rats with chronic pain [48][49][50]. Fgfr2 encodes receptor 2 of fibroblast growth factors, which are key players in peripheral nerve regeneration after injury [48]. Gphn encodes gephyrin, an assemble protein that clusters with glycine or the GABA receptor to form inhibitory postsynaptic density [51]. A loss in gephyrin may reflect disinhibition in the contralateral S1 under a chronic pain condition. Disruption of Gphn expression has been shown to alter postsynaptic clustering of gephyrin, leading to larger postsynaptic density (PSD)-95 clusters and larger presynaptic glutamatergic terminals [52] that collectively contribute to pain hypersensitivity. The restoration of these genes in the ipsilateral S1 might be a compensatory effect in response to the reductions in Gphn in the contralateral S1. The other three genes (Tma7-ps, Gm15459, and Gm6989) that were increased in both the ipsi-and contralateral S1 were pseudogenes with unknown functions.
The ACC is another important brain region in the perception of emotional and motivational components of pain. Although there was no global change in the levels of 5mC or 5hmC in the ACC, changes in DNA methylation may occur at individual gene levels but still be indetectable by global methylation markers. Thirty-three DEGs were found in contralateral ACC, indicating that they may be regulated by other mechanisms. The detected DEGs were mainly involved in mitochondrial function, oxidative phosphorylation signaling pathways, and inflammatory responses ( Figure 4A). Similarly, an alteration in genes was also found to be enriched in inflammation-related pathways in the ACC of a rat spinal nerve ligation (SNL) model [53].

EA Regulated Different Pathways in S1 and ACC for Analgesia
EA is a promising alternative treatment for NP and has been increasingly practiced clinically over the past decades. Its efficacy has been investigated in both patients and animal models. For example, eight weeks of EA treatment can remarkably decrease pain intensity in NP patients and be accompanied by noticeable psychosocial and functional improvement [54]. Application of 2 Hz EA has been shown to significantly improve the mechanical and thermal hyperalgesia in rats or mice with chronic constriction injury of the sciatic nerve. [55,56]. Application of 100 Hz EA has also been shown to relieve chronic pain arising from the spared nerve injury in rats [57].
EA exerts its analgesic effects by regulating the neuronal activity in the S1. For example, fibromyalgia patients receiving EA experienced great reductions in pain severity and showed increased S1 leg to anterior insula connectivity contributing to pain relief [58]. Several other labs have reported that EA stimulation can reduce the activity in the S1 cortex of rats with peripheral nerve injury [16]. Furthermore, the ipsilateral and contralateral S1 seems to have different responses to EA in pain patients. EA treatments in patients with carpal tunnel syndrome and median nerve injury can activate neuron activity in the ipsilateral S1 while deactivating the contralateral S1 [59]. However, no studies have investigated how S1 contributes to the analgesic effect of EA. Here, we found that EA exerted its regulatory effects differently on the ipsi-and contralateral S1. Specifically, in the contralateral S1, the DEGs induced by TNI were mostly enriched in plasticity-related pathways, while EA mainly regulated inflammation-related pathways ( Figure 3A,C). Only two genes were shared between the two comparisons ( Figure 3B). Specifically, 9330111N05Rik was reduced in TNI mice and restored by EA, while 9330111N05Rik was found to be upstream to the Nr2f1 gene that is known to be proneurogenic and contribute to neurological diseases [60]. The other is RN7sk, which was downregulated by TNI and further reduced by EA treatment. RN7sk is a long non-coding RNA (lncRNA) that is part of the 7SK nuclear ribonucleoprotein. RN7SK is typically overexpressed in the cortex and mediates global transcriptional regulation. A strong upregulation of RN7SK was observed during neuronal differentiation both in vivo and in vitro [35]. During the reprogramming process of converting fibroblasts to neurons, RN7SK was shown to increase chromatin accessibility and facilitate neuronal gene activation [61].
Different from the contralateral S1, 19 genes were found to be upregulated by TNI and reversed by EA in the ipsilateral S1. This unexpected finding suggested a possible role of the ipsilateral S1 in the analgesic effects of EA ( Figure 3D). A canonical pathway analysis suggested that these DEGs were involved in molybdenum cofactor biosynthesis, IL-15 production, the STAT3 pathway, and tRNA charging, most of which have been extensively demonstrated to be associated with chronic pain conditions. T-cell cytokine IL-15 was increased in both spinal cord and lesioned nerves in a chronic constriction injury (CCI) model and could last for several weeks [62,63]. Blocking IL-15 activity can reduce macrophage and T-cell infiltration, suggesting that IL-15 is essential to the proinflammatory state in the spinal cord after a peripheral lesion that generates NP [63]. Additionally, active and phosphorylated STAT3 was shown to accumulate in dorsal spinal cord of rats with spinal nerve ligation, while blocking the STAT3 pathway could attenuate both mechanical allodynia and thermal hyperalgesia [64]. These studies supported our findings that EA may facilitate the normalization of these aberrant pathways in the contralateral S1 to induce an analgesic effect.
Prior studies have shown that regulation of EA on the ACC in NP models is associated with the phosphorylated ERK, CREB, 5-HT, and BDNF signaling pathways [15,65]. In our study, the oxidative phosphorylation, GP6, HOTAIR, mitochondrial dysfunction, and HIF1a pathways were regulated by EA ( Figure 4A). Ten DEGs were oppositely regulated by TNI and EA. Among them, Col6a1 and Col6a4, which encode the alpha 1 and alpha 4 chain of collagen type VI, were found to be reduced in TNI mice and restored after EA. Studies found that a lack of collagen VI caused functional defects of peripheral nerves and delayed response to pain stimuli [65], while the other eight genes were all increased after TNI and reversed by EA ( Figure 4B). Gpnmb encodes glycoprotein nonmetastatic melanoma B (GPNMB), which plays a pivotal role in neuronal survival and neuroprotection. GPNMB was increased in spinal cord after sciatic nerve injury in rats, and inhibiting GPNMB with siRNA alleviated mechanical allodynia and thermal hyperalgesia and also reduced proinflammatory factors such as TNF-α, IL-1β, and IL-6 [66]. MMP-12, which encodes matrix metalloproteinase 12, is thought to maintain peripheral neuroinflammation. It is elevated in pain conditions along with concomitant macrophage infiltration, demyelination, and elastin fiber loss at the lesion site [67]. Mitochondrially encoded cytochrome C oxidase II and mitochondrially ATP synthase membrane subunit 6 (mt-Co2 and mt-Atp6, respectively) were both increased in TNI mice but reversed by EA.
Taken together, we found different regulation patterns of EA in the S1 and ACC and observed a remarkable number of overlapping DEGs that can serve as potential candidates that contribute to the analgesic effects of EA on NP.

EA Regulates Gephyrin and Mitochondrial Function
Our immunofluorescent staining verified an increased level of gephyrin in the ipsilateral S1 that was reversed by EA treatment. The expression of Gphn can be epigenetically regulated. Although gephyrin is shown to be involved in the homeostatic regulation of inhibitory synaptic strength [68,69], its role in NP is unknown. Our results thus suggested the involvement of gephyrin in mediating NP and EA effects, and we speculated that the epigenetic regulation of Gphn may contribute to the observed transcriptional changes.
Mitochondria are the energy-production center of cells. Abundant and connected mitochondria are typically indicative of active and functioning cells, while scarce and fragmented mitochondria are commonly associated with dormant or repairing cells. Moreover, mitochondria are known to contribute to neuron functions such as excitability and action-potential firing [70,71]. Mitochondrial dysfunction has been shown to contribute to peripheral neuropathy induced by chemotherapy or various disease conditions such as diabetes and HIV [72]. However, few studies exist that directly assessed the correlation between mitochondrial functions and NP. A recent work by Dai et al. suggested that a Drp1-mitochondrial fission-reactive oxygen species (ROS) cycle contributes to NP by modulating the phosphorylation of Drp1 [73] and calcium homeostasis [74]. In addition, traumatic peripheral nerve injury can induce persistent mitochondrial and bioenergetic dysfunction, supporting their possible involvement in pain pathogenesis. Pharmacological agents that normalize or promote mitochondrial function were hypothesized to be beneficial for pain treatment [75]. Pertinent to this work, we observed consistent upregulation of mitochondria content in neurons and microglia after EA that may contribute to meeting the increased energy needs after TNI and help to normalize neuron functions [75].

Conclusions
In summary, we found that one week of EA treatment at GB34 significantly alleviated the mechanical hyperalgesia in a TNI model of NP. TNI resulted in higher DNA methylation of both the contra-and ipsilateral S1, while EA reduced contralateral S1 methylation without affecting that of the ipsilateral S1. We identified several differentially expressed genes in the ipsilateral and contralateral S1 as well as the contralateral ACC that were related to energy metabolism, inflammation, neural plasticity, and repair. EA normalized the majority of the up-or downregulated genes in both the S1 and ACC regions. Particularly, an increased expression of gephyrin in the ipsilateral S1 but not the contralateral S1 after TNI was normalized by EA. In contrast, mitochondria were further enhanced by EA treatment in both the contralateral and ipsilateral ACC. We concluded that NP is associated with differential regulations of gene expression in the ACC and S1 and that the analgesic effect of EA may be mediated by both normalizing and compensatory changes in altered cortical gene expression. These findings provide novel insights toward understanding the underlying cortical molecular mechanisms of the analgesic effect of EA on NP.
Author Contributions: X.J. and C.Y. participated in the research design and conceptualization; X.P. and J.X. conducted the experiments and performed the data analysis; all authors contributed to writing the manuscript. All authors have read and agreed to the published version of the manuscript.

Informed Consent Statement: Not applicable.
Data Availability Statement: The datasets generated and/or analyzed in this study are available from the corresponding author upon reasonable request.