GWAS, MWAS and mGWAS provide insights into precision agriculture based on genotype-dependent microbial effects in foxtail millet

The study investigates how host genetic variation and root-associated microbiota jointly shape complex agronomic traits in foxtail millet. While GWAS has identified genotype-phenotype associations in many crops, genetic variants alone often fail to explain quantitative variation in growth and yield. Root-associated microbes influence plant development, stress tolerance, and defense, and their composition is modulated by host genotype. In foxtail millet, key yield components correlate with growth indices, but loci for growth or yield remain poorly characterized. The authors aim to integrate GWAS of plant traits with MWAS of the rhizoplane microbiota and mGWAS of microbe abundance to map a genotype–microbiota–phenotype network, identify microbial biomarkers of agronomic traits, uncover host genetic determinants of microbiome assembly, and experimentally validate genotype-dependent microbial effects on plant growth.

Prior crop GWAS have elucidated trait loci in maize, rice, sorghum, cotton, and soybean, but isolating genes controlling complex yield traits remains challenging. Plants engage root microbiota that confer nutrient acquisition and stress resilience; host genotype affects microbiome composition, providing an indirect path to modulate phenotype via environment. Specific factors such as plant immunity, nutrient stress, and specialized metabolites influence plant–microbe associations. Beneficial microbial effects are often cultivar-specific, complicating broad application. For foxtail millet, previous work profiled root zone microbiota and its correlation with yield. Few plant MWAS exist, though analogous approaches in humans linked microbiota to complex traits, and mGWAS in Arabidopsis and sorghum revealed host loci shaping microbiomes. These contexts motivate integrating host genetics and microbiome data to understand and predict agronomic traits.

Study design and samples: 1127 foxtail millet cultivars were genotyped; after QC and population structure filtering, 827 cultivars from China were retained for analyses. All 827 were grown in a field trial (Yangling, China) and 12 agronomic traits were measured: growth traits (TSLL, TSLW, MSH, MSW, MSPD, FNL) and yield traits (MSPL, PGW, MSPW, HKW, MSSN, SGN). Genotyping and SNP processing: RAD-seq was performed. Reads were quality filtered (SOAPnuke), aligned to the Setaria italica cv. Zhanggu v2.3 reference (BWA-MEM), processed with SAMtools and Picard, realigned around indels (GATK), variants called with HaplotypeCaller and filtered (coverage and quality thresholds), biallelic SNPs retained, imputed (BEAGLE), and additional QC applied (PLINK) including IBD filtering, MDS outlier removal, MAF ≥ 0.01, and call-rate thresholds. 161,562 high-quality SNPs were retained. LD half-decay distance was ~9 kb. GWAS of plant traits: Linear mixed models (GEMMA) with kinship and top 10 PCs as covariates were used. The suggestive genome-wide significance cutoff was set to 2.01e-5 based on 49,512 effective SNPs. Candidate genes were defined within 20 kb or high-LD (R^2>0.4) windows around significant SNPs; GO/KEGG enrichment was performed. Root microbiome data and MWAS modeling: Rhizoplane 16S rRNA V4-V5 data yielded an OTU table. A common sub-community of 1004 OTUs (≥70% occurrence) representing on average 61.3% abundance was defined and evaluated for stability (AVD index) and network topology. For each trait, OTUs significantly correlated (adjusted P<0.05) were selected as candidate microbial predictors. Linear regression models predicting traits from either SNPs or OTUs included 10 genetic PCs as covariates. Five-fold cross-validation was repeated 30 times; models were selected by AIC on training sets and evaluated by R^2 in testing sets. Mixed models combining SNPs and OTUs and their interactions were also fit similarly. Marker OTUs were defined from the best-performing combined models. Network analysis and OTU categorization: OTUs were categorized as abundant (AT), moderate (MT), or rare (RT) by prevalence and relative abundance thresholds. Co-occurrence networks were built on Spearman correlations (|r|>0.4, adjusted P<0.05) and centrality metrics computed (igraph). Hub OTUs had degree>400 and closeness>0.5. Heritability of OTUs: OTU counts were CSS-normalized, log-transformed, and SNP-based heritability estimated using GCTA GREML; BH adjustment across 1004 OTUs identified highly heritable microbes. mGWAS of microbial abundance: OTU abundances were rank-based inverse normal transformed (Blom method) and tested against SNPs using GEMMA LMM with kinship and 10 PCs. Associations with P<2.01e-5 were considered suggestive. Candidate host genes were mapped and enriched (GO/KEGG). Bacterial isolation and strain–phenotype assays: Roots from 2000 cultivars were processed to isolate culturable bacteria (NA medium). 644 strains were obtained; 257 with distinct full-length 16S sequences were retained. Strains were mapped to marker OTUs at ≥97% 16S identity to identify 24 representative strains. Ten strains (six positive and four negative markers) with top regression coefficients were chosen for validation. Germ-free foxtail millet (cv. Huagu12) seeds were inoculated on sterile plates (2 ml OD600=0.5 culture normalized) or grown in sterilized field soil, with water controls. After 7 days (plates) or 14 days post-inoculation (soil), shoot height and root length were measured (t-tests with BH adjustment for soil assays). Genotype-dependent effects: Reference and allele cultivars at FaQR and SUZ12 loci were inoculated with the corresponding marker strains (Kita594 for FaQR; Baci173 for SUZ12). Differences from controls were analyzed by one-way ANOVA with LSD test; Wilcoxon tests compared effects between genotype classes; PERMANOVA (adonis) quantified genotype and strain effects. Transcriptomics: Seedlings inoculated with selected strains (Kita594, Baci299, Acid550, Baci173) and controls (sterilized soil experiment) were pooled (n=5 per treatment) for RNA-seq (BGISEQ, PE100). Reads were filtered (SOAPfilter), aligned (hisat2) to the Setaria reference, counted (featureCounts), quantified (FPKM; ballgown), and DEGs identified (DESeq). Pathway enrichment used Fisher's exact test. Data and code availability: Genomic and RNA-seq data are in NCBI PRJNA873890 and CNGBdb CNP0001521; reference genome CNPhis0000549; analysis code archived at https://zenodo.org/badge/latestdoi/424864991.

• Population and genotyping: 827 cultivars analyzed; 161,562 high-quality SNPs post-QC; LD half-decay at ~9 kb; three phylogenetic groups identified.
• Trait heritability and GWAS: Of 12 traits, 11 were significantly heritable; HKW showed H^2=0.006 (P=0.15). Growth traits were more heritable (e.g., MSPD H^2=0.46) than yield traits (e.g., PGW H^2=0.16). GWAS identified 86 significant SNP loci and 91 associations across 10 traits (P<2.01e-5), with notable clusters on chromosomes 2 (15 SNPs), 4 (16), 6 (11), 7 (10), and 9 (16). Candidate genes near signals included regulators of growth/development (ATG8C, ERF1B/ERF1, PRR37, Cyclin-like F-box), drought/stress response (PP2C, ARR12, NPF1.2/NPF4.6, WDR26, CPK2a, PIP5K1, SAPK9, APX, thioredoxin fold proteins), immunity/defense (RPP13, RGA2, RPS2, LRR-RLKs, EF-Tu, SYP22, NB-ARC, WAK2), and nutrient uptake (IRT1/IRT2, phosphate transporters).
• Predictive modeling (MWAS framework): Using SNP-only predictors, average explained variance in test sets was 32.82% (TSLW), 28.55% (MSPD), 47.27% (MSW), 15.02% (MSPW), 38.89% (PGW), 64.60% (MSPL). Using microbiota-only predictors (1004 common OTUs filtered to trait-associated), test-set R^2 averaged 32.47% (TSLW), 17.43% (MSPD), 56.06% (MSW), 30.36% (MSPW), 35.17% (PGW), 12.61% (MSPL). Combined SNP+OTU models significantly improved performance for all traits except MSPL, with average R^2 of 46.50% (TSLW), 59.08% (MSPD), 65.69% (MSW), 38.45% (MSPW), 43.04% (PGW), 44.31% (MSPL). Best retained models achieved 53.42%, 63.73%, 70.54%, 50.16%, 55.88%, and 54.82% for TSLW, MSPD, MSW, MSPW, PGW, and MSPL, respectively. Microbiota contributed more strongly than genotype to MSPD and MSPW.
• Marker microbes: 257 rhizoplane marker OTUs were identified across 86 genera and 15 phyla; top phyla were Proteobacteria (68 OTUs), Actinobacteria (54), Bacteroidetes (36), Acidobacteria (35), Firmicutes (33). 145 markers associated with growth traits and 128 with yield traits; 17 markers were shared between growth and yield; none was shared across all six focal traits. Network analysis showed abundant marker OTUs exhibited higher degree, closeness, and betweenness centrality than moderate or rare markers (P<0.05).
• Experimental validation: Ten representative strains (six positive, four negative markers) corresponding to marker OTUs were isolated and tested. Multiple strains promoted or suppressed seedling growth on plates and in sterilized soil. Growth-promoting strains induced an expansin gene linked to cell wall loosening and shoot/root growth. The suppressing strain Baci173 uniquely induced 39 genes, including auxin biosynthesis/response components (amidase, ALDH), auxin-responsive IAA proteins, L-glutamine synthetase, and BCAT, consistent with root growth inhibition.
• mGWAS and host–microbe genetics: mGWAS identified 2,108 significant SNP loci associated with the abundance of 838 OTUs. Genes near associated SNPs were enriched in plant–pathogen interaction (ko04626), MAPK signaling (ko04016), steroid biosynthesis (ko00100), among others, including PRRs (e.g., FLS2-related), R genes (RPM1, RPS2), PTI regulators (PTI1, PTI6), CALM, and WRKY25, suggesting defense genes shape root microbiota. Among 257 marker OTUs, 219 were associated with host genetic variation; 77 were highly heritable (H^2≈0.15–0.32). Notable associations included immune-related FLS2, bHLH35, and WAK2 with Acidobacteria markers, and nutrient/metabolite/stress genes (e.g., magnesium transporter, achilleol B synthase, BGLU12, RXW8) with Actinobacteria, Bacteroidetes, and Proteobacteria markers. Genes like SUZ12 and WAT1 were linked to Firmicutes and Proteobacteria markers with opposite growth effects. Positive marker Acidovorax OTU_46 associated with an EREBP-like factor; Kitasatospora OTU_8 associated with FaQR. RNA-seq showed FaQR, VWA, SUZ12, and an EREBP-like factor were differentially expressed upon inoculation with corresponding marker strains.
• Genotype-dependent microbial effects on phenotype: At the FaQR locus, the growth-promoting strain Kita594 enhanced shoot growth significantly in FaQR allele cultivars but not in reference cultivars (ANOVA-LSD, adjusted P<0.05). At SUZ12, the suppressing strain Baci173 inhibited root growth, with stronger suppression in allele cultivars versus reference cultivars. PERMANOVA detected significant genotype-by-strain effects: Kita594 genotypes R^2=13.048, P<0.001; Baci173 genotypes R^2=0.07, P<0.001. Multiple ANOVAs for plant height and root length across genotypes and treatments were highly significant (e.g., F up to 19.73, adjusted P<2e-16).

Integrating GWAS with MWAS and mGWAS reveals that root microbiota, shaped in part by host genetics, substantially contribute to variation in key agronomic traits. For some traits (e.g., MSPD, MSPW), microbiota explained more variance than genotype alone, while for others (e.g., MSPL), genotype dominated, consistent with known loci such as PRR37 influencing panicle length. The identification of 257 marker OTUs and their validation through culture-based assays demonstrates that MWAS can discover functionally relevant rhizobacteria that promote or suppress growth. Transcriptomic responses to marker strains implicate hormone biosynthesis/signaling, immune pathways, and nutrient uptake, aligning with the genetic enrichment of plant–pathogen interaction and MAPK signaling in mGWAS. The observation that 219 of 257 marker OTUs associate with host genetic variants, and that marker effects are genotype-dependent at loci such as FaQR and SUZ12, supports a causal chain in which host genes shape the microbiota, which in turn modulate plant phenotypes. This genotype–microbiome interplay suggests new strategies for precision microbiome management and breeding, tailoring microbial inoculants to specific host genotypes to maximize agronomic benefits.

This work maps a reciprocal network among host genetic variation, rhizoplane microbiota, and agronomic traits in foxtail millet by combining GWAS, MWAS, and mGWAS across 827 cultivars. It identifies 257 microbial biomarkers predictive of growth and yield traits, uncovers host genes (notably immunity- and signaling-related) associated with microbial assembly, and experimentally validates genotype-dependent effects of marker strains on plant growth. The combined SNP+microbiome models substantially improve trait prediction beyond genetics alone, highlighting the importance of including microbiome data in crop improvement pipelines. Future work should validate newly implicated host genes and microbial taxa mechanistically, expand to multi-environment trials to resolve genotype-by-environment effects on microbiome assembly, and develop precision microbiome engineering approaches that match beneficial consortia to host genotypes for sustainable yield gains.

Analyses were conducted in a single field environment, which may limit generalizability and does not fully capture genotype-by-environment interactions; the authors note that multi-environment GWAS/mGWAS will be needed. The SNP density was relatively low, leading to broader candidate gene windows (±20 kb or high-LD regions) and potentially reduced mapping resolution. Functional validation of many candidate host genes and most microbial markers remains outstanding; only a subset of marker OTUs had culturable representatives tested. Validation assays focused on early seedling growth in plates and sterilized soil rather than field yield outcomes. Direct GWAS associations between SNPs linked to marker OTUs and plant traits were often below suggestive significance, indicating partially independent or mediated effects.