Peripheral ERK modulates acupuncture-induced brain neural activity and its functional connectivity

Acupuncture has been widely used as a therapeutic intervention, and the brain network plays a crucial role in its neural mechanism. This study aimed to investigate the acupuncture mechanism from peripheral to central by identifying how the peripheral molecular signals induced by acupuncture affect the brain neural responses and its functional connectivity. We confirmed that peripheral ERK activation by acupuncture plays a role in initiating acupuncture-induced peripheral proteomic changes in mice. The brain neural activities in the neocortex, hippocampus, thalamus, hypothalamus, periaqueductal grey, and nucleus of the solitary tract (Sol) were significantly changed after acupuncture, and these were altered by peripheral MEK/MAPK inhibition. The arcuate nucleus and lateral hypothalamus were the most affected by acupuncture and peripheral MEK/MAPK inhibition. The hypothalamic area was the most contributing brain region in contrast task PLS analysis. Acupuncture provoked extensive changes in brain functional connectivity, and the posterior hypothalamus showed the highest betweenness centrality after acupuncture. After brain hub identification, the Sol and cingulate cortex were selected as hub regions that reflect both degree and betweenness centrality after acupuncture. These results suggest that acupuncture activates brain functional connectivity and that peripheral ERK induced by acupuncture plays a role in initiating brain neural activation and its functional connectivity.

Acupuncture works by inserting and/or manipulating a needle at a specific point called the acupuncture point (acupoint). It produces therapeutic effects in various diseases, including pain-related disease 1,2 and neurodegenerative disease [3][4][5] . Although the clinical effectiveness of acupuncture has been well identified 2,5 , and its fundamental mechanisms are being actively investigated 6,7 , the underlying mechanism of acupuncture action remains unclear.
In determining the mechanism of acupuncture treatment, it is essential to distinguish between the specific and nonspecific effects, also known as the placebo effect of acupuncture. Since the placebo effect is mainly due to psychological mechanisms, it is necessary to study the mechanisms that result from the physical stimuli applied to specific acupoints to distinguish the specific effect of acupuncture 8 . Several brain imaging studies have shown that specific changes in brain neural activation reflect the specific effect of acupuncture [9][10][11] . Recently, it has been known that brain neuronal connectivity, a useful method to clarify the functional network of the brain, plays a crucial role in explaining the neural mechanism of acupuncture 12,13 . It has been used in imaging-based approaches 14 and electroencephalography analysis to detect coordinated activity across distributed and spatially remote brain regions 15,16 . As the demand for more direct analysis of the neural activity in the brain, c-Fos, an immediate early gene, has been used as a surrogate marker for neuronal activity. Increased neuronal activity will produce an up-regulation of c-Fos 17,18 , and subsequent computation of interregional correlations can be analysed to identify a network that is active after acupuncture.
The clinical effect of acupuncture is manifested through several stages of biologic changes. These changes are first initiated in the peripheral skin tissues and/or muscles around the acupoints that were stimulated with Identifying the importance of ERK signaling in the acupuncture mechanism. Previously, we found that peripheral ERK was activated prominently in the early part of stimuli-transmission after acupuncture treatment 22 . Although it can be assumed that the initially expressed ERK plays an important role in initiating the acupuncture effect, the role of other molecules expressed by acupuncture should also be considered. Accordingly, we screened the expression of other molecules after the acupuncture treatment, and in particular, we observed the expression of molecules related to ERK activation. We performed two-dimensional electrophoresis  to derive acupuncture-induced protein candidates in skin tissue. About 500 spots were changed 5 or 30 min after the acupuncture treatment, and 11 spots were remarkably detected. Protein identification of the 11 spots was performed by peptide mass fingerprinting (PMF). Keratin, apolipoprotein A-1, beta-actin, serum albumin, heat shock protein (HSP), transferrin, transketolase, and peroxiredoxin-1 (PRDX1) were topscored or highly scored protein candidates (Supplementary Table 1). Of these proteins, we focused on HSP and PRDX1, a protein associated with not only ERK signaling but also cell signaling transduction. The degree of protein expression was measured at each time point after acupuncture treatment. We found that the activation of phospho-ERK (p-ERK) was increased significantly 5 min after acupuncture treatment compared to the CON group (F 4,10 = 20.04, P < 0.05), and it decreased from 10 min later ( Supplementary Fig. 1A). The activation of phospho-HSP27 (p-HSP27) was also increased 5 min after acupuncture treatment compared to the CON group (F 4,10 = 3.107, P = 0.1394), and it decreased from 10 min later ( Supplementary Fig. 1B). In contrast, the activation of PRDX1 decreased 5 min after acupuncture treatment than that in the CON group (F 4,10 = 0.8736, P = 0.5513), and it was restored gradually ( Supplementary Fig. 1C). Based on these results, we focused on the activation of p-ERK and p-HSP27 after acupuncture as an essential part of initiating the acupuncture effect.
Brain neural activation after acupuncture treatment was regulated by peripheral ERK signaling. To investigate the scientific mechanism of acupuncture from peripheral to central, we first identified the changes in brain neural activity after acupuncture treatment. Then, we identified the correlation between acupuncture-activated peripheral ERK signaling and brain neural activity. Acupuncture was performed on the GB34 acupoints, and the c-Fos positive cells were activated in thirty-four brain regions of the cortex, cerebral nuclei, thalamus, hypothalamus, hippocampus, midbrain, and medulla were analysed. The details of acupuncture treatment and the location of 34 brain regions are shown in Fig. 2. After acupuncture treatment, the number of c-Fos-positive cells was significantly increased in the primary somatosensory cortex (S1; F 2, 20 = 9.008, P < 0.01), secondary motor cortex (M2; F 2, 21 = 18.98, P < 0.001), cingulate cortex area 2 (Cg2; F 2, 21 = 11.34, P < 0.001), insular cortex ( After acupuncture treatment, the number of c-Fos-positive cells was significantly decreased in the CA1 of the hippocampus (CA1; F 2, 22 = 5.004, P < 0.05) and paraventricular thalamic nucleus posterior (PV; F 2, 25 = 7.549, P < 0.05) compared to the CON group, and these changes were inhibited after peripheral U0126 administration (F 2, 22 = 5.004, P < 0.05 in CA1; F 2, 25 = 7.549, P < 0.01 in PV) compared to the ACU group (Fig. 3A,B).
In some brain regions, increased c-Fos activation after acupuncture was not altered by U0126 administration. The number of c-Fos positive cells was significantly increased in the cingulate cortex area 1 (Cg1; F 2, 21 = 7.019, P < 0.01), paraventricular hypothalamic nucleus (PVN; F 2, 18 = 4.250, P < 0.05), raphe magnus nucleus (RMg; F 2, 20 = 6.762, P < 0.01) compared to the CON group. However, these signals was not altered by U0126 administration (F 2, 21 = 7.019, P < 0.05 in Cg1; F 2, 18 = 4.250, P < 0.05 in PVN; F 2, 20 = 6.762, P < 0.05 in RMg when compared between CON vs. ACU + U) ( Table 1). The c-Fos expression in the whole brain regions and the results of one-way ANOVA are shown in Table 1. The activated p-ERK after acupuncture treatment was inhibited by U0126; however, it was not affected by CPZ, a TRPV1 inhibitor and DPCPX, an A1R inhibitor. The activated p-HSP27 were inhibited by U0126, CPZ, and DPCPX (each n = 3). The western blot image represent the cropped blots. Full length blots are presented in Supplementary Figure S2. CON control, CON + U U0126 administration, ACU acupuncture treatment, ACU + U U0126 administration followed by acupuncture treatment, ACU + C CPZ administration followed by acupuncture treatment, ACU + D DPCPX administration followed by acupuncture treatment, N-ACU acupuncture treatment at non-acupoint on the hips. *P < 0.05, ***P < 0.001 compared to the CON group, # P < 0.05, ### P < 0.001 compared to the ACU group, $ P < 0.05, $$ P < 0.01 compared to the ACU + U group. One-way ANOVA was followed by the Newman-Keuls post-hoc test. Data are expressed as the mean ± SEM. www.nature.com/scientificreports/ To identify the brain areas that were most affected by acupuncture and peripheral U0126 administration, the c-Fos changes in the ACU and ACU + U groups were sorted in order of change. LH, ARC, Insul, M2, substantia nigra reticular (SNR), S1, and Cg1 were ranked above the 80th percentile in c-Fos increase after acupuncture treatment compared to the control. ARC, dorsal lateral part of striatum (ST-DL), DMPAG, LH, LPAG, Sol, and posterior hypothalamic area (PH) were ranked above the 80th percentile in the c-Fos decrease of the ACU + U group compared to the ACU group (Fig. 3C,D). Then, the overlap between brain regions ranked above the 80th percentile for acupuncture increased, and U0126 decreased-c-Fos activation was selected. ARC and LH in the hypothalamic area were the most affected brain regions by acupuncture and peripheral U0126 administration, followed by acupuncture (Fig. 3E).
Contrast task PLS analysis after acupuncture. For identification of the brain regions that differ in c-Fos expression between control and acupuncture treatment or between acupuncture and peripheral U0126 administration followed by acupuncture, the contrast task PLS analysis was applied. Then, two pairs of CON vs. ACU and ACU vs. ACU + U were compared. The results showed strongly differentiated patterns in the CON vs. ACU (Fig. 4A) and ACU vs. ACU + U (Fig. 4B), and these two pairs were similar. These distinct patterns were associated with c-Fos activation in multiple regions of the cortex, cerebral nuclei, hippocampus, thalamus, hypothalamus, midbrain, and medulla, and mainly driven by c-Fos activation in the hypothalamus in both comparisons.
Generation of functional connectivity after acupuncture. Next, we computed inter-regional pairwise correlations for c-Fos activation in the 34 brain regions. These analyses revealed that acupuncture produces a higher correlation between the brain regions, mainly in the cortex area, compared to the control. In some Figure 2. Schematic diagram of experimental schedule for functional connectivity study. Acupuncture was performed at bilateral GB34 (Yangneungcheon, Yanglingquan) for 10 min and inhibitor was administrated 15 min before acupuncture treatment. The brains were removed 90 min after acupuncture treatment for identifying the brain neural responses and the brain functional connectivity. The c-Fos positive cells were quantified in each of the 34 brain regions, and a set of inter-regional correlations were computed. Finally, functional networks were generated. Abbreviations for each brain region are shown in Table 1  The overlap between brain regions ranked above the 80th percentile for acupuncture increased, and U0126 decreased-c-Fos activation. CON control, ACU acupuncture treatment, ACU + U U0126 administration followed by acupuncture treatment. Abbreviations for each brain region are shown in Table 1. *P < 0.05, **P < 0.01, ***P < 0.001 compared to the CON group. # P < 0.05, ## P < 0.01, ### P < 0.001 compared to the ACU group. One-way ANOVA was followed by the Newman-Keuls post-hoc test. Data are expressed as the mean ± SEM. The number of animals in each group is shown in Table 1. , and medulla (raphe pallidus nucleus; RPa, RMg, and Sol) also showed high correlation with other brain regions. Peripheral U0126 administration mostly diminished the increased correlation by acupuncture treatment. However, VP in the thalamus, PVN, and anterior hypothalamic area (AH) in the hypothalamus, RPa and RMg in the medulla, showed an increased correlation in the ACU + U group (Fig. 5A). For better optimization of brain functional connectivity, we constructed brain networks by thresholding square of Pearson's correlations in order to consider both positive and negative correlations, and the threshold was 0.288 which minimizes the number of edges when all brain networks have no isolated node (Fig. 5B). Then, Table 1. Abbreviations of brain regions and the number of c-Fos-positive cells in each brain region. Thirtyfour brain regions for c-Fos quantification. Brain regions are categorized by major brain subdivisions. *P < 0.05, **P < 0.01, ***P < 0.001 compared to the CON group. # P < 0.05, ## P < 0.01, ### P < 0.001 compared to the ACU group. CON control, ACU acupuncture treatment, ACU + U U0126 administration followed by acupuncture treatment, n number. One-way ANOVA was followed by the Newman-Keuls post-hoc test. Data are expressed as the mean ± SEM.

No Brain region Abbreviation
Number of c-Fos positive cells
The cortical areas showed a high correlation with other brain regions, including the hypothalamus, thalamus, hippocampus, and cerebral nuclei. In particular, the extensive hypothalamic area (PH, ARC, AH, and LH), Sol of the medulla, CA1 of the hippocampus, and Cg2 of the cortex showed a strong connection with other brain regions. PH, LH, and AH regions of the hypothalamus showed strong connectivity with Sol and RMg of the medulla. In contrast, ARC showed strong connectivity with Cg2 of the cortex. Sol showed strong connectivity with not only the hypothalamus but also the central medial thalamic nucleus (CM) of the thalamus and CA1 of the hippocampus. Additionally, CA1 showed strong connectivity with the SNR of the midbrain, RPa, and Sol of the medulla and S1 of the cortex. The correlation patterns of the ACU + U group were significantly different from those of the ACU group. Most of the increased correlation between the brain regions by acupuncture was diminished by peripheral U0126 administration. However, a correlation between the midbrain and cortical area was observed, and the correlation between the thalamus and cerebral nuclei still existed in the ACU + U group. Additionally, RPa in medulla showed a strong correlation with AH and PVN in the hypothalamus and ST-DM in the cerebral nuclei in ACU + U group (Fig. 5C).
Network clustering in force atlas format. Then, we generated a force atlas format to show the degree and the betweenness centrality between the brain regions. In all results, the nodes were connected by edges, and the weighted number of edges was indicated as a degree. The betweenness centrality for each node is the weighted number of shortest paths that pass through the node from all shortest paths of all possible pairs of nodes. A high degree indicates that one brain region is connected to many other brain regions, and a high betweenness centrality indicates that one brain region is functionally correlated with other brain regions. After acupuncture, many edges between the regions were created. The number of edges was increased in Sol, Insul, S1, M2, secondary somatosensory cortex (S2), CA1, and Cg2 after acupuncture treatment. PH in the hypothalamus showed the highest betweenness centrality after acupuncture treatment, and Cg2, primary motor cortex (M1), ARC, Sol, RMg, and SNR also showed higher betweenness centrality after the acupuncture treatment. The ACU + U group showed an entirely different pattern from the ACU group. The number of edges was increased in RPa, S1, SNR, RMg, Sol, ventral tegmental area (VTA), and PV in ACU + U group, however, these values were lower than those in the acupuncture treated group. Various cortex areas showed high betweenness centrality in the ACU + U group. Among them, S1 showed the highest betweenness centrality (Fig. 6A,B). . Salience scores (right) identified that the brain regions which differed in c-Fos expression between the experimental conditions were intensely similar between the ACU vs. CON and ACU vs. ACU + U. In both comparisons, c-Fos expression in multiple brain regions contributed to the contrast, and the hypothalamus contributed the most strongly to this contrast. CON control, ACU acupuncture treatment, ACU + U U0126 administration followed by acupuncture treatment. Abbreviations for each brain region are shown in Table 1. These results indicate that the Sol is most connected to other brain regions, and the hypothalamic area containing PH and ARC played the most important role in the acupuncture-increased functional connectivity of the brain.
Identification of hub brain regions. Identification of the putative hub regions that work as an important brain region mediating the brain neural mechanism of acupuncture, hub identification was performed. Brain regions were ranked in descending order for the degree and betweenness of all groups. In the ACU group, the regions of Sol, Insul, S1, M2, S2, CA1, and Cg2 were ranked above the 80th percentile for degree, and PH, Cg2, M1, ARC, Sol, RMg, and SNR were ranked above the 80th percentile for betweenness. Then, the overlap of ranked regions was selected. The two putative hub regions of Sol and Cg2 were selected as hub regions changed Figure 5. Matrices of inter-regional correlations and network construction for brain c-Fos expression within each group. (A) Colour-coded matrices showing inter-regional correlations for c-Fos activation between the 34 brain regions. (B) The threshold square of inter-regional correlations for c-Fos activation. Acupuncture produces a high correlation between the brain regions of the cortex, and these correlations were mostly altered by U0126 administration. Red colour indicates high correlation and blue indicates low-correlation (scale, right). (C) A circular layout grouped by major brain subdivision to show the connectivity between the brain regions. Nodes are connected by the edges of super-threshold inter-regional correlations. CON control; ACU acupuncture treatment, ACU + U U0126 administration followed by acupuncture treatment. Abbreviations for each brain region are shown in Table 1 www.nature.com/scientificreports/ after acupuncture treatment. In the ACU + U group, RPa, S1, SNR, RMg, Sol, VTA, and PV were ranked above the 80th percentile for degree, and S1, Pir, Insul, PV, S2, PH, and M1 were ranked above the 80th percentile for betweenness. The two putative hub regions of S1 and PV were selected as hub regions that reacted with peripheral U0126 administration followed by acupuncture treatment (Fig. 6C).

Discussion
Here, we found that peripheral ERK expression by acupuncture plays a triggering role in initiating acupunctureinduced peripheral proteomic changes. And we also found that the brain neural activation and its functional connectivity by acupuncture is regulated by peripheral ERK expression by acupuncture. To the best of our knowledge, this is the first study to explore the brain functional connectivity of acupuncture and to identify the peripheral molecular signaling initiating the central neural activity. Brain regions were ranked in descending order for degree and betweenness centrality of the CON, ACU, and ACU + U groups. Brain regions ranked above the 80th percentile for degree and betweenness centrality were indicated by a red-coloured box. (C) Venn diagram shows the overlap between brain regions ranked above the 80th percentile for degree and betweenness centrality in the ACU and ACU + U groups. Two putative hub regions (Sol and Cg2 in the ACU group; S1 and PV in the ACU + U group) were identified as hub brain regions in each group. CON control, ACU acupuncture treatment, ACU + U U0126 administration followed by acupuncture treatment. Abbreviations for each brain region are shown in Table 1. www.nature.com/scientificreports/ Since acupuncture works by inserting and/or manipulating a needle at a specific acupoint, the peripheral biomolecular changes after acupuncture may play an important role in initiating the therapeutic effect of acupuncture. Previously, several studies reported that peripheral biomolecules such as TRPV1 19,20 , adenosine A1 receptor 21 , or mast cell degranulation 23,24 change after acupuncture stimulation. We previously identified that peripheral ERK signaling, a major molecule activated on the skin tissue after acupuncture treatment, plays a key role in mediating the analgesic effect of acupuncture 22 . In this study, peripheral U0126 administration blocked the expression of acupuncture-induced p-HSP27 activation, while the TRPV1 inhibitor and A1R inhibitor did not affect the acupuncture-induced peripheral ERK expression (Fig. 1D,E). Based on these results, we confirmed that peripheral ERK expression acts as an upstream regulator of peripheral biomolecular expression due to acupuncture treatment.
Analysis of brain functional connectivity is used as a useful method to identify the specific effects of acupuncture. Several studies have reported that acupuncture stimulation at different acupoints induces differential brain functional connectivity in functional MRI (fMRI) 25,26 . Recently, altering functional properties of brain networks using graph theory 27,28 and hub configurations in some brain regions 29 were investigated in advance. These studies have commonly used the degree and mean value to perform hub identification. Despite the similarities between structural and functional networks 30 , it is difficult to detect direct neurophysiological evidence among large-scale brain networks because of their variability in statistical measurements. In an effort to observe more direct brain neuronal responses, the expression of c-Fos activity as an immediate early gene has been utilized to quantify brain neural activity at a single-cell unit 18,[31][32][33] . Previous studies have indicated that acupuncture induced c-Fos activation in some brain regions, including the ARC, ventrolateral periaqueductal grey (VLPAG), Sol, and RMg [34][35][36][37] . However, these studies only showed a partial indication to central nervous system changes, which limits the understanding of the broad spectrum of neural mechanisms regarding acupuncture. Therefore, it is necessary to observe whole-brain level analysis of neural activity and its functional connectivity after acupuncture treatment to improve our knowledge of the mechanisms of acupuncture.
In our results, the number of c-Fos-positive cells in S1, M2, Cg2, Insul, Pir, CA1, PV, LH, ARC, DMPAG, LPAG, and Sol were changed significantly after acupuncture treatment, and these changes were all inhibited by peripheral U0126 administration (Fig. 3). These results indicate that peripheral proteomic changes, especially ERK expression after acupuncture treatment, provoke brain neural changes. When we sorted the c-Fos changes by acupuncture and peripheral U0126 administration followed by acupuncture in order of change, ARC and LH were determined as the most affected brain areas by acupuncture and peripheral U0126 administration followed by acupuncture. Next, we performed salience analysis using the contrast task PLS method to identify the brain regions that are more affected by acupuncture than control or peripheral U0126 administration followed by acupuncture (Fig. 4). Analysis showed that CON vs. ACU and ACU vs. ACU + U showed strongly differentiated patterns in the contrast task PLS analysis, and the two patterns were quite similar. These similar patterns indicate that the brain areas altered by acupuncture, such as the cortex, cerebral nuclei, hippocampus, thalamus, hypothalamus, midbrain, medulla, and especially the hypothalamic area, coincide with the brain areas affected by peripheral U0126 administration. Among these brain areas, the hypothalamus showed the most differentiated pattern. These results are consistent with the previous results of this study, where the brain area most significantly altered by acupuncture and peripheral U0126 administration was the hypothalamus.
Next, we computed inter-regional pairwise correlations for c-Fos activation in the 34 brain regions and constructed brain networks in a circular format to explore the functional connectivity between each brain region altered by acupuncture (Fig. 5). We found that acupuncture produced a higher correlation between the brain regions. The correlation between most cortical areas increased, and the cortical areas had a high correlation with the hypothalamus, thalamus, hippocampus, and cerebral nuclei. Among them, PH, ARC, AH, and LH in the hypothalamus, Sol in the medulla, CA1 in the hippocampus, and Cg2 in the cortex showed a strong correlation with other brain regions. The ACU + U group showed different correlation patterns from that of the ACU group. Most of the acupuncture-increased correlation between the brain regions was diminished by peripheral U0126 administration, while some correlations between the midbrain and cortical area and the thalamus and cerebral nuclei were still observed. Unlike the ACU group, RPa in medulla showed strong correlation with AH and PVN in the hypothalamus and ST-DM in cerebral nuclei. These results suggest that the increased specific functional connectivity between the brain regions by acupuncture may be primarily involved in the clinical effects of acupuncture. Both our previous and current studies show that peripheral U0126 administration blocks the analgesic effects of acupuncture, suggesting that the peripheral biomolecular changes and the central neural activity and its functional connectivity are closely linked. Other brain regions that maintained an increase by peripheral ERK blockade might be associated with the mechanical stimulation of acupuncture rather than its therapeutic effects and is also thought to be related to the non-specific effects of acupuncture.
A force atlas format showed that acupuncture produced prominent functional brain connectivity compared to the control (Fig. 6A,B). The control group had fewer edges between the brain regions and showed a scale-free network. The ACU group had many edges unlike the control group, indicating that the connectivity between the brain regions increased after acupuncture. The extensive area of the cortex was highly connected between the regions, and they were closely connected with some key areas of ARC in the hypothalamus, Sol and RPa in the medulla, and CA1 in the hippocampus. The PH in the hypothalamus showed the highest betweenness centrality after acupuncture. These results are also consistent with our previous results of this study, with the hypothalamus being the most important region mediated by acupuncture. The ACU + U group showed lower values of degree and betweenness centrality than the ACU group, and showed distinct patterns compared to the ACU group in the force atlas format. In the ACU + U group, the high connectivity between cortical areas was diminished, and instead of the hypothalamic area, S1 of the cortex showed the highest betweenness centrality. Thus, these results suggest that acupuncture induces changes in brain functional connectivity patterns, and these changes are mediated by peripheral ERK activation. Additionally, the hypothalamus may be closely related to the www.nature.com/scientificreports/ therapeutic effect of acupuncture, and the S1 region maintained in the ACU + U group is considered to be related only to acupuncture stimulation. These stimuli-related brain regions might be associated with the non-specific effect of acupuncture rather than its therapeutic effects. In our findings, the hypothalamic area has been derived as the major brain region that mediates the acupuncture effect. In brain neural activity analysis, hypothalamic area of ARC, and LH were the most affected brain regions by acupuncture and peripheral U0126 administration. When we analysed distributed patterns between the CON vs. ACU and ACU vs. ACU + U in contrast PLS analysis, the two pairs of patterns were quite similar. The hypothalamic area contributed most significantly to the distinct pattern changes in the brain region in both comparisons. Additionally, PH in the hypothalamus showed the strongest betweenness centrality after acupuncture. From these results, we concluded that the hypothalamus is not only the brain region where brain neuronal activities are most remarkably changed, but also the brain region that is most closely connected to other brain regions.
The hypothalamus, a crucial brain region that modulates the endocrine, autonomic, and central nervous system 38 , is located under the thalamus and anatomically linked to the anterior cingulate cortex, amygdala, hippocampus, and midbrain area. Based on these anatomical connections, the interconnected prefrontal-limbic network is strongly related to depressive disorders 39 , controlling food intake 40,41 , hypertension 42 , and pain 43,44 . Additionally, the hypothalamus is a representative brain region that is activated by acupuncture and is known to be closely associated with the various therapeutic effects of acupuncture, such as pain, chronic stress, and alcohol dependence [45][46][47] . The hypothalamic area is where closely related orexigenic neuropeptide and melaninconcentrating hormone (MCH) has been recently investigated in the role of acupuncture effects. We previously revealed that hypothalamic MCH mediated the neuroprotective action of acupuncture in a Parkinsonian mouse model 48 and also promoted the analgesic effect in a persistent neuropathic mouse model 49 . Chen et al. 50 also reported that median nerve stimulation, which is similar to PC6 acupuncture treatment, activated hypothalamic orexin neurons, and it mediates acupuncture analgesia through a cannabinoid 1 receptor-dependent cascade. Therefore, acupuncture-activated hypothalamic neuronal activation is closely related to the therapeutic effect of acupuncture in various neuronal diseases such as pain or neurodegenerative disease.
Finally, we derived the hub region, representing both high degree and high betweenness centrality. Recent brain network studies using graph theory have typically used degree and mean values to identify hub brain regions [27][28][29] . A degree indicates how many different brain regions are connected to a particular brain region. This can determine the extent of connections between brain regions, but it is difficult to determine how closely each brain region is connected. Since the betweenness centrality is a measure of how closely connected each brain region is 51 , we included not only degrees but also betweenness centrality for hub identification. Hub identification indicated that two regions of Sol and Cg2 have a major role in acupuncture-induced functional connectivity (Fig. 6C). Sol is an important primary nervous center for receiving somatosensory processes in the medulla and has extensive connections with many other structures in the CNS 35,37 . Sol is also the parasympathetic centre of the vagal pathway and transmits various sensory information from the viscera to the brain via afferent vagus fibers. Then, Sol integrates information from various brain areas and sends a signal to regulate the gastrointestinal (GI) tract via efferent vagal nerves 52,53 . Fang et al. demonstrated that electroacupuncture (EA) on ST36 and ST37 induced significant neuronal activity in Sol, which participates in different EA effects by acupoint and frequency selection 52 . Other previous studies also suggested that Sol contributed to the therapeutic effect of acupuncture on inflammation, regulating the GI tract, and relieving pain [54][55][56][57] . The cingulate cortex receives neuronal information from the thalamus and the neocortex, and then projects them to other brain areas of the prefrontal cortex, amygdala, hippocampus, hypothalamus, and midbrain 58,59 . These interconnected prefrontallimbic networks are strongly involved in the formation and processing of emotion, learning, and memory, and these are related to acupuncture signaling transduction, controlling pain, and depression 38,39,60 . The cingulate cortex is known to be related to the mechanism of action of acupuncture in Parkinson's disease, ischemic stroke, depressive disorders, and pain 46,48,[61][62][63][64][65] .
Based on our results, it is thought that acupuncture can contribute to the treatment of numerous diseases by modulating major brain areas such as the hypothalamic area, Sol, and cingulate cortex. Brain regions that were not affected by peripheral U0126 administration were thought to be associated with the non-specific effect of acupuncture. From the results of this study, we suggest a hypothesis of acupuncture-induced peripheral and central mechanisms (Fig. 7). Acupuncture treatment produces peripheral biomolecular changes triggered by ERK activation, and peripheral signaling is received mainly by peripheral nerve endings, which are then transmitted to the afferent sensory nervous tract. Sol and RMg receive mechanical signals from the spinal cord and are functionally transmitted to the hippocampus, thalamus, PAG, ST, and hypothalamus. Then signals are sent to Cg, SC, MC, and Insul. These signals are then sent back to Sol and RMg and might activate the descending inhibitory pathway to provoke the acupuncture effect. The activation of c-Fos-positive cells, which are not affected by peripheral U0126 administration, might be due to the mechanical stimulus of needle insertion.
This study aims to investigate the discursive brain neuronal activity and brain functional connectivity after acupuncture and identify that these brain neural changes were exerted by peripheral molecular signaling. However, the study still has the limitation that no whole-brain areas have been observed and no more samples have been obtained for more rigorous functional connectivity analysis. In addition, there is a limitation that the control group for acupuncture stimulation, such as non-acupoint stimulation, was not used for brain functional connectivity analysis. For further study, it is necessary to include more brain regions to identify the most important brain region within the acupuncture-induced functional network. It is also necessary to investigate the functional network in the disease model with non-acupoint treated controls to determine the central therapeutic mechanism of acupuncture. Additionally, it is necessary to verify whether the brain areas derived from this study actually play an important role in mediating the therapeutic effect of acupuncture. These further studies could elicit a hypothesis of whole-brain functional connectivity related to the therapeutic effect of acupuncture.

Scientific Reports
| (2021) 11:5128 | https://doi.org/10.1038/s41598-021-84273-y www.nature.com/scientificreports/ In conclusion, acupuncture-induced peripheral ERK expression plays a triggering role in producing acupuncture-induced peripheral biomolecular signals. The neuronal activity and functional connectivity were noticeably changed after acupuncture, and these changes were altered by peripheral U0126 administration. The hypothalamic area was the most affected by acupuncture and peripheral U0126 administration. As a result of hub identification, Sol and Cg2 were derived as hub regions that reflect both degree and betweenness centrality after acupuncture. From these results, we have concluded that peripheral ERK expression after acupuncture plays a triggering role in initiating the acupuncture-induced peripheral biomolecular signaling pathway and also plays an important role in mediating brain neural activation and its functional connectivity.

Materials and methods
Experimental animals. Male C57BL/6 mice (8-9 weeks of age, weighing 20-25 g; Samtaco, Seoul, Korea) were used in all experiments. The mice were maintained under a 12 h:12 h light/dark schedule with free access to food and water. The mice were randomly divided, and three to four mice were placed in each cage (room temperature 20-25 °C and humidity of 40-60%) in a clean room. This study was carried out in accordance with the guidelines for the Care and Use of Laboratory  , DPCPX (A1R inhibitor, 3 μg/10 μl; Tocris Bioscience, Bristol, BS110 QL, UK) were used to inhibit molecular activation in the skin layer. Each solution was injected into the GB34 acupoint 15 min before acupuncture treatment according to the previously described 66 . Briefly, mice were slightly anesthetized with ether, then each inhibitor was injected into the GB34 acupoint (i.d., 10 µL/each acupoint) using a 0.5 mL syringe (BD Biosciences, San Jose, CA, USA). The control mice were injected with an equal volume of vehicle.
Acupuncture treatment. Acupuncture was performed at bilateral GB34 (Yangneungcheon) or non-acupoint on the hips for 10 min. GB34, located in the depression anterior and distal to the head of the fibula, is known to have an effect on motor dysfunction and pain (Fig. 1) [67][68][69] . The mice were immobilized in a clear acrylic holder and acupuncture needles (8 mm in length and 0.18 mm in diameter; Haenglim-seoweon Acuneedle Co., Gyeonggi-do, Republic of Korea) were inserted to a depth of 3 mm. Then, needles were turned at a rate of two spins per second for 10 s (one spin consists of 180° clockwise rotation and 180° counter-clockwise Figure 7. Schematic diagram of local molecular signaling and brain function after acupuncture treatment. The peripheral mechanism (right) and central mechanism (left) after acupuncture treatment are illustrated. Acupuncture-induced ERK1/2 activation acts as a triggering molecule that induces multiple signaling transduction around the acupoint. These local signals are received by the peripheral nerve ending and then transmitted to the afferent sensory nervous tract. The brain regions of Sol and Cg are mainly involved in acupuncture-mediated brain functional connectivity, and these brain regions are functionally connected with variable brain regions of SC, MC, HIP, Insul, HyTH, SN, and RMg. These central signals might activate the descending inhibitory pathway that provokes the acupuncture effect. The solid line is the mechanism of the findings in this study, and the dotted line is the potential prediction mechanism. Skin preparation and western blotting. For protein analysis, skin samples around the GB34 acupoints were extracted 5, 10, 30, or 60 min after acupuncture treatment (n = 3 in each group). Briefly, the hairs of legs were removed one day before the experiment, and the mice were anesthetized with ether after acupuncture treatment. The bilateral skin tissues around the GB34 acupoint (5 mm in diameter, consisting of epidermis and dermis) were collected, and were homogenized in a lysis buffer (CyQUANT; Invitrogen, Eugene, OR, USA) containing protease inhibitor and phosphatase inhibitor cocktail tablets. After that, the samples were centrifuged for 15 min at 12,000 rpm, then the supernatants were collected. The protein concentration was calculated using a BCA assay. Western blot analysis was performed according to the previous study 22 . Briefly, 10 µg of total protein were separated using a 10% sodium dodecyl sulfate-polyacrylamide gel electrophoresis, then transferred to a PVDF membrane (Millipore, Billerica, MA, USA). The membranes were blocked with 5% non-fat dry milk in Tris buffered saline for 1 h, then incubated with the primary antibodies such as rabbit p-ERK, total ERK (t-ERK), p-HSP27, total HSP27 (t-HSP27), PRDX1 (Cell Signaling Technology, Beverly, MA, USA) and β-actin (Sigma, St. Louis, MO, USA) overnight at 4 °C. Then, membranes were incubated with secondary antibodies (horseradish peroxidase-conjugated goat anti-rabbit antibody or anti-mouse antibody; Pierce, Rockford, IL, USA) at room temperature. The protein expression was visualized using a chemiluminescence kit (Super Signal West Pico; Pierce, Rockford, IL, USA). 2D gels were run for 1700 Vh then stained with coomassie G250 as previously described 70 . Digitized images were analysed with PDQuest (version 7.0, BioRad) software according to the manufacturer's protocols. Each spot was normalized by total spot intensity, then the protein expression level over two-fold compared with control sample was selected. Protein spots were excised, digested with trypsin (Promega, Madison, WI), and mixed with α-cyano-4-hydroxycinnamic acid in 50% acetonitrile / 0.1% TFA. Then they were subjected to MALDI-TOF analysis (Microflex LRF 20, Bruker Daltonics) for identification by peptide mass fingerprinting as previously described 71 . Spectra were collected from 300 shots per spectrum over m/z range 600-3000 and calibrated with a two-point internal calibration with trypsin auto-digestion peaks (m/z 842.5099, 2211.1046).

Two-dimensional electrophoresis and peptide mass fingerprinting (PMF
The peak list was produced by Flex Analysis 3.0 (threshold: 500 for minimum resolution of monoisotopic mass; 5 for S/N), and the protein identification was performed by peptide mass fingerprinting with the search program MASCOT (http://www.matri xscie nce.com). The parameters were used for the database search as follows: a maximum of one missed cleavage, iodoacetamide (Cys) as a complete modification, trypsin as the cleaving enzyme, monoisotopic masses, oxidation (Met) as a partial modification, and a mass tolerance of ± 0.1 Da.
Immunohistochemistry. Mice were deeply anesthetized 90 min after acupuncture, then perfused transcardially using phosphate buffered saline (PBS) followed by 4% paraformaldehyde (PFA). Brains were removed, post-fixed overnight, then cryoprotected in 30% sucrose for three days. Brains were cut into 40 μm using cryocut microtome (Leica Microsystems Inc., Tokyo, Japan) and prepared for free-floating sections. The brain sections were incubated in 1% hydrogen peroxide for 15 min, incubated with 5% normal goat serum in PBS for 1 h. Then the brain sections were incubated with primary antibody (rabbit c-Fos antibody; Santacruz Biotechnology, CA, USA) overnight at 4 °C, then incubated with secondary antibody (biotinylated anti-rabbit IgG; Vector Laboratories Inc., Burlingame, CA, USA) for 1 h. Then, they were incubated with ABC reagent (Vector Laboratories Inc., Burlingame, CA, USA) followed by 0.02% diaminobenzidine and 0.003% hydrogen peroxide in 1 M Tris-buffered saline (pH 7.5). After that, the sections were dehydrated in ethanol (95% and 100%), incubated in xylene, and mounted. Protein expression was pictured using microscope (BX51; Olympus, Tokyo, Japan). The number of c-Fos positive cells was analysed in 34 brain regions of cortical, cerebral nuclei, hippocampal, thalamic, hypothalamic, midbrain, and medulla regions 72 . The c-Fos positive cells were quantified bilaterally for 2-3 times in each 2.5 mm 2 tissue area using Image Pro Plus 6.0 software (Media Cybernetics, Inc., Silver www.nature.com/scientificreports/ Spring, MD, USA), then a mean value was computed. All procedures were performed blind to minimize the observer bias. To find the brain regions most affected by acupuncture or local ERK inhibition, all controls were converted to 100%, and then the c-Fos increases in ACU groups were arranged in order. In the same way, all ACU groups were converted to 100%, and the c-Fos decreases in ACU + U groups were arranged in order. Then, the overlap between the brain regions ranked above the 80th percentile was selected.
Contrast task PLS analysis (salience analysis). The partial least squares (PLS) analysis, a multivariate statistical technique was used to identify optimal patterns of functional activity that differentiate experimental conditions 73 . Contrast task PLS was used to analyse the brain region activity to describe the relationship between central c-Fos activity and sets of contrasts which represent specific experimental conditions (CON, ACU, and ACU + U). There were three experimental groups in this study, but only two pairs of ACU vs. CON and ACU vs. ACU + U were compared. Brain measurement were stored in an n × m matrix denoted X where the rows are observations and the columns are brain regions. Matrix X is made up of two priori sub-matrices, with n 1 , n 2 being the number of observations in two sub-matrices. The set of contrasts were stored in an n × l matrix Y . The contrasts are pairwise orthogonal and the sum of squares of a given contrast is equal to one. The matrix Y is made up of the same two priori sub-matrices as X, with l contrasts that code for different aspects of the experimental design. Because there are two experimental groups in each pair, the maximum number of orthonormal contrasts is one, i.e. l = 1 . The relationships between X and Y , denoted R, were stored in a cross-product matrix: The singular value decomposition (SVD) of R decomposes it into three matrices: where U, V, and are the normalized left singular vectors, the normalized right singular vectors, and diagonal matrix of the singular values, respectively. The singular vectors U and V are called saliences. The left singular vectors U represent the experimental design that best characterize R , whereas the right singular vectors V represent the brain measurement that best characterized R. In this paper, two saliences representing brain measurement, V of ' ACU vs CON' and ' ACU vs ACU + U' , were compared. These saliences identified the brain regions showing the activation profile across tasks, indicating the maximally expressed brain areas in differences between pairs of experimental groups.
Network generation. Within each of the three experimental groups (CON, ACU, and ACU + U), all possible pairwise correlations between the c-Fos expression in the 34 regions were determined by computing Pearson correlation coefficients. Each complete set of correlations was computed from each vector and were color-coded correlation matrices using MATLAB 2018a (version 9.4, http://www.mathw orks.com/). Networks were constructed by thresholding square of correlations in order to consider both positive and negative correlations. In this procedure, an optimization algorithm was applied to minimize the number of edges when all networks had no isolated nodes. Then, the graph data were visualized using Gephi (version 0.9.2, https ://gephi .org/). Functional connectivity analysis. Two measurements of network in this paper are usually used in brain connectivity analysis, and these are used to detect central brain regions in many previous studies.
The degree represents the number of nodes to which a focal node is connected. It measures the involvement of the node in the network. This is the most basic network measurement, and almost all other measurement are related to degree.
The betweenness centrality is the number of shortest paths that pass through that node from all shortest paths of all possible pairs of nodes 18,74,75,76 . The betweenness centrality of a node ' i' is given as: where 'N' is the number of all nodes, ' ρ st (i) ' is the total number of shortest paths between a source node 's' and a target node 't' that pass through node 'i'. ' ρ st ' is the total number of all shortest paths between the node 's' and the node 't' 77 . All betweenness centralities were normalized.
The degrees and the betweenness centralities of every nodes in each network were measured. To show the connectivity between the brain regions effectively, a circular layout grouped by major brain subdivision was used 78 . And a force atlas format was used to show the degree and the betweenness centrality. Then, we organised degrees and betweenness centralities of each network to distribution tables for hub identification.
Hub identification. Hub regions play important roles in the function of a network. Some previous researchers used the intersection of degree and betweenness centrality to identify hub regions in electrophysiology studies or large-scale brain imaging 77 . In this study, we selected brain regions that are represented as above 80% of both degree and betweenness centrality.

Statistical analysis.
All the data were expressed as the mean ± standard error of the mean (SEM), and analysed using GraphPad Prism 5 software (GraphPad Software Inc., San Diego, CA, USA) for the statistical analysis. For the behaviour test, western blotting, and immunostaining, data were analysed using one-way ANOVA R = Y T X.

Data availability
The datasets used in this manuscript and all analysed data from this study are available from the corresponding author upon reasonable request.