Transcriptome Analysis of Multi-leaf and Standard Clover of F1 Hybrids of Mengnong Clover No.1 and White Clover

Caucasian clover (Trifolium ambiguum ) is native to the cold Caucasus mountains and also grows in alpine meadows, relatively arid low-elevation grass steppes, eastern Turkey and northern Iran (Yin, 2004; Zhou, 2017; Ren, 2014). Caucasian clover is a perennial forage grass of the genus Trifolium (Sleuch, 2000). It has a well-developed root system, with the main root penetrating deeply into the soil and it exhibits a strong reproductive ability (Yang, 2008). Extensive experiments and observations carried out over numerous years have offered in-depth understanding of its biological characteristics, genetic traits, production performance and adaptability. Caucasian clover shows a high forage yield, along with drought resistance and cold tolerance (Cui, 2006; He, 2008; Wan, 2013; Li, 2015). As a result, it is currently the only clover species capable of overwintering and being cultivated and used in inner Mongolia.
       
In contrast, white clover (T. repens), a shallow-rooted perennial legume herb native to Europe, is now globally distributed in temperate and subtropical high-altitude regions. As one of the most important legume forages worldwide, it exhibits strong natural regeneration through seed persistence (remaining viable for ~8 years). Additionally, its symbiotic relationship with rhizobacteria enables biological nitrogen fixation, providing essential nutrients for itself and neighboring plants.
       
In 1978, Elizabeth Williams successfully regenerated hybrid progeny plants through the in vitro culture of hybrid embryos from Caucasian clover and white clover. These plants could be propagated outdoors in large quantities (Williams, 2012). In 1995, Meredith used in vitro culture of hybrid ovules to backcross white clover with the progeny of Caucasian clover and white clover, successfully obtaining one seed (Meredith, 1995). In other studies, Caucasian clover has demonstrated excellent cold resistance. Its drought resistance, however, warrants further evaluation and study. Since 2005, Chinese researchers have introduced Caucasian clover and successfully obtained hybrid F1 plants through techniques such as hybrid embryo isolation culture and tissue regeneration. Using Mengnong clover No.1 as the maternal parent and white clover as the paternal parent, they have successfully expanded the nutrient base to develop a lot of plant population (Yang, 2008; Luan, 2009; Yang et al., 2009). Moreover, genetic and molecular marker (ISSR) analysis has confirmed that the hybrid F1 is a genuine hybrid progeny (Huang et al., 2012). White clover and mengnong clover no.1 are distantly related and the hybrid F1 generation shows high sterility. Hence, after years of research, researchers have successfully obtained F1 generation plants from the cross between Mengnong clover No.1(T. ambiguum Bieb. cv. Mengnong clover No.1) and white clover. The progeny of this cross can be preserved for a long time in both laboratory and field environments by taking advantage of its strong clonal propagation ability through histoculture succession and crossbreeding (Li et al., 2015).
       
Chen Lili and others analyzed the genetic differentiation between the parents and their F1 generation. They found that the genetic similarity coefficient between the F1 generation and the parent White Clover was 0.3465, while the coefficient with the parent Mengnong clover No. 1 was 0.5842. This indicated that the F1 generation had obtained a greater amount of genetic information from Mengnong clover No.1. (Chen et al., 2011). Li Xiaoxiu carried out a study on the morphological differences between the hybrid progeny of Mengnong clover No.1 and white clover. The results showed that these differences in the F1 progeny were determined by karyotype variations rather than chromosome number (Li, 2015). Huang Fan and others compared the photosynthetic physiological characteristics of Mengnong clover No. 1, white clover and their hybrid F1 generation. The results indicated that the hybrid F1 generation showed a relatively higher net photosynthetic rate and a lower transpiration rate. Moreover, its water-use efficiency was found to be significantly greater than that of both parent species, suggesting a high yield potential (Huang, 2012).
       
During the breeding of the F1 generation hybrid between Mengnong clover No.1 and white Clover, multi-leafed plants (with 4-5 leaflets) emerged. After further research, a multi-leafed strain was formed. Wu Na and others investigated the characteristics of the four multi-leafed strains of the F1 generation and found that the actual plant height of the multi-leafed strains was greater than that of the parent plants. Additionally, although the leaf width and leaf area of the multi-leafed strains were smaller, the leaf length was greater than that of the normal-type three-leafed strains. Regarding photosynthetic characteristics, the transpiration (Tr) values of the normal-type F1 generation hybrid strains were higher than those of the multi-leaf type strains. Although the Tr values of the five selected hybrids did not show significant differences from those of the parents, they exhibited significant differences from the parents in other traits. Nutrient analyses of the parents, leafy lines and common lines indicated that the conventional nutrient contents in the stems and leaves of the F1 generation hybrids were similar to those of the parent, Mengnong clover No.1. Moreover, the crude fat content in the stems of the leafy lines was higher than that in the common lines. The water use efficiency of the common lines was lower than that of the leafy lines and the parents. However, differences in the stem crude protein content among multi-leaf type strains at the initial flowering stage were not significant. The neutral detergent fiber content in the common type strain was lower than that in the leaves of the leafy type strain during the same growth cycle. Additionally, the neutral detergent fiber content in the leaves of Mengnong clover No.1 was higher than that of the parent and white clover (Wu, 2017).
       
Mengnong clover no.1 is a cultivated variety of Caucasian clover that was painstakingly selected and developed by Inner Mongolia agricultural university over a number of years. Mengnong clover No.1 possesses strong overwintering ability and remarkable drought resistance in most areas of Inner Mongolia, making it suitable for planting. However, it has limited nitrogen fixation ability. By employing in-vitro hybrid embryo culture and various tissue regeneration techniques, a large number of hybrid F1 generation plants have been produced. Previous studies have observed clover plants with leafy traits displaying 4-5 leaflets, yet these exhibited low fruiting rates. Intriguingly, both trilobed and multi-leafed traits were observed when using seed planting or root-tiller propagation. Furthermore, both multi-leaved and trilobed traits were observed when plants were cultivated in greenhouses or experimental fields. The advancements in transcriptome sequencing technologies have expanded the research scope in this field, enabling its application across a broader range of species. Transcriptome studies offer a novel theoretical framework for understanding gene function and structure. By leveraging high-throughput sequencing, it is possible to quickly and accurately obtain comprehensive sequence information for the majority of transcripts within specific tissues or organs under particular treatment conditions. This technology has been applied in fundamental research on crop breeding and various other fields.
       
To further explore the factors influencing the development of multi-leaf traits and uncover the molecular mechanisms governing the stable inheritance of multi-leaf clover, we isolated mRNA from the aboveground organs (stems and leaves) of both multi-leaf and standard clover. Subsequently, we constructed a cDNA library and sequenced the transcriptome using RNA-Seq technology. The genes associated with multi-leaf expression were identified and screened and bioinformatics analyses were conducted to establish a foundation for understanding the stable inheritance of multi-leaf traits. The RNA-Seq library facilitated transcriptome sequencing, gene screening related to leaf expression and bioinformatics analyses, which collectively provided a framework for comprehending the stable inheritance of leafy traits.
       
To elucidate the factors contributing to changes in multi-leaf clover traits and explore the molecular mechanisms underlying the stable inheritance of these traits, we employed the Illumina HiSeq150 high-throughput sequencing platform for transcriptome sequencing based on RNA-Seq technology. Subsequently, we utilized relevant bioinformatics methods to screen and analyze differentially expressed genes (DEGs), thereby establishing a theoretical foundation for the future cultivation of novel clover varieties.

The experiment  was conducted at vocational and technical college of inner mongolia agricultural university of inner mongolia autonomous region located at tumote right banner, baotou sity. This place is marked by a diverse continental semi-arid monsoon climate. It has an annual average temperature of 7.5°C. The highest average temperature in July reaches 22.9°C. The annual sunshine averages 3095 hours and the annual precipitation averages 346 mm. The frost-free period lasts for 130~140 days, with the first frost occurring in September and the late frost period ending in May. The soil in this area has a pH of approximately 7.5 and is classified as highly fertile sandy chestnut soil.
       
The experimental material consisted of six multi-leaf clovers (1-1) lines (D1, D2, D3, D4, D5 and D6) and nine standard three-leaf (2-2) lines. The details of the test materials taken and their grouping are given in Table 1. The test materials were surface sterilized by rinsing with 70% ethanol and then the plant material was dried. The stem-leaf mixed samples were stored in 2-mL screwcap cryovials, flash-frozen in liquid nitrogen for 3-4 hours and then transferred to -80°C for long-term storage.To guarantee sufficient RNA extraction, two multi-leaf plants were selected as samples.The numbering in T1, T2, T3, T4, T5 and T6 was used for follow-up test database construction and information analysis.


 
Functional annotation and differential expression gene analysis
 
Bowtie software was utilized to compare the sequencing reads with those in the unigene library. Based on the comparison results, the expression levels were estimated using RSEM. The FPKM values were employed to represent the expression abundance of the corresponding Unigene as follows,

 
DESeq2 software was employed to analyze the differential expression between the sample groups and identify differentially expressed genes between the two conditions. Genes with a p-value less than 0.05 or a fold change greater than 2 were classified as significantly differentially expressed. Subsequently, GO enrichment and KEGG pathway enrichment analyses were carried out on the differentially expressed genes using the hypergeometric distribution method (Guo, 2024).
 
Real-time fluorescence quantitative PCR
 
Sequences of the five target genes were extracted from the Unigene fa assembly file of the non-reference transcriptome. One internal reference gene was selected from the red clover GAPDH mRNA sequence and primers were designed correspondingly (Table 2).


       
A quantitative polymerase chain reaction (qPCR) system was established by using a real-time fluorescence qPCR instrument for each sample. For each gene in each sample, there were three replicate sets (Table 3).


       
The PCR reaction conditions were as follows: 94°C, 30 s; 94°C, 5 s; 58°C, 15 s; and 72°C, 10 s. The final three stages were performed for 45 cycles. Finally, the dissociation stage.

Analysis of differentially expressed genes between multi-leaf and standard clovers
 
The screening of DEGs involved the utilization of T1, T2 and T3 as test groups and T4, T5 and T6 as control groups. A total of 1985 DEGs were identified, including 794 upregulated genes and 1191 downregulated genes. In the process of screening for differential genes, customary thresholds were employed. Typically, empirical values such as |Fold Change| ≥ 2 and FDR < 0.05 were used (Table 4).


       
M-versus-A plots (MA plots) were generated based on the gene expression profiles (Fig 1). The MA plot offers a visual representation of the overall distribution of expression abundance and the differential ploidy of genes in the two sets of samples. In this representation, green indicates downregulated gene expression, red indicates upregulated gene expression and black dots denote genes with no significant differences in expression. Fig 1 also illustrates the differential expression patterns between standard and multi-leaf type clovers. Notably, there were 57,070 genes showing no significant differences, a quantity approximately 29 times greater than the number of DEGs.


       
Through a statistical analysis of the expression profiles of DEGs, we identified 68 genes that were specifically expressed in polyphyllous clover and were absent in standard clover. Conversely, we identified 121 specifically expressed genes in standard clover that were absent in polyphyllous clover.
 
DEG Gene ontology (GO) functional enrichment
 
Protein enrichment in differentially expressed up-regulated and down-regulated proteins was categorized using the GO secondary functions (Fig 2).


       
In biological processes, DEGs showed significant enrichment in various processes, including protein phosphorylation, defense response, cellular homeostasis, photosynthesis, response to oxidative stress, ATP metabolic process, DNA integration, hydrogen ion transmembrane transport, transmembrane transport, oxidant detoxification, seed development, osmotic stress response, cellular respiration, response to salt stress and cytoplasmic translation.
       
In cellular components, DEGs demonstrated significant enrichment in structures such as the nucleus, cytosol, mitochondria, ribosome, extracellular region, plasmodesma, cell wall, chloroplast subunit, thylakoid membrane, nucleolus, apoplast and cytosolic small ribosomal subunit.
       
Moving to molecular functions, DEGs were predominantly involved in oxidoreductase activity, protein binding, protein serine/threonine kinase activity, ADP binding, TPase activity, electron carrier activity, hydrogen ion transport across membranes, copper ion binding, rRNA binding, peroxidase activity and mRNA binding.
       
In this study, a GO enrichment analysis was performed on 1,985 DEGs, with a focus on biological processes. Among these, 94 DEGs were found to be involved in protein phosphorylation, including 24 upregulated genes and 70 downregulated ones. Defense response processes included 68 DEGs, of which 13 were upregulated and 55 were downregulated. Seventeen DEGs were associated with cellular homeostasis, with 5 upregulated and 12 downregulated. DNA integration processes involved 14 DEGs, consisting of 6 upregulated and 8 downregulated genes. Notably, five DEGs were found to participate in both protein phosphorylation and defense responses.When examining the DEGs from a molecular function perspective, remarkable enrichment was observed in the processes of protein serine/threonine kinase activity and oxidoreductase activity. Specifically, 77 DEGs were identified in the protein serine/threonine kinase activity process, comprising 19 upregulated genes and 58 downregulated ones.Under identical external environmental conditions, the composition and activity of hormone signaling pathways may differ between standard and multi-leaf clover, which can affect the expression of genes related to mitogen-activated protein kinase activity.
       
The activation of one hormone signaling pathway may inhibit another, leading to differences in gene expression between standard and multi-leaf clover.Furthermore, ATPase activity associated with transmembrane movement revealed 11 DEGs, consisting of 3 upregulated and 8 downregulated genes. Notably, three DEGs were common to both protein serine/threonine kinase activity and protein binding molecular functions, indicating enhanced activity in these two molecular functions. In contrast, only one DEG was shared among the combination of protein binding, ATP enzyme activity and transmembrane movement of substances. Additionally, six distinct DEGs were identified in both oxidoreductase activity and protein binding molecular functions. An examination of the GO annotation results for these categories and their 15 subcategories indicated a predominance of downregulated genes across all categories. This suggests that the metabolic processes associated with key pathways involving these genes are somewhat affected. The intricate network relationships among these DEGs play a crucial role in influencing plant tissue growth and development (Shang, 2016).
 
Annotation of DEG kyoto encyclopedia of genes and genomes (KEGG)
 
The KEGG annotation results of DEGs were classified into five main functional groups: Organismal systems, metabolism, genetic information processing, environmental information processing and cellular processes. The distribution of DEGs among these categories differed. KEGG annotations were divided into five functional categories and the proportion of the total number of annotated genes varied. Among the cellular processes, the endocytosis pathway had the largest number of genes (7), accounting for 3.83% of the total number of annotated genes. Among the environmental information processing pathways, the plant hormone signal transduction pathway had the largest number of genes (17), accounting for 9.29% of the total number of annotated genes. Among the gene information processing pathways, homologous recombination was annotated with 12 genes, accounting for 6.56%, followed by nucleotide excision repair, mismatch repair and protein processing in the endoplasmic reticulum. Eleven genes were annotated, accounting for 6.01% of the total number of annotated genes. This was followed by nucleotide excision repair, mismatch repair, protein processing in the endoplasmic reticulum and DNA replication. Among the metabolic pathways, starch and sucrose metabolism contained the largest number of annotated genes (15), accounting for 8.20% of the total number of annotated genes. Finally, in the organic system, there were 24 genes related to plant-pathogen interactions (the largest number of genes in this category), accounting for 13.11% of the total number of annotated genes (Fig 3).


       
The top 20 pathways with the most significant enrichment were meticulously screened and analyzed. These pathways were mainly enriched in various biological processes, including plant hormone signal transduction, mismatch repair, homologous recombination, DNA replication and nucleotide excision repair (Fig 4).


       
The plant-pathogen interaction signaling pathways exhibited the highest number of enriched DEGs at 24, whereas the lowest number, 11, was noted in the mismatch repair, DNA replication and nucleotide excision repair signaling pathways. Hongliang Zhang (Zhang, 2016) conducted a KEGG metabolic pathway analysis of two DEGs (Upadhyay, 2016), revealing that in the cold-resistant zucchini SS1, a total of 179 DEGs were annotated across 73 metabolic pathways in the KEGG database. Notably, DEGs were significantly enriched in the photosynthetic metabolic pathway (ko00195) and the phytohormone signal transduction pathway (ko4075). In contrast, the cold-sensitive zucchini exhibited 132 DEGs annotated to 59 metabolic pathways in the KEGG database, with a significant concentration in the phytohormone signaling pathway (ko04075) and the starch sucrose metabolism pathway (ko0050). These metabolic pathways encompass processes such as photosynthesis, transcription regulation, translation regulation, phytohormone signaling and other metabolic activities (Ramesh, 2016). A comparative analysis of the enriched metabolic pathways in DEGs across various zucchini varieties revealed an enrichment of the phytohormone signal transduction pathway following low-temperature stress. In the present study, we identified more DEGs enriched in plant-pathogen interaction metabolic pathways compared to other pathways (Zhang, 2016). This finding suggests that the specific expression of these genes may influence the ability of leafy clover to resist pathogen infection and invasion.
 
Important kyoto encyclopedia of genes and genomes (KEGG ) pathway analysis
 
In the homologous recombination pathway, the sites of the three enzymes were affected by 12 differential genes. The activity of RAD51 B was affected by the down-regulation of c84426.graph _ c1 and the activity of pol ∂  was affected by the down-regulation of c88352.graph _ c0. RPA activity was affected by the down-regulated expression of c68409.graph _ c0, c86105.graph _ c1, c90048.graph _ c2, 90796.graph _ c0, c99695.graph _ c1 and the up-regulated expression of c74406.graph _ c0, c86456.graph _ c0, c86571.graph _ c0, c97390.graph _ c0, c99817.graph _ c1 (Fig 5).


       
In the DNA replication pathway, the sites of three enzymes are affected by 11 differentially expressed genes. The activity of RPA is jointly affected by the up-regulated expressions of c68409.graph_c0, c86105.graph_c1, c90048.graph_c2, c90796.graph_c0 and c99695. graph_c1 and the down-regulated expressions of c74406.graph_c0, c86456.graph_c0, c86571.graph_c0, c97390.graph_c0 and c99817.graph_c1. The activity of DNA polymerase ∂ 2 is affected by the down-regulated expression of c88352.graph_c0. The activity of RFA1 is jointly affected by the up-regulated expressions of c68409.graph_c0, c86105.graph_c1, c90048.graph_c2, c90796.graph_c0 and c99695.graph_c1 and the down-regulated expressions of K07466 c74406.graph_c0, c86456.graph_c0, c86571.graph_c0, c97390.graph_c0 and c99817.graph_c1 (Fig 6).


       
The aforementioned results suggest that differential genes might be involved in the entire process of DNA replication by regulating the activity of related enzyme sites in the homologous recombination and DNA replication pathways, thus promoting the occurrence of multi-leaf traits.
       
It is widely accepted that phytohormone-triggered disease resistance occurs primarily in response to pathogenic bacterial infestations rather than in the absence of such pathogens (Danyang, 2018; Hardie, 2003). This observation suggests that the various phytohormone signaling pathways form a complex regulatory network (Mao et al., 2010) that governs plant growth, development and resistance adaptation. This intricate network of regulatory mechanisms enables plants to respond more rapidly to diverse pathogen stress conditions, facilitating quicker adjustments. Furthermore, it implies that clover will activate this complex regulatory network to adapt more swiftly to disease stress when infested by specific pathogens (Lorenzo et al., 2013).
 
Real-time qPCR (RT-qPCR) validation
 
RT-qPCR was employed for validation to ensure the accuracy of the transcriptome sequencing data and assess the expression levels obtained from RNA-Seq. GAPDH (NCBI Gene Bank accession number JF968420.1) was utilized as an internal reference gene, while five differentially expressed genes (c49410.graph_c0, c68409.graph_c0, c72784.graph_c0, c88352.graph_c0 and c93936. graph_c0) served as target genes. The RT-qPCR validation results, as depicted in Fig 7, revealed consistency with the transcriptome sequencing data. Notably, the differential gene c88352.graph_c0, which showed higher expression in standard clover compared to polyphyllous clover at the RNA-Seq level, was validated at the RT-qPCR level, aligning with the transcriptome sequencing results (Fig 7a). Similarly, the differentially expressed gene c93936. graph_c0, which was differentially expressed between polyphyletic and standard clover, exhibited higher expression in polyphyletic clover at the RNA-Seq level and RT-qPCR validation corroborated the transcriptome sequencing results (Fig 7b). Consistency between transcriptome sequencing and RT-qPCR validation was also observed for differentially expressed genes c49410.graph_c0 and c72784.graph_c0. Both genes were expressed significantly more in standard clover than in polyphyllotypes (Fig 7c and d). The differential gene c68409.graph_c0 was expressed significantly in multi-leaf type clover at higher levels than in standard-type clover and the results were consistent between RT-qPCR validation and RNA-Seq sequencing (Fig 7e).

Transcriptome data analysis revealed that there are 1985 differentially expressed genes between multi-leaf and standard clover. Specifically, it includes 794 up-regulated genes and 1191 down-regulated genes. Among these differentially expressed genes, 68 genes are specifically expressed in multi-leaf clover, while 121 genes are specifically expressed in standard clover.
               
GO functional annotation indicated that differentially expressed genes (DEGs) were mainly enriched in protein phosphorylation, defense response, nuclear process, cytoplasmic activity, oxidoreductase activity and protein binding. KEGG enrichment analysis demonstrated that in the homologous recombination pathway, the sites of three enzymes were affected by 12 differential genes. In the DNA replication pathway, the sites of three enzymes were influenced by 11 differential genes. This suggests that differential genes can promote the occurrence and development of multiple leaves by participating in important biological processes such as homologous recombination and DNA replication during cell meiosis. These findings contribute theoretical evidence to improve our understanding of molecular mechanisms governing multi-leaf clover.

The present study was supported by Inner Mongolia Natural Science Foundation ( 2025QN03062 ); Organic dry farming national key laboratory ( preparation ) construction project open fund project ( Z135050009017-K9).
 
Disclaimers
 
The views and conclusions expressed in this article are solely those of the authors and do not necessarily represent the views of their affiliated institutions. The authors are responsible for the accuracy and completeness of the information provided, but do not accept any liability for any direct or indirect losses resulting from the use of this content.
 
Informed consent
 
All animal procedures for experiments were approved by the Committee of Experimental Animal care and handling techniques were approved by the University of Animal Care Committee.

The authors declare that there are no conflicts of interest regarding the publication of this article. No funding or sponsorship influenced the design of the study, data collection, analysis, decision to publish, or preparation of the manuscript.