Targeted Metabolomics Provide Chemotaxonomic Insights of Medicago ruthenica, with Coupled Transcriptomics Elucidating the Mechanism Underlying Floral Coloration

Medicago ruthenica, a wild legume forage widely distributed in the Eurasian steppe, demonstrates high genetic and phenotypic variation. M. ruthenica with a purely yellow flower (YFM), differing from the general phenotype of M. ruthenica with a purple flower (PFM), was recently discovered. The similar characteristics of YFM with Medicago falcata have led to conflicting opinions on its taxonomy using traditional morphological methods. The lack of chemotaxonomy information about M. ruthenica species and the unclear flower coloration mechanisms have hampered their study. Here, we investigated M. ruthenica using targeted metabolomics based on the chemotaxonomy method and elaborated the floral coloration mechanisms using transcriptomics. The identified flavonoids were the same types, but there were different contents in YFM and PFM, especially the contents of cyanidin-3-O-glucoside (C3G), an anthocyanin that causes the purple-reddish color of flowers. The over-accumulation of C3G in PFM was 1,770 times more than YFM. Nineteen anthocyanin-related genes were downregulated in YFM compared with their expression in PFM. Thus, YFM could be defined as a variety of M. ruthenica rather than a different species. The loss of purple flower coloration in YFM was attributed to the downregulation of these genes, resulting in reduced C3G accumulation. The taxonomic characteristics and molecular and physiological characteristics of this species will contribute to further research on other species with similar external morphologies.


Introduction
Germplasms of wild relatives are urgently needed for forage crop breeding and increased agricultural food production [1]. Medicago ruthenica is a cross-pollinated, diploid (2n = 2x = 16), perennial legume forage crop and a close relative of alfalfa (Medicago sativa) [2]. This species is widely distributed in the steppes of Eurasia and has unique physiological mechanisms for responding to environmental stress compared with other Medicago species [3]. M. ruthenica shows high genetic variation, particularly for stressresistant genes, thus providing abundant resources for forage crop resistance breeding [4,5]. An alfalfa cultivar (M. sativa L. cv. Longmu No. 803), which was obtained by the hybridization of alfalfa (M. sativa L. zhaodong) and M. ruthenica, has been successfully cultivated and exhibits better survival than alfalfa in a cold environment and higher yield [6]. Therefore, the germplasms of wild M. ruthenica should be urgently collected for germplasm exploitation and breeding.
The program for the collection of germplasm of M. ruthenica species was established in the 1990s in China [4]. Ample natural germplasm resources are available for this species, and phenotypic characteristics vary widely among this germplasm. In the wild, researchers identified the species and obtained germplasm resources using traditional morphologybased methods. Generally, the corolla of M. ruthenica is purple-reddish on the outside and yellow-colored on the inside [2,7]. It is interesting to note that the purely yellow corolla types were found by Bu Ning in the second year, after seeds collected in the wild were planted in the field, which is similar to the corolla of Medicago falcata. Bu Ning was engaged in forage germplasm collection at the Grassland Research Institute of the Chinese Academy of Agricultural Sciences and collected seeds from the hilly terrain of Shaerqin Town, Hohhot, Inner Mongolia, in the 1990s. Plants always exhibit similar characteristics, and the use of morphological taxonomy gives rise to controversial conclusions [7]. It has also led to discordant opinions that the purely yellow-flowered type should be ascribed to M. falcata or it may be a new variant of M. ruthenica [8]. The traditional method is not adequate. However, the lack of chemotaxonomy information in M. ruthenica species and the unclear flower coloration mechanisms have hampered studies on M. ruthenica.
Therefore, in this study, we utilized chemotaxonomy to analyze the chemical difference between purple-flowered and purely yellow-flowered phenotypes of M. ruthenica and elaborated on the floral coloration mechanisms using transcriptomics. We aim to provide novel information for the discrimination and classification of M. ruthenica species.

Floral Color Phenotypes
During the investigation of M. ruthenica, two distinct petal color phenotypes were observed: one with a purple-reddish color outside and yellow inside the petals (PFM) and the other with purely yellow petals (YFM), as shown in Figure 1. The development and blooming of the flowers take seven days and can be divided into the following stages: initial floret separation and petals packaged in the calyx; the appearance of petals among the calyx lobes; petals emerging from the calyx; more than 2 mm of petals appearing from the calyx but the keel still wrapped by the vexil; and a ready to bloom flower followed by complete blooming at the final stage ( Figure 2). During flowering, the floral colors of PFM and YFM were stable. exhibits better survival than alfalfa in a cold environment and higher yield [6]. Therefore, the germplasms of wild M. ruthenica should be urgently collected for germplasm exploitation and breeding.
The program for the collection of germplasm of M. ruthenica species was established in the 1990s in China [4]. Ample natural germplasm resources are available for this species, and phenotypic characteristics vary widely among this germplasm. In the wild, researchers identified the species and obtained germplasm resources using traditional morphology-based methods. Generally, the corolla of M. ruthenica is purple-reddish on the outside and yellow-colored on the inside [2,7]. It is interesting to note that the purely yellow corolla types were found by Bu Ning in the second year, after seeds collected in the wild were planted in the field, which is similar to the corolla of Medicago falcata. Bu Ning was engaged in forage germplasm collection at the Grassland Research Institute of the Chinese Academy of Agricultural Sciences and collected seeds from the hilly terrain of Shaerqin Town, Hohhot, Inner Mongolia, in the 1990s. Plants always exhibit similar characteristics, and the use of morphological taxonomy gives rise to controversial conclusions [7]. It has also led to discordant opinions that the purely yellow-flowered type should be ascribed to M. falcata or it may be a new variant of M. ruthenica [8]. The traditional method is not adequate. However, the lack of chemotaxonomy information in M. ruthenica species and the unclear flower coloration mechanisms have hampered studies on M. ruthenica.
Therefore, in this study, we utilized chemotaxonomy to analyze the chemical difference between purple-flowered and purely yellow-flowered phenotypes of M. ruthenica and elaborated on the floral coloration mechanisms using transcriptomics. We aim to provide novel information for the discrimination and classification of M. ruthenica species.

Floral Color Phenotypes
During the investigation of M. ruthenica, two distinct petal color phenotypes were observed: one with a purple-reddish color outside and yellow inside the petals (PFM) and the other with purely yellow petals (YFM), as shown in Figure 1. The development and blooming of the flowers take seven days and can be divided into the following stages: initial floret separation and petals packaged in the calyx; the appearance of petals among the calyx lobes; petals emerging from the calyx; more than 2 mm of petals appearing from the calyx but the keel still wrapped by the vexil; and a ready to bloom flower followed by complete blooming at the final stage ( Figure 2). During flowering, the floral colors of PFM and YFM were stable.   (11-12 am) and stored at 4 °C to keep them fresh; then, the samples were inspected under a stereomicroscope and photographed. Flower face (left), side (middle), and bottom (right) images are shown. During the investigation of M. ruthenica, two distinct flower color phenotypes were observed, with one having the purple-reddish outside and yellow petal inside (PFM) and the other with the purely yellow petal (YFM). and stored at 4 °C to keep them fresh. The samples were then inspected under a stereomicroscope and photographed. It takes seven days for the development and blooming of the flower tissue. The developmental process from the stage of floret includes initial separating and calyx packaging the petals, the stage of petals appearing among the calyx lobes, the stage of the petals exceeding the calyx, the stage of the petals exceeding the calyx more than 2 mm but with the vexil still wrapped by the keel, the stage of the flower ready to bloom, and finally the flower in full bloom, shown from left to right. The floral colors of PFM and YFM were stable.

Sample Quality Control and the Discrimination of Samples in the Validation Set
To identify the potential metabolites related to petal color in M. ruthenica species, we used the UPLC-MS-based targeted metabolomics approach for qualitative and quantitative analyses of the flavonoids. PCA and OPLS-DA are widely used in metabolomics for multivariate statistical analysis. PCA is used to analyze all data, highlighting specific samples, whereas OPLS-DA is used for sample grouping and focuses on analyzing differences in grouped samples. Principal component 1 (PC1) and principal component 2 (PC2) explained 96% and 2.3% of the variance, respectively, which separated PFM and YFM (Figure S1a). The OPLS-DA value indicated the complete separation of PFM and YFM, suggesting a significant difference in the flavonoid content of PFM and YFM ( Figure S1b). Therefore, PCA and OPLS-DA are accurate and could be used for further screening of differentially accumulated metabolites (DAMs), particularly flavonoids.

Composition of Anthocyanins and Flavonoids and Concentrations' Discrepancy in YFM and PFM
The flavonoid metabolites were identified using UPLC-MS, the metabolite data acquisition of samples was performed in both positive and negative ionization modes (Figure 3a), and 48 DAMs were screened based on PCA and OPLS-DA (Figure 3b). There are two metabolic phenotypes in PFM and YFM (Figure 3c). Compared with those in PFM, 25 The developmental process from the stage of floret includes initial separating and calyx packaging the petals, the stage of petals appearing among the calyx lobes, the stage of the petals exceeding the calyx, the stage of the petals exceeding the calyx more than 2 mm but with the vexil still wrapped by the keel, the stage of the flower ready to bloom, and finally the flower in full bloom, shown from left to right. The floral colors of PFM and YFM were stable.

Sample Quality Control and the Discrimination of Samples in the Validation Set
To identify the potential metabolites related to petal color in M. ruthenica species, we used the UPLC-MS-based targeted metabolomics approach for qualitative and quantitative analyses of the flavonoids. PCA and OPLS-DA are widely used in metabolomics for multivariate statistical analysis. PCA is used to analyze all data, highlighting specific samples, whereas OPLS-DA is used for sample grouping and focuses on analyzing differences in grouped samples. Principal component 1 (PC1) and principal component 2 (PC2) explained 96% and 2.3% of the variance, respectively, which separated PFM and YFM ( Figure S1a). The OPLS-DA value indicated the complete separation of PFM and YFM, suggesting a significant difference in the flavonoid content of PFM and YFM ( Figure S1b). Therefore, PCA and OPLS-DA are accurate and could be used for further screening of differentially accumulated metabolites (DAMs), particularly flavonoids.

Composition of Anthocyanins and Flavonoids and Concentrations' Discrepancy in YFM and PFM
The flavonoid metabolites were identified using UPLC-MS, the metabolite data acquisition of samples was performed in both positive and negative ionization modes (Figure 3a), and 48 DAMs were screened based on PCA and OPLS-DA (Figure 3b). There are two metabolic phenotypes in PFM and YFM (Figure 3c). Compared with those in PFM, 25 metabolites were downregulated in YFM, including cyanidin-3-O-glucoside (C3G) and C-glycosylflavones, 14 were upregulated in YFM, and 9 were not significantly different between YFM and PFM ( Figure 3c and Figure S2). These DAMs could be considered potential chemical markers for distinguishing the germplasms of PFM and YFM. Forty-eight DAMs were screened based on PCA and OPLS-DA (b). In this figure, each peak represents a kind of anthocyanin or flavonoid. The contents of anthocyanins and flavonoids were calculated using the area of related peak. (c) Heatmap of differentially accumulated flavonoid metabolites (DAMs) in PFM and YFM. A heatmap was drawn with TBtools using log 2 fold-change values of the concentrations of these 48 DAMs. Each line presents the color group representing related metabolites, which are listed on the right. The colors reflect the contents of metabolites, cells with orangecolor represent high accumulation in samples, and the color scale from blue to red represents low to high accumulation in samples, respectively. PFM−1-3 and YFM−1-3 represent the three biological independent PFM and YFM samples, respectively. In this study, PFMs served as the control group.
As summarized in Table 1, these forty-eight DAMs included twenty flavonoids, six dihydroflavonoids, eleven flavonols, four flavanols, three chalcones, three isoflavones, and one anthocyanin. The types of anthocyanins and flavonoids were identical in YFM and PFM. C3G was the major anthocyanin in YFM and PFM (Table 1). Interestingly, C3G accumulation was reduced in YFM. The over-accumulation of the C3G in PFM was 1770 times more than YFM (Table 1), implying changes in the regulation of anthocyanin synthesis-related genes. was used for signal acquisition and these separated components' identification. The signals of these DAMs were acquired and identified in both positive ion mode and negative ion mode (a). Forty-eight DAMs were screened based on PCA and OPLS-DA (b). In this figure, each peak represents a kind of anthocyanin or flavonoid. The contents of anthocyanins and flavonoids were calculated using the area of related peak. (c) Heatmap of differentially accumulated flavonoid metabolites (DAMs) in PFM and YFM. A heatmap was drawn with TBtools using log 2 fold-change values of the concentrations of these 48 DAMs. Each line presents the color group representing related metabolites, which are listed on the right. The colors reflect the contents of metabolites, cells with orange-color represent high accumulation in samples, and the color scale from blue to red represents low to high accumulation in samples, respectively. PFM−1-3 and YFM−1-3 represent the three biological independent PFM and YFM samples, respectively. In this study, PFMs served as the control group.
As summarized in Table 1, these forty-eight DAMs included twenty flavonoids, six dihydroflavonoids, eleven flavonols, four flavanols, three chalcones, three isoflavones, and one anthocyanin. The types of anthocyanins and flavonoids were identical in YFM and PFM. C3G was the major anthocyanin in YFM and PFM (Table 1). Interestingly, C3G accumulation was reduced in YFM. The over-accumulation of the C3G in PFM was 1770 times more than YFM (Table 1), implying changes in the regulation of anthocyanin synthesis-related genes. FC represents the fold-change and the multiple relationship of the substance after the comparison between the two groups. Rt represents retention time. Related contents were calculated using the area of mass spectrum peaks. Forty-eight differentially accumulated flavonoid metabolites (DAMs) were identified, comprising twenty flavonoids, six dihydroflavonoids, eleven flavonols, four flavanols, three chalcones, three isoflavones, and one anthocyanin. The signals of these DAMs were acquired and identified in both positive ion mode and negative ion mode. In column polarity, "1" represents positive ion mode and "−1" represents negative ion mode. The contents of 25 metabolites were lower in YFM than in PFM, especially that of cyanidin-3-O-glucoside.

Identification of 3319 Differentially Expressed Genes (DEGs) Using Transcriptome Sequencing
We used transcriptome sequencing to analyze DEGs in purple-reddish and purely yellow-colored M. ruthenica flowers to further investigate the molecular basis of the flower color variation. The total RNAs extracted from the petals of both plants were sequenced, yielding 42.11 Gb of clean data (Table S1). More than 84.36% of the clean reads were mapped to the M. ruthenica reference genome [9]. Mean Q20 and Q30 values were ≥98.48% and ≥95.25%, respectively. The GC content was ≥42.48%, indicating high-quality reads that could be used for differential gene expression analysis. PCA of PFM and YFM samples showed that the PC1 contribution was 76.1%, which showed a clear separation of PFM and YFM (Figure 4a).

Identification of 3319 Differentially Expressed Genes (DEGs) Using Transcriptome Sequencing
We used transcriptome sequencing to analyze DEGs in purple-reddish and purely yellow-colored M. ruthenica flowers to further investigate the molecular basis of the flower color variation. The total RNAs extracted from the petals of both plants were sequenced, yielding 42.11 Gb of clean data (Table S1). More than 84.36% of the clean reads were mapped to the M. ruthenica reference genome [9]. Mean Q20 and Q30 values were ≥98.48% and ≥95.25%, respectively. The GC content was ≥42.48%, indicating high-quality reads that could be used for differential gene expression analysis. PCA of PFM and YFM samples showed that the PC1 contribution was 76.1%, which showed a clear separation of PFM and YFM (Figure 4a).   (Table S2).
According to NR annotation, the highest homology matched for PFM and YFM was that for Medicago truncatula (85%) (Figure 4b). According to GO classification, 36,566 unigenes were classified into 53 functional terms; 18 terms were categorized in the biological process, 20 in the cellular component, and 15 in the molecular function (Figure 5a). Among them, the groups of the metabolic process (16,979) and cellular process (16,952), cell (17,257) and cell part (17,257), and binding (18,995) and catalytic activity (17,203) were the highestranked categories for the biological process, cellular component, and molecular function, respectively (Table S3).  (Table S2).
According to NR annotation, the highest homology matched for PFM and YFM was that for Medicago truncatula (85%) (Figure 4b). According to GO classification, 36,566 unigenes were classified into 53 functional terms; 18 terms were categorized in the biological process, 20 in the cellular component, and 15 in the molecular function (Figure 5a). Among them, the groups of the metabolic process (16,979) and cellular process (16,952), cell (17,257) and cell part (17,257), and binding (18,995) and catalytic activity (17,203) were the highest-ranked categories for the biological process, cellular component, and molecular function, respectively (Table S3).   For the annotation of potential gene functions, we divided these unigenes into 25 categories based on the KOG database, among which 1195 unigenes were annotated into the group of "secondary metabolites biosynthesis, transport, and catabolism" (Figure 5b). Moreover, 20,285 unigenes were matched to 137 KEGG pathways (Table S4). The most enriched pathway focused on metabolism, particularly the phenylpropanoid biosynthesis pathway. Among them, 491, 158, 4, 70, and 126 key genes were predicted to be involved in the phenylpropanoid, flavonoid, anthocyanin, flavone and flavonols, and isoflavonoid biosynthesis pathways, respectively (Table S4). In total, 3319 DEGs were identified (Table S5). Among them, 1261 and 1621 genes were upregulated (Table S6) and downregulated when comparing YFM with PFM, respectively (Table S7). According to the KEGG annotations, DEGs involved in flavonoid biosynthesis (19 unigenes) and isoflavonoid biosynthesis (10 unigenes) were downregulated (Figure 6a), and genes involved in flavone and flavonol biosynthesis (9 unigenes) and isoflavonoid biosynthesis (15 unigenes) were upregulated (Figure 6b). synthesis, transport, and catabolism".
For the annotation of potential gene functions, we divided these unigen categories based on the KOG database, among which 1195 unigenes were anno the group of "secondary metabolites biosynthesis, transport, and catabolism" (F Moreover, 20,285 unigenes were matched to 137 KEGG pathways (Table S4). enriched pathway focused on metabolism, particularly the phenylpropanoid bi pathway. Among them, 491, 158, 4, 70, and 126 key genes were predicted to be in the phenylpropanoid, flavonoid, anthocyanin, flavone and flavonols, and iso biosynthesis pathways, respectively (Table S4). In total, 3319 DEGs were identif S5). Among them, 1261 and 1621 genes were upregulated (Table S6) and down when comparing YFM with PFM, respectively (Table S7). According to the KE tations, DEGs involved in flavonoid biosynthesis (19 unigenes) and isoflavono thesis (10 unigenes) were downregulated (Figure 6a), and genes involved in fla flavonol biosynthesis (9 unigenes) and isoflavonoid biosynthesis (15 unigenes) regulated (Figure 6b).  number of downregulated genes common to the corresponding KEGG annotation. In total, the number of downregulated genes involved in categories of flavonoid biosynthesis and isoflavonoid biosynthesis was 19 and 10, respectively. (b) KEGG analysis of upregulated genes from the RNAseq data. The Y-axis on the left represents the categories of the KEGG annotation, and the X-axis indicates the number of upregulated genes that are common to the corresponding KEGG annotation. The number of upregulated genes involved in categories of flavone and flavonol biosynthesis and isoflavonoid biosynthesis was 9 and 15, respectively.

Differential ABP Gene Expression and Their Quantitative Real-Time PCR (qRT-PCR) Validation
We further explored ABP DEGs to elucidate the floral coloration mechanism in PFM and YFM. Anthocyanin synthesis was regulated by the following structural genes-CHS, CHI, F3H, F3 H, DFR, ANS, and UFGT. CHS, CHI, F3H, and F3 H were upstream genes, and DFR, ANS, and UFGT were downstream genes. A total of 23 DEGs were considered candidate genes ( Table 2). Compared with those in PFM, a total of 19 genes (4 CHS, 2 CHI, 3 F3H, 2 F3 H, 4 DFR, 2 ANS, and 2 UFGT) were downregulated, and 4 UFGT genes were upregulated in YFM. Thus, different ABPs were inferred in PFM and YFM; a schematic model to better understand the formation of pale-yellow M. ruthenica flowers is presented in Figure 7. To verify the transcriptome data, the expression levels of these 23 candidate genes were examined using qRT-PCR ( Figure 8); primers used for qRT-PCR are listed in Table S8. Most of the genes responsible for anthocyanin synthesis were significantly downregulated in YFM, whereas three UFGT genes showed significantly higher expression levels in YFM than in PFM, which was consistent with the results of DEG analysis. To verify the transcriptome data, the expression levels of these 23 candidate gen were examined using qRT-PCR ( Figure 8); primers used for qRT-PCR are listed in Tab S8. Most of the genes responsible for anthocyanin synthesis were significantly downre ulated in YFM, whereas three UFGT genes showed significantly higher expression leve in YFM than in PFM, which was consistent with the results of DEG analysis.

Discussion
Botanical taxonomy is the basis for understanding and classifying plants, exploring plants' diversity, and utilizing and conserving plant resources. The categorization and phylogeny of plants relying on traditional external morphological characteristics are known to be conflicting [10]. DNA molecular markers have successfully assisted in plant classification [11]. The limitation of this technique is that it can still only partially reveal the relationship among plant species by gene fragments rather than the complete genome [12]. Therefore, other methods need to be employed.
Notably, the secondary metabolite components are often similar in plants within a taxonomical unit [13]. This indicates that phytochemistry can be used to provide additional evidence for plant taxonomy [14]. With the availability of metabolomics, it has been used generally in Vicia and Siegesbeckiae Herba species for plant species identification [15,16]. Among the wide variety of secondary metabolites, flavonoids, and anthocyanins are potential chemotaxonomic markers [17]. Flavonoids are natural pigments, and the combinations of flavonoid compounds differ significantly in plant species, thus being responsible for the yellow or white color of flowers [12,18]. Anthocyanins originate from an essential branch of flavonoids' biosynthesis, and more than 500 anthocyanins can be produced by modifications such as hydroxylation, glycosylation, or methylation [19,20]. The types and content of anthocyanins can also result in flower color polymorphisms that cause variation from pink to red, blue, or purple [21,22]. For example, flowers rich in cyanidin, pelargonidin, or delphinidin display purple, red, or blue petals, respectively [23].

Discussion
Botanical taxonomy is the basis for understanding and classifying plants, exploring plants' diversity, and utilizing and conserving plant resources. The categorization and phylogeny of plants relying on traditional external morphological characteristics are known to be conflicting [10]. DNA molecular markers have successfully assisted in plant classification [11]. The limitation of this technique is that it can still only partially reveal the relationship among plant species by gene fragments rather than the complete genome [12]. Therefore, other methods need to be employed.
Notably, the secondary metabolite components are often similar in plants within a taxonomical unit [13]. This indicates that phytochemistry can be used to provide additional evidence for plant taxonomy [14]. With the availability of metabolomics, it has been used generally in Vicia and Siegesbeckiae Herba species for plant species identification [15,16]. Among the wide variety of secondary metabolites, flavonoids, and anthocyanins are potential chemotaxonomic markers [17]. Flavonoids are natural pigments, and the combinations of flavonoid compounds differ significantly in plant species, thus being responsible for the yellow or white color of flowers [12,18]. Anthocyanins originate from an essential branch of flavonoids' biosynthesis, and more than 500 anthocyanins can be produced by modifications such as hydroxylation, glycosylation, or methylation [19,20]. The types and content of anthocyanins can also result in flower color polymorphisms that cause variation from pink to red, blue, or purple [21,22]. For example, flowers rich in cyanidin, pelargonidin, or delphinidin display purple, red, or blue petals, respectively [23]. Thus, in this study, we discussed the differences in chemical compounds between purple-flowered and yellow-flowered phenotypes of M. ruthenica; we provided information for the classification of M. ruthenica species, and elaborated on the floral coloration mechanisms in M. ruthenica species.
In total, 48 flavonoid metabolites were identified in M. ruthenica, with no difference in the types of metabolites between YFM and PFM (Table 1). These flavonoids could be considered potential chemical markers for discriminating M. ruthenica from other species. Among the Medicago species, more than 28 flavonoid metabolites in M. truncatula were different from those in M. ruthenica, including formononetin-7-O-glucoside, liquiritigenin, and purpurin [24]. A total of 13 flavonoids in alfalfa were different from those found in M. ruthenica such as malvidin 3-O-glucoside, petunidin 3-O-glucoside, and robinin [25]. Six flavonoid metabolites in M. falcata differed from those in M. ruthenica, particularly salidroside, laricitrin, and daidzein [26]. These unique-species flavonoid metabolites were not detected in either YFM and PFM. Hongmei et al. [27] found that the chromosome set was highly similar in yellow-flowered M. ruthenica and purple-flowered M. ruthenica based on karyotype analysis, which suggests a close relationship between them. The study of inter simple sequence repeats also revealed that M. ruthenica with purely yellow flowers was closely related to M. ruthenica with purple-reddish flowers, rather than M. falcata [5]. Hence, M. ruthenica with purely yellow flowers could be defined as a variety of M. ruthenica.
Among these detected flavonoids, C3G was the major anthocyanin in the M. ruthenica species. The over-accumulation of the C3G in PFM was 1,770 times more than YFM (Table 1), which implies changes in the regulation of anthocyanin pathway genes. The relationship between phenotypic variation, metabolic products, and the regulation of anthocyanin pathway genes has been described through studies on colorful flowers [28]. Either a change in gene expression or loss-of-function of the ABP genes could be responsible for the failure of anthocyanidin accumulation, which further impacts the floral color. Compared with PFM, anthocyanin biosynthesis in YFM is blocked by low expression levels of most structural genes, including 4 CHS, 2 CHI, 3 F3H, 2 F3 H, 4 DFR, 2 ANS, and 2 UFGT genes, which might disrupt anthocyanin accumulation, leading to the loss of the purple color and the presentation of pale-yellow petals. (Figure 7). Single or multiple ABP genes are downregulated, resulting in reduced anthocyanin accumulation, which has been found in lotus (Nelumbo nucifera) with yellow flowers and alfalfa with white flowers [29][30][31]. Thus, the loss of purple flower coloration in YFM was attributed to the downregulation of these genes, resulting in reduced C3G accumulation.

Plant Materials
Germplasms of PFM and YFM were provided by China National Forage Germplasm Bank and planted in the experimental station (with a random plot size) at the Institute of Grassland Research, China Academy of Agricultural Sciences, in 2016. The floral color of these two M. ruthenica phenotypes was stable for the entire duration of the study. Petal samples were collected during the flowering and bud phases for targeted metabolomic and transcriptomic analyses, respectively. A total of three biological samples, each from three independent plants, were used in the experiment.

Metabolite Extraction and Detection
To determine total flavonoids, extraction and targeted metabolite profiling were performed by Ollwe gene Technologies Co., Ltd. (Nanjing, China). Briefly, 100 mg powder of each freeze-dried petal sample was extracted using 3 mL of 75% methanol (containing 1% acetic acid). After 30 s of vortexing, the mixture was sonicated at 4 • C for 30 min and centrifuged at 10,000× g for 15 min. Finally, the supernatant solution was filtered (0.22 µm pore size) and used for further analysis by ultra-performance liquid chromatography (UPLC). In addition, a quality control sample was prepared by mixing equal aliquots of the supernatants from all of the samples [32]. ical characteristics will contribute to further research on forage crop species with similar external morphologies and their breeding.