A total of 40 Institute of Cancer Research (ICR) female mice (weight: 35 ± 5 g) were obtained and raised at the experimental animal center of the Xiamen University. The animal study was reviewed and approved by the Xiamen University’s Animal Ethics Committee (Permit Number: SCXK2018-0003). All procedures were conducted in accordance with the regulations of the “International Council for Laboratory Animal Science”. The experimental mice were randomly assigned to the control, POF, electro-acupuncture at the acupoints (EA), and electro-acupuncture at the non-acupoints (EN) groups (n = 10 in each group). All mice were fed and drank freely during the experimental process. The POF mice were established via intraperitoneal injection of cisplatin (2 mg/kg) daily for 2 weeks, except for the control group (12).

CV4 and SP6, as the acupoints, were selected in the EA group. SP6 was located at 0.5 cm above the medial malleolus of the hind limb, whereas CV4 was located 1 cm below the navel (the navel was located at the lower third of the line between the xiphoid process and perineum). Correspondingly, the non-acupoints were at 3-mm horizontal distance to CV4 and 3 mm higher than SP6 ( Figure 1A ). Generally, the non-acupoints did not belong to any known meridians. In the process of treatment, a breathable and opaque headgear was used to completely cover the head of the mice so as to help the mice remain calm during the treatment. Both acupoints and non-acupoints were inserted at 0.8 mm with acupuncture needles. The frequency of the EA instrument was adjusted until a slight beating of the skin was observed. Stainless steel acupuncture needle (0.25 mm × 13 mm; Suzhou Medical Supplies Factory Co., Ltd, Suzhou, China) was used ( Figures 1B, C ). All mice in the EA and EN groups were treated for 30 min, daily, for 3 weeks.

Into the vagina of the mice, 20 μl of 0.9% NaCl solution was dripped from 09:00 AM to 10:00 AM daily. Then, 0.9% of the NaCl solution was mixed gently and repeatedly dropped in the mice vagina with a pipette gun. Next, the solution was sucked out with a pipette gun and dropped on a slide. Nucleated epithelial, cornified epithelial, and leukocytes of mice in each group were observed under a light microscope to determine the estrous cycle.

All mice were anaesthetized with 10% chloral hydrate after treatment. One side of the ovary from the mice was collected surgically and placed in 4% paraformaldehyde. Next, these dehydrated ovaries were embedded in paraffin and sectioned into 5-μm-thick slices with a freezing microtome (CM1950, Leica Biosystems Division of Leica Microsystems Inc., Germany). Subsequently, the slices were dehydrated in alcohol and stained with hematoxylin–eosin (HE). The pathological morphology of the ovary was observed under the light microscope.

The apoptosis of granulosa cells in the ovary was detected by terminal deoxynucleotidyl transferase dUTP nick-end labeling (TUNEL). All processes were conducted as per the instructions of the TUNEL Assay kit (G1501, Wuhan Servicebio Technology Co., Ltd., China). The apoptotic granulosa cells were labeled green with a fluorescent reagent and observed under the fluorescence microscope. To avoid any pathological differences between the samples, three ovarian sections were prepared for each sample. The percentage of TUNEL-positive granulosa cells was calculated and analyzed by the Image-Pro Plus 6.0 software.

Blood samples were collected from the orbital artery and sacrificed after treatment. FSH, luteinizing hormone (LH), estradiol (E₂), and anti-Mullerian hormone (AMH) in the serum were detected by using the enzyme-linked immunosorbent assay (ELISA) kit (0555M2/44039M2/0546M1/44204M2, Jiangsu Meimian Industrial Co., Ltd, China). All procedures were conducted as per the manufacturer’s instructions.

The relative expression of estrogen receptor-α (ER-α), receptor-β (ER-β), and G protein-coupled ER (GPR30) in the ovary were detected by quantitative real-time (qPCR). Total RNA was extracted from the ovary by the Trizol method and reversed transcribed into cDNA. Subsequently, cDNA was used as the template for amplification. All data were analyzed by the software of Quantity One, and the relative expressions were calculated according to the 2^−△△CT method.

The liver and kidney were removed by surgery after all the mice were sacrificed. The metabolites of the liver and kidney were performed with the ¹HNMR spectrometer (Bruker AVANCE-III 600MHz, Switzerland) at 298 K. The samples of the liver and kidney were weighed (300 mg), and the homogenate was mixed with 600 μl of methanol and 300 μl of double-distilled water. All the samples were placed and allowed to stand for 10 min on an ice bath. Next, the samples were centrifuged (10,000 rpm/10 min at 4°C) and dried with nitrogen. Then, 600 μl of D2O was mixed, and 500 μl of the solution was extracted into a nuclear magnetic tube. For all samples, 64 FIDs were collected into 64K data points over a spectral width of 12,000 Hz with a relaxation delay of 6.5 μs.

MestReNova software (Mestrelab Research, Santiago de Compostela, Spain) was used to calibrate and optimize the ¹HNMR spectra. Specifically, the ¹HNMR spectra of each sample were processed through baseline calibration, phase correction, and water peak removal. To accurately analyze the ¹HNMR spectrum, the chemical shift interval (0.5–9) was selected for piecewise integration, and the resulting data were imported into the SIMCA-P software (Umetrics, Sweden). Next, the regression model was established through partial least squares discriminant analysis (PLS-DA) along with the discriminant analysis. In addition, to improve the effectiveness of data analyses, orthogonal partial least squares discriminant analysis (OPLS-DA) was performed to correct the orthogonal transformation based on PLS-DA. OPLS-DA was applied to construct the relationship model between the metabolite expression and the sample category. The detected characteristic metabolites by ¹HNMR spectroscopy were identified with reference to the National Center for Biotechnology Information database (https://www.ncbi.nlm.nih.gov/) and Human Metabolome Database (http://www.hmdb.ca/). Furthermore, the metabolic pathway of the characteristic metabolites was established according to the MetaboAnalyst 5.0 database (https://www.metaboanalyst.ca/).

To compare the differences in weight, all mice were weighed daily. The weight after modeling decreased significantly (P < 0.01). Although the weight of POF mice increased after EA (P < 0.05), it was still lower than that of the control group (P < 0.01). The weights of the POF group and the EN group were not significantly different (P > 0.05) ( Figure 2A ).

On the basis of the different types of vaginal exfoliated cells, the estrous cycle of mice was identified as proestrus, estrus, metestrus, and diestrus. The estrous cycle of normal mice was 4–6 days. In proestrus, the vaginal exfoliated cells were mainly composed of nuclear epithelial cells and keratinocytes ( Figure 2B ). As the estrous cycle changed, most keratinocytes were observed in the estrous period ( Figure 2C ). In the metestrus period, nucleated epithelial cells, keratinocytes, and leukocytes were observed ( Figure 2D ). In the diestrus period, mainly nuclear epithelial cells and leukocytes were observed ( Figure 2E ). The vaginal exfoliated cells were mainly composed of nuclear epithelial cells and leukocytes after POF modeling was established ( Figure 2F ). Interestingly, the estrous cycle of mice in the EA group was changed regularly as the proestrus–estrus–metestrus–diestrus cycle progressed after treatment. On the contrary, all vaginal exfoliated cells in the mice of the EN group were still mainly composed of nuclear epithelial cells and leukocytes.

As shown in Figure 3A , atresia follicles increased, which was observed in the ovaries of POF mice. The number of mature follicles and primary follicles increased, and atresia follicles decreased in the EA group. Although a few antral follicles were observed in the EN group, no significant difference was found between the antral follicles in the EN group and the POF group ( Figures 3A ).

The positive rate of granulosa cell apoptosis increased in the POF group (P < 0.01). Compared with the POF group, the positive rate of granulosa cell apoptosis in the EA group was reduced (P < 0.05). However, the positive rate of the EA group was still higher than that of the control group (P < 0.01). Although the positive rate of the EN group was lower than that of the POF group, it was higher than that of the EA group (P < 0.01) ( Figures 3B ).

Compared with the control group, the serum levels of FSH and LH in the POF group increased (P < 0.05). Meanwhile, the levels of E₂ and AMH decreased after POF was established (P < 0.05). The levels of FSH and LH between the control group and the EA group were approached (P > 0.05). On the contrary, although the levels of FSH, LH, E₂, and AMH in the EN group were also regulated, they were significantly different from the control group. In addition, no significant difference was found in the levels of AMH among the EA group, EN group, and POF group (P > 0.05) ( Figures 4A–D ).

To determine the effect of EA on ovarian ERs, the relative expression of ER-α, ER-β, and ERP30 was assessed. The results showed that the expression of ER-α, ER-β, and ERP30 was lower than those in the control group after the POF modeling was established (P < 0.01). We found that ER-α, ER-β, and ERP30 increased after EA; however, they were lower than those in the control group (P < 0.01) ( Figures 4E–G ).

To determine the characteristic metabolites underlying the mechanism of EA on POF, the ¹HNMR profiles of extracts from the liver and kidney were performed and analyzed. The characteristic metabolites in the obtained ¹HNMR profiles of the liver and kidney were difficult to identify intuitively ( Figure 5 ). Therefore, the PLS-DA and orthogonal projections to latent structures discriminant analysis (OPLS-DA) were performed according to ¹HNMR profiles.

The four groups of ¹HNMR profiles were merged and analyzed. As shown in Figures 6A, B , the control group and the POF groups in the liver and kidney were distinctly dispersed. Concurrently, a good dispersion was present between the POF group and the EA group. Although the dispersion between the POF group and the EN group was observed intuitively in the liver, no significant difference was found in the metabolites between them. Poor dispersion was observed between the POF group and the EN group in the kidney. Therefore, the metabolism of the liver and kidney in POF mice was abnormal. The metabolites in the EN group and POF group were similar.

To determine the effect of EA on liver and kidney metabolism in POF, the metabolic profile of each group were analyzed in pairs. The analysis further proved that the dispersions between the control group and the POF group in the liver and kidney was good ( Figures 6C, D ). This finding indicated that the metabolites in the liver and kidney were abnormal after POF modeling. Similarly, the metabolites in the liver and kidney of POF mice were changed distinctly due to EA. The OPLS-DA of the liver and kidney showed that the metabolites of the EA group were greatly different from the POF group ( Figures 7A, B ). The result showed that the metabolites of the liver and kidney in POF can be regulated by EA. Therefore, S-plot and T-test were performed to identify the characteristic metabolites between the POF group and the EA group ( Figures 7a, b ). In contrast, OPLS-DA showed that no significant difference was present in the metabolites of the liver and kidney between the POF group and EN group. The result showed that EN could not improve the metabolic disorder of the liver and kidney in POF mice.

In the liver, the relative levels of lysine, glucose, phosphocholine, taurine, glycine, and glycerol increased after EA. On the contrary, the levels of glutamate, and creatine decreased. In the kidney, the relative levels of lactate and glycerol increased by EA. Furthermore, compared with the POF group, the relative levels of glutamate, creatine, threonine, serine, alanine, and lysine decreased. Creatine, glutamate, glycerol, and lysine were the overlapping metabolites between the liver and kidney after being treated by EA. Interestingly, the levels of overlapping characteristic metabolites in the liver and kidney showed consistent trends after EA.

All the characteristic metabolites of the liver and kidney in the EA group were regulated ( Figures 8 , 9 ). Compared with the POF group, the level of characteristic metabolites in the EA group was closer to that in the control group. The characteristic metabolites in the EN group and the POF group were similar. All metabolic pathways containing the above characteristic metabolites in the liver and kidney were established based on the MetaboAnalyst 5.0 database ( Figure 10 ).

To determine the relationship between POF and metabolites, the correlation of POF with characteristic metabolites was obtained by performing the Pearson analysis.

In the liver, FSH and LH were positively correlated with creatine, glutamate, glycerol, lysine, and taurine. FSH and LH were negatively correlated with glucose, glycine, and phosphocholine. AMH, E₂, ER-α, ER-β, and GPR30 were positively correlated with glucose, glycine, and phosphocholine. AMH, E₂, ER-α, ER-β, and GPR30 were negatively correlated with creatine, glutamate, glycerol, lysine, and taurine ( Figure 11 ).

In the kidney, FSH and LH were positively correlated with alanine, creatine, glutamate, glycerol, lysine, serine, and threonine. FSH and LH were negatively correlated with lactate. AMH, E₂, ER-α, ER-β, and GPR30 were positively correlated with lactate. AMH, E₂, ER-α, ER-β, and GPR30 were negatively correlated with alanine, creatine, glutamate, glycerol, lysine, serine, and threonine ( Figure 12 ).