Genome-wide association studies provide genetic insights into natural variation of seed-size-related traits in mungbean

Although mungbean (Vigna radiata (L.) R. Wilczek) is an important legume crop, its seed yield is relatively low. To address this issue, here 196 accessions with 3,607,508 SNP markers were used to identify quantitative trait nucleotides (QTNs), QTN-by-environment interactions (QEIs), and their candidate genes for seed length (SL), seed width, and 100-seed weight (HSW) in two environments. As a result, 98 QTNs and 20 QEIs were identified using 3VmrMLM, while 95, >10,000, and 15 QTNs were identified using EMMAX, GEMMA, and CMLM, respectively. Among 809 genes around these QTNs, 12 were homologous to known seed-development genes in rice and Arabidopsis thaliana, in which 10, 2, 1, and 0 genes were found, respectively, by the above four methods to be associated with the three traits, such as VrEmp24/25 for SL and VrKIX8 for HSW. Eight of the 12 genes were significantly differentially expressed between two large-seed and two small-seed accessions, and VrKIX8, VrPAT14, VrEmp24/25, VrIAR1, VrBEE3, VrSUC4, and Vrflo2 were further verified by RT-qPCR. Among 65 genes around these QEIs, VrFATB, VrGSO1, VrLACS2, and VrPAT14 were homologous to known seed-development genes in A. thaliana, although new experiments are necessary to explore these novel GEI-trait associations. In addition, 54 genes were identified in comparative genomics analysis to be associated with seed development pathway, in which VrKIX8, VrABA2, VrABI5, VrSHB1, and VrIKU2 were also identified in genome-wide association studies. This result provided a reliable approach for identifying seed-size-related genes in mungbean and a solid foundation for further molecular biology research on seed-size-related genes.


Background
Mungbean (Vigna radiata (L.) R. Wilczek) is a basic source of protein and carbohydrate, as it contains approximately 20% protein and 75% carbohydrate, and is a traditional and important legume in Asia (Somta et al., 2007). Due to its short life cycle (60-75 days), relative drought tolerance, and the ability to restore atmospheric nitrogen in association with Rhizobium/Bradyrhizobium bacteria, mungbean plays a crucial role in cropping systems and soil improvement (Somta et al., 2007;Alam et al., 2014).
The crop is generally grown as a cash crop in cerealbased farming systems. However, the major constraint in mungbean production is low seed yield. The average seed yield of mungbean is only approximately 700 kg per ha (Islam et al., 2015). Therefore, improving seed yield is the main goal in mungbean breeding. Understanding the genetic basis underlying seed-size-related traits is critical for the genetic improvement of mungbeans. In mungbeans, the ideotype of high-yielding cultivars are generally characterized by a large seed size, a short and synchronous maturity, a low sensitivity or insensitivity to day length, and the resistances to insects and disease (Fernandez et al., 1988). However, the knowledge on genes related to seed size has been limited. Moreover, the genes involved in the pathway of seed developments are not yet fully known.
Seed weight is the most important yield component and directly proportional to seed yield per plant in mungbean. To date, there have been seven studies of QTLs for seed weight in mungbean. Most of these studies are based on bi-parental segregation populations derived from interspecific crosses between cultivated and wild (V. radiata var. sublobata) mungbeans, and only two studies have evaluated seed size in more than one environments. The number of QTLs identified in those studies ranged from 3 to 11. Humphry et al. (2010) reported 11 loci for seed weight using SSR-marks, and Mei et al. (2009) identified a major QTL associated with both bruchid resistance and seed mass. Nonetheless, no candidate gene was identified for this trait.
Knowledge regarding seed development pathway is also a valuable source for transgenic strategies to improve crop production. As reported, there are several signaling pathways that control seed size, including the G-protein signaling, ubiquitin proteasome pathways, mitogen-activated protein kinase (MAPK) signaling, auxin pathways, and some transcriptional regulators (Li et al., 2019). In Arabidopsis, GPA1, AGB, and AGG3 were involved in G-protein-signaling pathways. DA1, DA2, SOD2, UBP15, EOD1, and SAMBA were involved in ubiquitin proteasome pathways. In addition, ABA2, ABI5, SHB1, MINI3, IKU2, and CKX were involved in the HAIKU (IKU) pathway. Additional genes were found to be related to seed size developments, but their pathways are uncertain, such as KIX8, BES1, MES1, and KLU (Orozco-Arroyo et al., 2015;Li et al., 2019). However, some reports have been focused on genetic foundation and molecular mechanism of seed developments in mungbean.
Genome-wide association studies (GWASs), along with multi-omics analysis, have been frequently used to mine candidate genes for most important agronomic traits in crops. Integrating GWAS with comparative genomics, transcriptome analysis, and molecular experiments, genes have been identified to be associated with complex traits (Liu et al., 2020c). For example, Gong et al. (2022) conducted a GWAS with high-quality single nucleotide polymorphism (SNP) data and seed-size traits, and found that Cla97C05G104360 and Cla97C05G104380, which are involved in abscisic acid metabolism, played important role in regulating the seed size in watermelon. Duan et al. (2022) identified GmST05 to be associated with soybean seed size through the GWAS of 1800 soybean germplasm resources, and GmST05 differed significantly at the transcriptional level. Liu et al., 2022a,c used GWASs and biological experiments to identify a pleiotropic gene GmPDAT for seed size-and oil-related traits in soybean, and a salt-stress-tolerance gene VrFRO8 in mungbean. Nonetheless, the related genes responsible for seed-size-related traits remained unknown in mungbean.
To address the above issues, 196 mungbean accessions with 3,607,508 SNP markers were used to conduct GWAS for seed length (SL), seed width (SW), 100-seed weight (HSW) using 3VmrMLM (Li et al., 2022b), efficient mixed-model association expedited (EMMAX) (Kang et al., 2010), genome-wide efficient mixed-model association (GEMMA) (Zhou and Stephens, 2012), and compressed mixed linear model (CMLM) (Zhang et al., 2010) methods. Candidate genes around quantitative trait nucleotides (QTNs) and QTN-by-environment interactions (QEIs) for the three traits were predicted by transcriptomics and comparative genomics. Key candidate genes were verified by RT-PCR analysis. Moreover, genes in seed-development-regulation pathway were also mined by comparative genomics. It should be noted that VrEmp24/25 and VrKIX8 were found to be associated with SL and HSW, and a major gene VrPAT14 (LOD = 61.95, r 2 = 5.80%) was identified in QEI detection via 3VmrMLM.

Plant materials and treatments
A diverse set of 196 mungbean accessions including 20 wild and 176 cultivated accessions from 23 countries, were used in this study (Supplementary Data Set 1). All the accessions were planted in a randomized complete block design with two replicates in an experimental field of Kasetsart University, Kamphaeng Saen Campus, Nakhon Pathom, Thailand in 2018 and 2020. In each replicate, each accession was planted in a single row 2.5 m long with 12.5 cm intra-row spacing (ca. 20 plants/row) and 50 cm inter-row spacing. Cultural practices were performed according to Park (1978). SW (mm), SL (mm), and HSW (g) were measured. At maturity. The SL and SW traits for each accession were averaged based on 20 seeds and 100SW for each accession was averaged based on three replicates.

Whole-genome resequencing
The young leaves of the above 196 mungbean accessions were collected 1 week after planting. The DNA was extracted in 2018, using the CTAB method (Smith et al., 2005). Short reads sequenced by an Illumina HiSeq 4000 platform (Illumina, San Diego, CA, United States), and mapped to scaffolds using Burrows-Wheeler-Alignment Tool (BWA) (Version 0.7.15) 1 (Li and Durbin, 2009). Genome Analysis Toolkit (GATK) was used to select SNP and indel 2 (McKenna et al., 2010). Sulv 1 genome was selected as the reference genome in the GATK analysis (Yan et al., 2020). High-quality SNPs and Indel variations were obtained as the following steps. (a) Retaining concordant sites both identified by GATK and VCFtools were retained (Danecek et al., 2011). (b) Filtering out SNP with quality value below 30, removing SNPs with an average coverage depth < 8× and with minor allele frequency (MAF) less than 5%. (c) Deleting insertions and deletions (InDels) with length less than10 bp were deleted. A total of 3,607,508 SNPs were identified.
Genome-wide association study for seed width, seed length, and 100-seed weight Only the SNPs with MAF ≥ 0.05 and missing rate < 10% were used in GWAS (Pongpanich et al., 2010). The lines with more than 95% missing for trait were filtered out (Liaw and Wiener, 2002). SW, SL, and HSW, and the above SNP markers in 196 mungbean accessions were used to conduct GWAS using four different methods, including 3VmrMLM (Li et al., 2022b) via software IIIVmrMLM (Li et al., 2022a), EMMAX (Kang et al., 2010), GEMMA (Zhou and Stephens, 2012), and CMLM (Zhang et al., 2010). The probability threshold for significant QTNs was set at 1/m = 2.77e-07 (m = 3,607,508) for all the four methods (Xu et al., 2018;Zhang Y. M. et al., 2019;Zhang Y. M. et al., 2019), and the LOD score threshold for suggested QTNs was set at LOD ≥ 3.0 for 3VmrMLM (Li et al., 2022b). Heatmaps of the linkage disequilibrium was generated by LDheatmap package (Shin et al., 2006), haplotype analysis was conducted by LDheatmap package (Barrett et al., 2005). The averages for those traits measured in 2018 and 2020 were used in GWAS.

Candidate gene identification
Candidate genes for salt tolerance were mined in the follow steps. (a) All the genes between the 30 Kb around regions for each of the significantly QTN were mined, where the LD-value was about 20 Kb in mungbean, (b) mined the Arabidopsis, rice and soybean homologous genes of those candidate genes, which were reported related to seed developments, seed production, phytohormone signaling pathways and carbohydrate metabolism pathways, etc. (Li et al., 2019), as the candidate genes. (c) The selected genes showing different expression between two groups of mungbean accessions contrasting in seed size (large seed vs. small seed) (see below) were considered as candidate genes.
In the key candidate gene identification, the extracted RNA in two large-seeded accessions at 10 and 25 DAF were treated with RNase-free DNase I (Promega, Madison, WI, United States). After reverse transcription, the cDNA was used as a template for RT-qPCR using the Takara Bio TB Green Premix Ex Taq (Tli RNase H Plus). The detail progress was described by Liu et al. (2022b). Reactions were run on a Bio-Rad CFX96 system. EVM0007380 (homologous of At3g18780) was used as the CK in this experiment. Primers were designed by NCBI and tested by RCR of tubulin. The t-test was adopted in the hypothesis testing, P < 0.05, P < 0.01, and P < 0.001 indicated significant probability levels at 0.05, 0.01, and 0.001, respectively. Information of the primers used is presented in Supplementary Table 1.

Genome-wide association studies for seed-size-related traits in mungbean
Detection of main-effect quantitative trait nucleotides for seed-size-related traits in each environment After removing the SNPs with an average coverage depth < 8× and with a MAF less than 5%, we identified more than 3.6 million SNP markers. In the single-environment analysis, the phenotypic observations for each trait in 196 accessions measured in 2018 and 2020 were used to associate with 3,607,508 SNPs using 3VmrMLM, EMMAX, GEMMA, and CMLM under the situations of five subpopulations and polygenic background control (kinship matrix) (Supplementary Data Set 3). As more than 10,000 QTNs were identified by GEMMA for HSW in 2018, the relevant results were not used in the subsequent analysis. As a result, 208 significant QTNs were identified for the above traits. Thirteen significant QTNs were simultaneously identified in two environments by two GWAS methods (Supplementary Table 3; Supplementary Data Set 4), some significant QTNs are presented in Figure 2. For example, Chr10-25206533-25223155 The frequency distributions of seed-size-related traits. Frequency distributions of HSW (A) (g), SL (B) (mm), and SW (C) (mm) in 196 mungbean accessions, which were measured in 2018 (brown bar) and 2020 (black bar). SD, standard deviation. The associations of HSW with SW and SL, the average dates of those traits measured in 2018 and 2020 were used in the partial correlation analysis (D).

Detection of quantitative trait nucleotides for seed-size-related traits in multiple environments
To detect more stable QTNs, three seed-size-related traits of 196 mungbean accessions measured in 2018 and 2020 were used to associate with 3607508 SNP markers using twoenvironment 3VmrMLM joint analysis. As a result, 32, 33, and 18 significant QTNs were identified for HSW, SL, and SW, respectively (Supplementary Table 3), and had a 1.08% average proportion of total phenotypic variation explained by each QTN. Moreover, eight significant QTNs were identified (Supplementary Table 4). For example, Chr1-8161305-8347626 (LOD = 24.09∼36.33) and Chr10-25222572-25223133 loci (LOD = 29.75∼37.89) were detected to be associated with HSW and SL, respectively (Supplementary Tables 3, 4).
Based on all the above main-effect QTNs in singleand multiple-environment analysis, five stable QTNs across various methods and/or two environments were found (Supplementary Table 5 Detection of quantitative trait nucleotide-by-environment interactions for seed-size-related traits in multiple environments All the above datasets in GWAS were used to detect QEIs using 3VmrMLM. As a result, 5, 10, and 5 significant QEIs were The P1-values were calculated by CMLM, EMMA, and 3VmrMLM, The P2-values were calculated using paired t-test from the average FPKM values at three stages between two high seed weight (n1 = 2) and tow seed weight (n2 = 2) mungbeans, and their significances were marked by * (0.05 level); FC and NA represent fold change and no expression, respectively. found to be associated with HSW, SL, and SW, respectively (Supplementary Figure 2; Table 2). Among these QEIs, 5 had zero dominant-by-environment interaction effects, and 7 had zero additive-by-environment interaction effects. For example, the two loci Chr4-26262890 and Chr4-31677341 for HSW had only additive-by-environment interaction effects of 0.12 ( Supplementary Figures 2A-C In addition, among these QEIs, the QEI locus Chr9-24007163 for SW had large effect, and r 2 was 5.8% (Supplementary Figure 2B, LOD = 61.95). The additive and dominant effects in environment 1 were −0.14 and −0.098, respectively.
Using the same approach described above, among 65 genes around 20 QEIs, four were homologous to previously reported seed development related genes in rice and A. thaliana ( Table 2), although new experiments are necessary to explore these novel GEI-trait associations. The four genes were described as below. VrFATB was linked to the locus Chr4-30176682 (Supplementary Figure 2A). As described in Bonaventure et al. (2003) and Sun et al. (2014), FATB is  The P-values were calculated using paired t-test from the average RPKM values at three stages between two high seed weight (n1 = 2) and tow seed weight (n2 = 2) mungbeans, and their significances were marked by * (0.05 level); FC and NA represent fold change and no expression, respectively.

FIGURE 3
The expression of eight key candidate genes. The expression profiling of eight key candidate genes significantly associated with seed-size-related traits. The expression profiling of eight key candidate genes between two high-seed-weight and two low-seed-weight mungbeans (A). Real-time PCR analysis of the eight key candidate genes; the t-test was used to test the significant differences of genes expression between two high-seed-weight mungbeans at 10 DAF and 25 DAF (B). DAF, days after flowering.
a major determinant of saturated fatty-acid synthesis, and increases FATB activity at low temperature during seedling establishment caused high saturated fatty-acid content in plant. VrGSO1 was linked to the locus Chr4-42563100 (Supplementary Figure 2A). As observed in Creff et al. (2019), GSO1 was a stress signal-pathway-related gene, and stress-associated MPK6 protein acted downstream of GSO1 in developing embryo. VrPAT14 was linked to the locus Chr9-24007163 (Supplementary Figure 2B). In Zhao et al. (2016), PAT14 was involved with NPR1-dependent salicylicacid signaling. VrLACS2 was linked to the locus Chr8-32848165 (Supplementary Figure 2C), in which VrLACS2 was essential for normal cuticle development in Arabidopsis (Schnurr et al., 2004) and CrLACS2 suppression resulted in 50% less oil, yet with a higher amount of chloroplast lipids under N-deprivation (Bai et al., 2022).

Haplotype analysis of the main candidate genes
Two DEGs, VrEmp24/25 and VrKIX8, were detected in the single-and multi-environment analyses (Figures 4A,B), and verified by RT-qPCR. Their haplotypic analyses were described as below.
In the haplotype analysis of VrEmp24/25, five SNP markers were found to be within VrEmp24/25 and the promoter region (Supplementary Data Set 8), and the two SNP markers in VrEmp24/25 were used to consist of three haplotypes ( Figure 4D). Among the three haplotypes, hap 1 (5.17 g) had significantly higher HSW than hap 2 (1.58 g) and hap 3 (4.50 g; P = 2.11E-29) (Supplementary Table 7). Thus, hap 1 is elite haplotype. And the elite haplotypes TT made up more than 90.9% (160/176) in the cultivated mungbeans. VrEmp24/25 with elite haplotype frequencies less than 45% in wild mungbeans (Supplementary Table 7; Figure 4) can be exploited for the improvement of mungbean cultivars.
Based on these results, we deduced that these two SNP and six SNP cause the difference expression of the VrEmp24/25 and VrKIX8 gene, respectively. The discovery of VrEmp24/25 and VrKIX8 two domestication/improvement genes can accelerate breeding selections and facilitate ideal crop designs.

Expression patterns of seed development pathway genes in mungbean
As seed development pathway genes were largely unknown in mungbean, we mined seed development pathway genes by comparative genomics and transcriptomics analysis. As a result, 54 genes in seed-development pathway were identified in this study (Figure 6; Supplementary Data Set 9), such as two GPA1, one AGB, and one AGG3. In the ubiquitin proteasome pathways, two DA1, one DA2, one SOD2, one EOD1, and one UBP15 rather than SAMBA were identified. In the auxin pathways, two ABA2, one ABI5, three SHB1, five IKU2, and three CKX2 rather than IKU1 and MINI3 were identified (Figure 6A). Five transcription factors including three BES1, and two SOD7 were identified. Moreover, 16 genes for seed size developments were found to be with uncertain pathways, including three KIX8, five MES1, and one KLU (Figure 6A; Supplementary Data Set 9). Among the 54 genes, 13 genes were significantly differentially expressed (P-value < 0.05, t-test) between two low-seed-weight (nos. G169 and G171) and two high-seed-weight (no. G141 and G143) accessions in the 196 mungbean accessions using the transcriptome data at 10, 15, and 25 DAF (Figure 6B; Supplementary Data Set 8). Moreover, almost 90% of the 54 genes (48/54) had higher expressions in the early stage of seed development (10 and 15 DAF) than in the late maturation stage (25 DAF), including VrKIX8 (EVM0032114), which was commonly identified in the GWAS by 3VmrMLM for HSW and SL. And EVM0010067/VrABA2, EVM0033315/VrSHB1, EVM0028440/VrABI5, and EVM0030447/VrIKU2 were also identified in the GWAS by 3VmrMLM, within 100 Kb region of significant QTNs (Table 1).
We also did the PPI analysis among the seed development pathway genes, and found five pairs of PPIs were larger than the medium confidence value of 0.40 (Supplementary  Table 7), indicating the existence of significant PPIs, i.e., EVM0013794.1 (VrAGG3) and EVM0006667.1 (VrDA2) Genetic analysis of VrKIX8. Local Manhattan plots for HSW in multi-environments. LOD ≥ 3.0 for the 3VmrMLM as the significant QTN (A). LD heatmaps surrounding Chr1-8161305-8347626 loci (B). Genes around the significant QTN region, shown at the bottom (C). Haplotype analysis of VrKIX8 (D). Wil, the wild accessions. Cul, the cultivated accessions. symbol "**" means omit the same sequence part.

Discussion
The high-yield and efficiency breeding progress of mungbeans have been limited by the lack of ideal yieldrelated genes. At present, few QTNs or QTLs of yield-related traits in mungbeans have been reported (Kang et al., 2014). This study provided a genetic analysis of seed-size-related traits in mungbeans, to improve the accuracy of significant QTNs, we used multiple genome-wide M0017 association studies combined with multi-omics analysis to mine candidate genes associated with yield-related traits. Firstly, a total of 98 QTNs and 20 QEIs were identified using 3VmrMLM, while 95 and 15 QTNs were identified using EMMAX, and CMLM, respectively. Then, in the identification of candidate genes, 12 key candidate genes were mined, and seven of them including VrKIX8, VrEmp24/25, and VrPAT14 were evidenced by transcriptome analysis and RT-qPCR analysis. Lastly, through haplotype analysis, the thirtieth amino acid of VrEmp24/25 in the elite haplotype was changed from Ile to Phe. And there were six SNP in the promoter and 5 UTRs of VrKIX8, however, the amino acid sequence of VrKIX8 in the elite haplotype was not changed. The results provided the theoretical basis for both the functional identification of seed-size-related genes and for quality improvements in mungbean breeding.
Multiple genome-wide association studies methods combined with multi-omics analysis in mining candidate genes In the GWAS, how to identify candidate genes around significant QTNs has been a challenge. Liu et al. (2020c), Zhang et al. (2021), and Gong et al. (2022) selected the 100kb interval upstream and downstream of the significant QTN as the candidate interval in watermelon and soybeans. Usually, the interval has been chosen according to the LD decay values.
Notably, 3VmrMLM showed more powerful ability in the detection of significant QTN than GEMMA, EMMAX, and CMLM, as it found more differentially expressed key candidate genes than other methods. The combination of 3VmrMLM and multi-omics analysis in the genetic analysis of complex traits was helpful.
Besides, in maize, the loss function of EMP24 and Emp25 would impair embryo and endosperm development (Xiu et al., 2020). EMP24 was required for the splicing of nad4 (Ren et al., 2019), and the lack of either Nad4 or Nad5 blocked the assembly of complex I holoenzyme in Arabidopsis (Ligas et al., 2019). The loss of the steady-state level of mitochondrial nad5 mature mRNA blocked the assembly of complex I and caused an arrest in endosperm development . Lastly, the elite haplotypes of VrEmp24/25 (TT) made up the main proportion of more than 90.9% in cultivated mungbeans, 45% in wild mungbeans ( Figure 4E). The HSW in hap 1 haplotypes accessions was significantly higher than that in hap 2 and hap 3 (P = 2.11E-29). It was reported that a single amino acid completely prevented the appearance of the enzyme in the medium, and we inferred that the related variation could lead to the change in enzyme activity (East et al., 1990;Alfson et al., 2018).
There have four evidences to take VrKIX8 as another important seed-size trait gene. Firstly, VrKIX8 associated with Chr1-8161305-8347626 (LOD = 24.09∼36.33) for HSW and SL were detected in multi-environment by 3VmrMLM ( Figure 5A; Supplementary Table 5). Secondly, VrKIX8 (LOD = 24.09∼36.33) had significantly differentially expressed between high-and low-HSW accessions ( Figure 3A). Then, in Arabidopsis, the disruption of KIX8/9 and PPD1/2 could cause large seeds due to increased cell proliferation and cell elongation in the integuments (Liu et al., 2020a). In soybeans, the loss of the function GmKIX8-1 showed a significant increase in the size of seeds and leaves. In addition, the increase in organ size was due to the increased cell proliferation, rather than cell expansion. GmKIX8-1 showed negatively regulated cell proliferation in plants (Nguyen et al., 2021). Lastly, the elite haplotypes of VrKIX8 (ATCGAA) made up the main proportion of more than 73% in cultivated mungbeans, 40% in wild mungbeans. Moreover, there are four SNPs in the promoter and of VrKIX8, and two SNPs in the CDS region, however the amino acid sequence did not change between the elite haplotypes and the other haplotypes ( Figure 5C). The HSW in hap 1 haplotypes accessions was higher than that in hap 2 to hap 4 (P = 1.19E-21). We supposed that the mutations may have influenced the translation efficiency of VrKIX8 and caused low expression in cultivated accessions during mungbean domestication.
We also analyze the possible correlation between the main seed development pathways. Among the 54 genes, five genes (VrAGG, VrDA2, VrAGB, VrGPA1-1, and VrGPA1-2) consisted of five pairs of significant PPIs. Interestingly, four pairs PPIs were found to be in the G-protein-signaling pathway, and one pair of PPIs was found to be in the G-protein-signaling and the ubiquitin proteasome pathways (Figure 6; Supplementary Table 6). Ubiquitin proteasome pathway is an important pathway for the selective degradation of proteins and seed development (Smalle and Vierstra, 2004), and the G-protein-signaling pathway is a ubiquitous cell transmembrane signal transduction pathway in eukaryotes (Huang et al., 2006). Moreover, mutations in GPA1 or AGB1 could cause short flowers (Lease et al., 2001;Ullah et al., 2001). The overexpression of AGG3 promoted seed and organ growth by increasing cell proliferation, and loss-offunction mutations in AGG3 caused small seeds and organs (Chakravorty et al., 2011;Li et al., 2012). The ubiquitin receptor DA1 could control seed size by restricting cell proliferation in maternal integuments (Li et al., 2008). DA1 functioned synergistically with DA2 to restrict seed growth, and DA2 physically interacted with DA1 in vitro and in vivo (Song et al., 2007;Xia et al., 2013). This interaction could mediate the interactions between the G-protein-signaling pathway and the ubiquitin proteasome pathway, which might offer an important clue in the mechanism analysis of seed development.
In addition, 48 genes had higher expressions in the early stage of seed development than in the late maturation stage of seed development, indicating that seed-development-related genes function primarily in the early stages of seed development, which was consistent with the findings of Zuo et al. (2022) in soybean.

Conclusion
This study conducted GWAS for seed-size-related traits in mungbeans. 98 QTNs and 20 QEIs were identified using 3VmrMLM, while 95, >10,000, and 15 QTNs were identified using EMMAX, GEMMA, and CMLM, respectively. A total of 12 key candidate genes were mined, which were homologous to known seed-development genes in rice and A. thaliana. VrEmp24/25 and VrKIX8 were identified as main candidate genes around two stable QTNs, the two candidate genes were further confirmed by RT-qPCR and haplotype analysis, and prevalent haplotypes of VrEmp24/25 and VrKIX8 may be useful in mungbean breeding.

Data availability statement
The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found below: The WGS sequencing data of 196 mungbean accessions was uploaded to NGDC, with subCRA011538, subSAM100395, and PRJCA010704 ID.