Non-targeted Metabolomics Reveals Metabolic Characteristics of Porcine Atretic Follicles

Follicular atresia is one of the main factors limiting the reproductive power of domestic animals. At present, the molecular mechanisms involved in porcine follicular atresia at the metabolic level remain unclear. In this study, we divided the follicles of Bama Xiang pigs into healthy follicles (HFs) and atretic follicles (AFs) based on the follicle morphology. The expression of genes related to atresia in granulosa cells (GCs) and the concentration of hormones in the follicular fluid (FF) from HFs and AFs were detected. We then used liquid chromatography–mass spectrometry-based non-targeted metabolomic approach to analyze the metabolites in the FF from HFs and AFs. The results showed that the content of estradiol was significantly lower in AFs than in HFs, whereas that of progesterone was significantly higher in AFs than that in HFs. The expression of BCL2, VEGFA, and CYP19A1 was significantly higher in HFs than in AFs. In contrast, the expression of BAX and CASPASE3 was significantly lower in HFs. A total of 18 differential metabolites (DMs) were identified, including phospholipids, bioactive substances, and amino acids. The DMs were involved in 12 metabolic pathways, including arginine biosynthesis and primary bile acid biosynthesis. The levels of eight DMs were higher in the HF group than those in the AF group (p < 0.01), and those of 10 DMs were higher in the AF group than those in the HF group (p < 0.01). These findings indicate that the metabolic characteristics of porcine AFs are lower levels of lipids such as phospholipids and higher levels of amino acids and bile acids than those in HFs. Disorders of amino acid metabolism and cholic acid metabolism may contribute to porcine follicular atresia.


INTRODUCTION
Reproductive ability, which is a base index, typically refers to productivity in the porcine industry. In mammals, the ovary produces generative cells. The number of follicles is an index that determines the reproductive ability of female animals, especially polyembryonic animals. However, when more than 90% of the follicles degenerate before ovulation, the condition is termed follicular atresia (1). It is well-known that follicular cells include not only oocytes but also follicular fluid (FF) and granulosa cells (GCs). Many factors affect follicular atresia; however, GC apoptosis is the main molecular mechanism (2,3). Additionally, the physiological activities of the hypothalamus-pituitary-gonadal axis occupy an influential position in follicular atresia (4). The molecules in FF, such as hormones [including melatonin (5) and growth hormone (6)] and cytokines [transforming growth factor beta (TGF-β) superfamily (7), insulin-like growth factor (IGF) (8), interleukin-6 (IL6) (9), and epidermal growth factor (EGF) (10)], also affect follicular atresia. Because of the transfer of plasma and the secretion of follicle cells, the FF component is highly complex. The concentration of FF components may be inconsistent in diverse follicular cell types, especially between atretic follicles (AFs) and healthy follicles (HFs). Nishimoto et al. (11) studied the difference between HFs and AFs in cattle and proposed that glucose and lactic acid may be markers of follicular atresia in FF. However, the consistency of these indices in Bama Xiang pigs is unknown.
Fortunately, sequencing technology has quickly developed with the development of biotechnology and computer technology. After RNA sequencing (RNA-seq) explained the difference at the transcription level, proteomics explained the difference at the translation level. Metabolomics is an emerging omics technology that explains differences at the metabolic level. We can perform qualitative and quantitative analyses of all small molecule metabolites (<1 kDa) in specific biological samples (such as blood, urine, and FF) by metabolomics. This technology is suitable for the study of physiological mechanisms and the discovery of biomarkers (12). At present, metabolomics based on FF samples have shown great potential in reproductive medicine. It is used to identify biomarkers that assess oocyte quality to predict the outcome of assisted reproduction (13). Additionally, metabolomic analysis of FF screened out biomarkers of multiple reproductive diseases such as polycystic ovary syndrome (14), endometriosis (10), and infertility (15) and characterized their global metabolic profile to provide a reference for the clinical diagnosis of these diseases and an in-depth study of their pathogenesis. This study indicated that this method can be used to identify biomarkers in animal husbandry and veterinary medicine.
Bama Xiang pigs are a Chinese indigenous pig and are famous for their small size, excellent meat quality, and early maturation (16). Bama Xiang boars produce sperm in the seminiferous tubules at the age of 29 days, and the sex maturation time is 76 days. For Bama Xiang sows, the sex maturation time is 110 days (17). However, the litter size of Bama Xiang pigs is lower than that of many Chinese indigenous breeds, such as Meishan pigs, Bamei pigs, and Erhualian pigs. Therefore, using metabolomics to identify biomarkers between HFs and AFs will be of great benefit in improving the understanding of Bama Xiang pig reproduction traits. Meanwhile, because of the parallels with humans in anatomy and physiology, Bama Xiang pigs could be used in human medical research (18)(19)(20)(21). Our study used a liquid chromatography-mass spectrometry (LC-MS)-based non-targeted metabolomic approach for the first time to analyze the FF of HFs and AFs in Bama Xiang pigs. The results will further increase the understanding of the mechanism of follicular atresia in Bama Xiang pigs and provide a reference for experimental animals to benefit human reproduction studies.

Follicle Classification
A total of 10 ovaries were collected from the right side of 6-month-old Bama Xiang pigs from the Bama Xiang Pig Agriculture and Animal Husbandry Co., Ltd. (Guangxi, China). After collection, the ovaries were rinsed with physiological saline, placed in a heat preservation pot containing physiological saline at 37 • C with high-pressure processing, and transported to the laboratory within 2 h. After rinsing with 75% ethanol, the ovaries were washed with physiological saline at 37 • C. Ophthalmic scissors were used to cut open and separate a single follicle (3 mm < n < 5 mm) from the ovarian tissue. The single follicle was repeatedly rinsed with 0.1 M phosphate buffer (pH 7.25). As previously mentioned (22), we divided the follicles into HF and AF groups based on follicle morphology. The morphology of the follicle mainly included the following two aspects: the appearance of the follicle (such as color, capillary, and transparency) (23) and state of cumulus-oocyte complex (COC) and GCs.

Collection of Granulosa Cells and Follicular Fluid
After morphological classification, we removed the oocytes, and the remaining liquids (including FF and GCs) were mixed. The FF and GCs were separated by centrifugation (4 • C, 1,200 rpm, 5 min). After washing 2-3 times with precooled phosphate buffered saline, the collected GCs were mixed with 1 ml TRIzol and stored at −80 • C. The FF was centrifuged again (4 • C, 8,000 rpm, 5 min) to collect the supernatant. Then, 30 µl of the supernatant was used to detect the concentration of estradiol (E2) and progesterone (PROG), and the rest was stored at −80 • C for LC-MS analysis. In the same individual, we collected only one HF and one AF sample. A total of 10 samples each were collected from the HF group and AF group.

Measurements of Estradiol and Progesterone
The hormone levels related to follicular atresia, such as E2 (24) and PROG (25) in FF, were measured using ELISA kits (Shanghai Enzyme Biotechnology Co., Ltd., China, ml002366, ml002422) to confirm the follicle classification. A total of 20 FF samples were analyzed using a microplate reader (Infinite M200 PRO, TECAN) according to the manufacturer's instructions. All assays were validated in our laboratory by showing parallels between serial sample dilutions and the provided assay standard curve (range: E2, 2.5-80 pg/ml; PROG, 0.625-20 ng/ml). Sensitivity of the E2 and PROG assays was 0.1 pg/ml and 0.1 ng/ml, respectively. The intra-assay coefficient of variation was <15% for each.

cDNA Synthesis and Quantitative RT-PCR for Related Genes
The TRIzol Reagent (Ambion, Life Technologies, New York, USA), chloroform, isopropanol, and ethanol were used to extract the total RNA. RNase-free water was used to dissolve total RNA from GCs. Oligo (dT) 23 VN (10 µM), 4× gDNA wiper Mix, and RNase-free ddH 2 O at 42 • C were used for 2 min to remove the DNA from the GCs. We used the 5× HiScript II

Follicular Fluid Sample Preparation
Prior to FF metabolite analysis, a 50-µl FF sample with 125 µl of methanol and 125 µl of acetonitrile was added to 1.5-ml Eppendorf (EP) tubes, vortexed for 30 s to precipitate protein, underwent a low-temperature ultrasound for 15 min, and then was centrifuged (13,000 rpm, 20 min, 4 • C). The supernatant fractions were collected and evaporated using a vacuum concentrator. The resulting dry residues were redissolved in 50 µl of acetonitrile:water (1:1) and centrifuged again (13,000 rpm, 20 min, 4 • C). Finally, the liquid supernatant was transferred to sampler vials for analysis using a ultraperformance LC-MS (UPLC-MS) system (27).

Data Analysis
The original data were processed using the Compound Discover 3.1 software for peak alignment of metabolites, to correct the retention time of metabolites, and to extract the peak area of metabolites. To identify the metabolite structures, accurate mass matching (<25 ppm) and matching with secondary spectra were used in the experiment. The peak areas of all metabolites were imported into the SIMCA-P program (version 14.1, Umetrics) for multivariate analysis, consisting of unsupervised principal component analysis (PCA) and supervised orthogonal partial least squares discriminant analysis (OPLS-DA). PCA was used to determine intragroup aggregation and inter-group separation tendencies, whereas OPLS-DA was performed to further determine inter-group differences. The OPLS-DA models were validated based on the interpretation of variation in Y (R2Y) and forecast ability based on the model (Q2) in cross-validation and permutation tests by applying 200 iterations. Univariate analyses were performed, including independent-samples t-test and fold change analysis. Significant differential metabolites (DMs) were screened using variable importance in projection (VIP) scores (VIP > 1) obtained from the OPLS-DA model and p-values (p < 0.05) using the t-test. In addition, the metabolic pathway analysis network tools MetPA (http://metpa.metabolomics.ca), MetaboAnalyst (https://www.metaboanalyst.ca/MetaboAnalyst/ home.xhtml), and the Encyclopedia of Genes and Genomes in the Kyoto Protocol (KEGG, http://www.genome.jp/kegg/) were used to identify potentially disordered metabolic pathways. The differences between groups for hormone concentration and the qRT-PCR results were analyzed using independent-samples t-test analysis in SAS 9.2, and the results are expressed as the mean ± standard error of the mean. Unless otherwise stated, differences were considered statistically significant at p < 0.05.

Appearance of Follicles
We identified one HF and one AF from the same ovary. The COCs, GCs, and FF were noticeably different between HFs and AFs. The cumulus, ooplasm, zona pellucida, GCs, and FF were damaged or deteriorated in AFs ( Table 3). In follicles of the same size, more blood vessels were observed in the HFs than in the AFs (Figure 1).

Determination of Estradiol and Progesterone in Follicular Fluid
Each sample in the HF and AF groups was analyzed three times. A total of 60 records of E2 and PROG were used to carry out the independent-samples t-test analysis in SAS 9.2. The concentrations of E2 in the HF and AF were 35.20 ± 8.30 and 16.43 ± 3.06 pg/ml, respectively. The concentrations of PROG in the HF and AF were 6.41 ± 2.45 and 11.89 ± 2.14 ng/ml, respectively. Compared with the HF, in the AF, the content of E2  in the FF was significantly reduced (p < 0.01) and that of PROG was significantly increased (p < 0.01) (Figure 2A).

Expressions of Genes Related to Atresia in Granulosa Cells
Similar to the content of E2 and PROG, the expression of genes related to atresia in GCs from the HF and AF was assessed three times in each sample. A total of 60 records of each gene from the 20 samples, including BCL2, VEGFA, CYP19A1, BAX, and CASPASE3, were used to carry out the independent-samples ttest analysis using SAS 9.2 software. While we used the gene expression of HF as 1, the expression levels of BCL2, VEGFA, CYP19A1, BAX, and CASPASE3 in AF were 0.22 ± 0.14, 0.53 ± 0.19, 0.53 ± 0.20, 2.26 ± 0.59, and 3.05 ± 1.76, respectively. The expression of BCL2, VEGFA, and CYP19A1 in HF was significantly higher than that in AF (p < 0.01). The expression of BAX and CASPASE3 in HF was significantly lower than that in AF (p < 0.01) (Figure 2B).

Metabolite Profiles of Follicular Fluid
To obtain an overview of the FF metabolic characteristics of Bama Xiang pig follicular atresia, non-targeted metabolic profiling was conducted using the UPLC-MS/MS method. In this experiment, 3,044 peaks in the positive and negative ion modes were obtained. After filtering using the mzCloud database, HMDB database, total ion chromatograms, molecular weight, retention time, peak area, and peak intensity, a total of 994 metabolites were matched. The PCA showed that the R2X model interpretation rates for the two groups under the positive and negative ion mode conditions were 0.62 and 0.679, respectively. Four quality control (QC) samples gathered together not only in the positive PCA score plot but also in the negative PCA score plot (Figure 3). The OPLS-DA showed that samples from the two groups were well-separated, and each tended to cluster in the same group (Figures 4A,B). In the positive ion mode of OPLS-DA, R2Y and Q2 were 0.986 and 0.517, respectively, whereas in the negative ion mode, R2Y and Q2 were 0.979 and 0.38, respectively. The intercepts of permutation tests in the positive ion mode and negative ion mode were −0.257 and −0.353, genes in the HF and AF samples, respectively. *** highly significant difference between the two groups (p < 0.01).
Frontiers in Veterinary Science | www.frontiersin.org  respectively (Figures 4C,D). A total of 18 metabolites, including L-glutamine, beta-lysine, testosterone cypionate, and L-arginine, were identified as DMs between the HF and AF groups (Table 4, Figure 5). The levels of eight DMs, including palmitamide, octadecanamide, and tetradecanamide, were higher in the HF than in the AF (p < 0.01), and those of 10 DMs, including Lglutamine, L-arginine, and guanine, were higher in the AF than in the HF (p < 0.01) (Figure 6). Using MetaboAnalyst 4.0 to analyze the 18 DMs showed that they were completely enriched in 12 metabolic pathways, such as arginine biosynthesis, alanine,  The VIP score plot for the two groups analyzed in the negative ion mode. p [1] was the load of each substance under the first principal component; p(corr) [1] was the correlation between each substance and the first principal component. The red plots were the differential metabolites.

DISCUSSION
First, we identified HFs and AFs based on the appearance of the follicles. We then measured the concentration of hormones in FF and gene expression in GCs to prove the identification accuracy of HFs and AFs. Similar to the research on VEGFA (28) in follicular atresia, the expression of VEGFA was approximately two times higher in the HF than in the AF in our study. Furthermore, HFs appeared to have more blood vessels than AFs.
This result indicated that the expression of VEGFA gene in HF plays a key role in blood vessel generation to maintain the health of the follicles. Meanwhile, the expression of BAX (29), BCL2 (29, 30), CASPASE3 (30), and CYP19A1 (31) and the concentration of E2 (24) and PROG (25) were significantly different, which is consistent with the results of prior research. All these studies illustrated that the grouping of our study was reasonable.
The results of this study suggest that the concentrations of bioactive phospholipid substances,   1-O-hexadecyl-lyso-sn-glycero-3-phosphocholine (Lyso-PAF), are lower in the AF than in the HF. In vivo, LPC and LPE can be converted into PC and PE by acylation, which is mediated by lysophospholipid acyltransferases (LPLATs). PC and PE are the most abundant phospholipids and important structural and functional components of the biological membrane in all mammals (32, 33). The physiological significance of LPE remains unknown, but LPC has been shown to participate in several physiological activities as an effective lipid mediator. LPC is beneficial for GC growth and oocyte maturation by activating the extracellular signal-regulated kinases and nitric oxide (NO) (34). The exogenous addition of LPC could reduce the inhibitory effect of zearalenone on the oocyte maturation process (35). Similar to our study, some recent studies based on FF metabolomics have found more evidence that LPC and LPE are involved in the regulation of follicle development (36,37). PAF is considered to be the most potent lipid mediator known to date, and the remodeling pathway mediated by LysoPAF is one of the important pathways of PAF biosynthesis (38). PAF has been reported to promote ovulation, facilitate oocyte maturation in vitro, enhance extracellular matrix formation around luteal cells, and promote the growth of vascular networks (39). This result is consistent with our findings, as the expression of VEGFA was higher in the HF than in the AF, and more blood vessels were observed in HFs. In contrast, estrogen has been proven to cause accumulation of PAF in the ovary by inhibiting the activity of PAF-acetylhydrolase, and the level of PAF is positively correlated with estrogen (40). More importantly, testosterone cypionate as a testosterone can transform to E2 by aromatase reaction (41), and testosterone cypionate level in the AF was significantly higher than that in the HF. Therefore, we hypothesized that testosterone cypionate in the HF was transformed into E2 to promote the accumulation of PAF, which benefited oocyte maturation and maintenance of HFs. However, further research is required to determine the relationship between these factors and follicular atresia. In our results, phospholipid levels were higher in the HF than in the AF. The level of phospholipids that promote follicle development may be achieved by the activation of extracellular signal-regulated kinases and production of NO; the lack of phospholipids may lead to damage of the zona pellucida and GCs, which further results in the deterioration of FF. An interesting finding of our study was that the levels of amino acids such as L-glutamine, L-pyroglutamic acid, and L-arginine in FF samples were increased by atresia. However, studies have shown that L-pyroglutamic acid can be converted into glutamate to synthesize glutathione (GSH). Reduced GSH is a powerful antioxidant that protects GCs against oxidative stress caused by reactive oxygen species, maintaining normal mitochondrial function and inhibiting cell apoptosis (42). L-arginine improves oocyte maturation and embryo development by increasing NO production (43). These amino acids appear to have a positive effect on follicle development; however, they are upregulated in the FF of AFs, which seems contradictory. However, the concentration of amino acids was not higher. A higher concentration of amino acids means that the balance of osmotic pressure is affected, and the disruption of normal amino acid metabolism may damage the balance of intracellular osmotic pressure and ammonia level in FF, which may injure the oocyte and GCs (44). Therefore, disorders of amino acid metabolism may be one of the factors contributing to follicular atresia in pigs.
Additionally, this study revealed differences in BA metabolism. We compared the FF from HFs and AFs and found that in the FF from AFs, cholic acid (primary BA, p < 0.1) and its synthetic precursor 7-hydroxy-3-oxy-4-cholesteric acid were significantly upregulated, while secondary BAs, including glycoursodeoxycholic acid (p < 0.1) and glycodeoxycholic acid (p < 0.1), were downregulated. BAs are cytotoxic molecules, and excessive BAs have been proven to disrupt the mitochondrial function of cells, induce oxidative stress, and ultimately damage the ovary (45). The aliphatic side chain of primary BAs is conjugated to an amide linkage (N-acyl amidation) with glycine or taurine to form secondary BAs and reduce BA toxicity (46). Our results indicated that cholic acid synthesis was increased and the degradation was decreased, which led to the accumulation of cholic acid in AFs. Finally, the disruption of cholic acid metabolism may damage the follicles by inducing oxidative stress.

CONCLUSION
In summary, this study analyzed the different metabolic profiles of FF between HFs and AFs in Bama Xiang pigs. Our results revealed that the metabolic characteristics of porcine AFs are lower levels of lipids (such as phospholipids) and higher levels of amino acids and BAs than those in HFs. Disorders of amino acid metabolism and cholic acid metabolism may contribute to follicular atresia in pigs.

DATA AVAILABILITY STATEMENT
The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation.

ETHICS STATEMENT
The animal study was reviewed and approved by the College of Animal Science & Technology, Guangxi University (Ethics approval reference number: GXU-2021-055).

AUTHOR CONTRIBUTIONS
JM, LS, and JC designed the study, conducted the experiments, and drafted the paper. YL, YW, and GQ conducted parts of the experiments. JL and GL revised the paper. All authors have read and agreed to the published version of the manuscript.

FUNDING
This study was supported by the National Modern Agricultural Industrial Technology System (nycytxgxcxtd-15-01) and the Science and Technology Major Project of Guangxi (Guike-AA17292002).