Wild rice locus enables microbiome re-domestication for enhanced nitrogen-use efficiency

Ensuring sustainable food production while reducing environmental impacts is a major challenge¹. Intensive fertilizer and pesticide use has increased crop yields but also caused environmental degradation^2,3,4. Sustainable strategies are therefore needed to improve crop nutrient-use efficiency (NUE) and yield. Enriching beneficial plant-associated microbes can enhance nutrient uptake and stress resilience, reducing reliance on synthetic fertilizers and pesticides^5,6.

Extending the microbiome rewilding hypothesis⁷, which proposes restoring microbial partners lost during domestication, we asked whether this process is genetically controlled by host loci. We hypothesized that beneficial microbiota and their functions can be restored in cultivated crops through the inheritance of host genetic loci from wild relatives, a process we term ‘microbiome re-domestication’. We therefore analysed wild rice microbiomes to identify key ‘ancestral’ microbiota and corresponding microbiome-associated genes (M-genes), defined as host genes involved in microbiome recruitment and function that were counterselected during domestication and breeding.

Modern varieties were primarily bred for above-ground agronomic traits under high-input conditions, potentially overlooking beneficial below-ground microbiota adapted to low-input or stress environments. Growing evidence shows that microbiome selection differs between wild accessions and cultivated plants^8,9,10. Studies using tomato (Solanum lycopersicum) and barley (Hordeum vulgare) recombinant inbred line (RIL) populations demonstrated associations between plant loci and rhizosphere microbiome composition^11,12,13. Despite these advances, previous studies have focused on statistical associations between host loci and microbiome composition under controlled conditions, without addressing whether beneficial microbiome member(s) from wild relatives can be stably reinstated in crop progeny or contribute to agronomic traits under real field conditions. In particular, it remains unclear whether core microbiome members can be retained following hybridization and genetic recombination in RIL populations, and whether specific host loci can restore the functions of these microbes associated with NUE.

We previously showed that wild rice (Oryza rufipogon) hosts a core rhizobacterial community that promotes nutrient acquisition under nutrient limitation and provides a valuable resource for M-gene-based crop improvement¹⁴. We further demonstrated that wild-rice rhizosphere exhibits higher abundance and activity of nitrogen fixation gene, whereas cultivated rice is enriched in genes for N₂O production and emissions¹⁵. This differentiation is partially reproduced by root exudates¹⁵, linking wild rice genetics and root chemistry to microbial recruitment supporting nitrogen acquisition and reduced greenhouse gas emissions.

In this study, we investigated whether beneficial rhizosphere microbiota lost during rice domestication can be genetically recovered in hybrid progeny through the introgression of wild-rice loci. We first assessed the rhizosphere microbiota assembly in two RIL populations derived from wild rice (O. rufipogon) and cultivated rice (O. sativa L. ssp. japonica and indica) and determined whether paternal wild-rice core microbiota could be reinstated in RILs grown under field conditions (Fig. 1). We then mapped quantitative trait loci (QTLs) associated with core taxa and overlapped microbiome-associated QTLs with domestication loci identified from whole-genome data of 382 rice accessions. Subsequent transcriptome analyses and microbiome profiling identified candidate genes, particularly the receptor-like kinase-encoding gene OsRLK (Os04g0659300), whose association with an Acinetobacter synthetic community (SynCom) was functionally validated using overexpression and knockout rice lines.

Collectively, by integrating genetics, genomics and rhizosphere microbiome profiling, we link OsRLK to the reinstatement of wild-rice-derived bacteria that improve nitrogen-related functions and support enhanced plant performance in domesticated rice. This study goes beyond previous association-based studies by connecting host genetic loci to functionally relevant core microbiota and nitrogen-use efficiency under field conditions.

Rhizobacterial diversity and core microbiome identification

To characterize the core rhizosphere bacteria of paternal wild rice, we performed 16S ribosomal (r)RNA amplicon sequencing on 424 samples, including bulk soil (n = 6), parental rhizosphere samples (n = 6 per genotype) and 2 RIL populations (RIL-J and RIL-I), each comprising 100 lines and 200 rhizosphere samples in total (2 replicates per line). Rhizosphere communities differed in alpha diversity (Chao1 index) and composition from bulk soil (Extended Data Fig. 1, permutational multivariate analysis of variance (PERMANOVA), R² = 0.17, P = 0.001).

To define the paternal wild-rice-specific core rhizosphere bacteria, we first identified core bacterial taxa in the wild rice, japonica and indica parental lines using generalized joint attribute modelling (GJAM) analyses and then excluded core bacterial taxa shared with the maternal japonica or indica parents. The remaining bacterial taxa were defined as paternal wild-rice-specific core members and tracked in the RIL populations (Extended Data Fig. 2). The RIL-J core comprised 45 amplicon sequence variants (ASVs), while the RIL-I core comprised 48 ASVs (Extended Data Fig. 2a). Of these, 27.40% (n = 20) of ASVs were shared, while 34.25% (25) and 38.35% (n = 28) were unique to RIL-J and RIL-I, respectively (Extended Data Fig. 2b).

The 20 core bacterial taxa (ASVs) were assigned to Acinetobacter, Asticcacaulis, Azospirillum, Bacillus, Gaiella, Kaistia, Lysobacter, Methylomonas, Nitrospira, Nocardioides, Novosphingobium, Paenarthrobacter, Pantoea, Pseudogulbenkiania, Pseudolabrys, Pseudomonas, Roseomonas, Sideroxydans, Streptomyces and Vogesella, and were present in different frequencies within the 2 RIL populations relative to the paternal wild rice (Fig. 2). To assess the influence of paternal wild rice on the 2 RIL populations, we calculated the percentage of each core taxon across 200 RIL accessions. Gaiella occurred in 88.5% of RIL-J accessions and 96.6% of RIL-I accessions (Fig. 2). Acinetobacter, Paenarthrobacter and Duganella reached 60.5%, 64.6% and 80%, respectively, in RIL-J while Bryobacter, Lysobacter, Anaeromyxobacter, reached 62%, 75.2% and 60.7%, respectively, in RIL-I, demonstrating that paternal core bacterial taxa are maintained and detectable in the RIL progeny.

Associations between core bacterial taxa and rice phenotypes

Partial Mantel tests were performed to assess correlations between the abundances of the core bacterial taxa and the phenotypes of the two RIL populations of rice, including plant height, fresh and dry shoot weights, fresh and dry root weights, main spike length, ground diameter, chlorophyll a (Ca), chlorophyll b (Cb), total chlorophyll (CT) and total nitrogen (TN) content (Extended Data Fig. 3).

In RIL-J, plant ground diameter and TN significantly correlated with the core bacterial taxa, while shoot and root biomass and CT did not (Extended Data Fig. 3). Conversely, in RIL-I, fresh and dry root weights, Ca, main spike length, CT and Cb significantly correlated with the core bacterial taxa (Extended Data Fig. 3).

Genetic loci associated with core bacterial taxa in rice

Linkage analysis of the 20 core ASVs common to wild, RIL-J and RIL-I populations identified 31 QTL peaks (logarithm of odds (LOD) > 2) associated with the abundance of 17 core ASVs. In RIL-J, 17 QTLs were associated with 13 core ASVs, while in RIL-I, 14 QTLs were linked to 12 core ASVs (Fig. 3 and Supplementary Table 1). Using wild-rice allele as reference, negative effects were relatively more associated with the cultivated rice allele and positive effects with the wild-rice allele. Of the lineages with multiple QTLs, Acinetobacter and Vogesella were wild allele enriched (Fig. 3b). The chromosomal locations associated with each core bacterial taxon differed between RIL-J and RIL-I. For example, Acinetobacter was linked to a QTL on chromosome (Chr) 4 in RIL-J but to a QTL on Chr 2 in RIL-I, while Sideroxydans was associated with a QTL on Chr 8 in RIL-J but with 2 QTLs on Chr 2 and 5 in RIL-I (Fig. 3c). Together, these analyses indicate that wild-rice loci consistently influence the abundance of specific core bacterial taxa, with Acinetobacter representing one of the microbial taxa across populations.

Rice NUE genes linked to core bacterial functions

Given the urgent need to improve NUE in rice systems, where only ~30–40% of applied nitrogen is utilized by plants, we focused on microbial functions related to nutrient cycling. Previous studies showed that nitrogen-cycling-related microbes are enriched in wild rice^14,15; however, their genetic basis in host control remains unclear. Using FAPROTAX¹⁶ for mapping taxonomy to validated functions, we inferred functions of the 20 core ASVs shared by wild, RIL-J and RIL-I. Ten ASVs in nitrogen metabolism were associated with genes on Chr 3, 4, 5, 8, 9, 10 in RIL-J and Chr 1, 2 and 6 in RIL-I (Extended Data Fig. 4). Seven ASVs including Acinetobacter (ASV6549), Novosphingobium (ASV6074), Sideroxydans (ASV6368), Streptomyces (ASV6088), Pantoea (ASV6098), Bacillus (ASV5438) and Pseudomonas (ASV6093) were predicted to be nitrogen fixers (Supplementary Table 2).

We identified 10 ASVs predicted to be involved in nitrogen metabolism. Among them, one ASV was associated with the previously reported rice NUE gene OsAMT1;3, which encodes an ammonium transporter¹⁷. Specifically, the gene was associated with Sideroxydans (ASV6368) in RIL-I (Supplementary Table 2). Notably, the remaining ASVs were linked to unreported NUE-related rice genes. For instance, 113 rice genes, including OsRLK (Os04g0659300) were significantly associated with Acinetobacter (ASV6549) abundance (Supplementary Table 3).

QTLs associated with core bacterial taxa overlap rice domestication loci

To explore whether the host and rhizosphere core bacterial taxa co-evolved, we sequenced the whole genomes of 382 rice accessions, including 117 wild-rice O. rufipogon accessions and 265 cultivated rice O. sativa accessions (Extended Data Fig. 5 and Supplementary Table 4). We estimate genetic diversity ratio (π_wild/π_cultivar) between wild and cultivated rice accessions and identify rice domestication loci (Fig. 4a,b). This analysis identified 21 domestication loci across all chromosomes except Chr 10 and 11 (Supplementary Table 5). Comparing QTLs for the 20 core ASVs shared by wild, RIL-J and RIL-I populations with these loci revealed two strong selective sweeps overlapping with QTLs associated with Acinetobacter ASV6549 and Pantoea ASV6098 (Fig. 4c,d). The locus associated with Acinetobacter showed higher additive effect value (ADD), phenotypic variance explained (PVE), along with stronger selective sweep signals (Supplementary Table 2), and was therefore selected for isolation and subsequent bioassays to assess their potential contribution to nitrogen acquisition.

Rice candidate genes associated with Acinetobacter

Nine Acinetobacter strains exhibiting high 16S rRNA gene similarity to A. baumannii, A. halotolerans, A. sichuanensis, A. schindleri, A. geminorum, A. defluvii and A. lactucae were isolated from wild-rice rhizosphere soil to make a SynCom (Extended Data Fig. 6c). Notably, nifH and nitrogenase protein abundance in Acinetobacter 7 was significantly higher than in the other 8 Acinetobacter strains, indicating potential functional heterogeneity in nitrogen-related traits among the strains (Extended Data Fig. 6a,b).

To identify plant genes associated with Acinetobacter, we compared root transcriptomes of paternal wild rice grown in non-sterilized soil (WN) versus sterilized soil (W), paternal wild rice grown in sterilized soil inoculated with the Acinetobacter SynCom (WA) versus sterilized soil (W). Kyoto Encyclopedia of Genes and Genomes (KEGG) comparative analysis revealed significant pathways ordered by proportion of differentially expressed genes (DEGs). In WN versus W, the DEGs were enriched in betalain biosynthesis, flavonoid biosynthesis and nitrogen metabolism (Fig. 5a), while in WA vs W, enriched pathways comprised nitrogen metabolism and biosynthesis of various plant secondary metabolites (Fig. 5b). Nitrogen metabolism was the only pathway common to both comparisons. Fourteen reported NUE-related DEGs, including OsNRT2.3a, NRT2.2, OsNAR2.1, OsNPF7.9, OsDREB1C and OsAMT1;3 were identified (Fig. 5c and Supplementary Table 6).

A total of 2,494 upregulated and 2,661 downregulated genes were shared between WN vs W and WA vs W comparisons (Fig. 5d). Within the Acinetobacter-associated 0.64 Mbp region on Chr 4, there were 16 upregulated genes and 3 downregulated genes. Among these, only the receptor-like kinase-encoding gene OsRLK in paternal wild rice met the criteria of high expression with log₂(fold change) > 2 (q-value < 0.05) in both comparisons (Fig. 5e).

To validate these findings, OsRLK from wild rice was cloned to generate overexpression lines for testing interactions with the Acinetobacter synthetic community. Southern blot analysis confirmed one or two copies of the OsRLK expression cassette in six transformants (Y165, Y166, Y168, Y172, Y174 and Y175; Supplementary Fig. 1), and two single-copy lines were selected for subsequent functional analyses.

Interaction of Acinetobacter SynCom with OsRLK rice lines

Wild-type indica rice (1504), OsRLK-overexpressing lines, OsRLK-knockout lines and wild rice, inoculated or not with Acinetobacter SynCom and inoculated with a heat-killed Acinetobacter SynCom were evaluated. Total nitrogen (TN) levels were significantly higher in SynCom-inoculated wild rice and OsRLK-overexpressing lines, whereas no significant differences occurred among these rice lines without and with a heat-killed SynCom (Fig. 6a, P < 0.05). Gene expression was significantly higher in SynCom-inoculated transgenic and wild rice compared with non-inoculated and with heat-killed SynCom plants (Fig. 6b). Furthermore, OsRLK expression was consistently higher in transgenic rice than in wild-type rice (1504) under all treatments. Similarly, plant dry biomass increased in SynCom-inoculated wild and OsRLK-overexpressing lines compared with wild-type and OsRLK-knockout lines, but did not differ among lines under a heat-killed SynCom treatment (Fig. 6c).

To further assess OsRLK effects on Acinetobacter SynCom function, we quantified Acinetobacter and nifH abundances by quantitative real-time PCR (qPCR), and quantified nitrogenase protein abundance by an enzyme-linked immunosorbent assay (ELISA) kit. Without SynCom, traits and growth did not differ among wild-type (1504), OsRLK-overexpressing lines, OsRLK-knockout lines and wild rice (Fig. 6d–f, P < 0.05). After 30 days with SynCom, Acinetobacter, nifH and nitrogenase protein abundance ranked OsRLK-1 and wild rice > OsRLK-2 > wild-type > OsRLK-knockout mutant lines (Fig. 6d–f, P < 0.05). Plant growth was comparable across all lines without SynCom inoculation (Fig. 6g,h) but was enhanced in OsRLK-overexpressing lines and wild rice grown with SynCom inoculation (Fig. 6i).

To compare agronomic traits between the wild-type and OsRLK overexpression lines, we conducted field trials and evaluated yield-related traits. OsRLK overexpression lines significantly increased effective panicle number, seed-setting rate, thousand-grain weight and grain yield per plant (Supplementary Table 7, P < 0.05). In contrast, no significant differences were observed in panicle length, grain length or grain width between the overexpression lines and the wild type.

To analyse natural variation of OsRLK and its role in Acinetobacter enrichment, we aligned OsRLK sequences from 382 rice accessions (including wild rice and cultivated rice), revealing significant differences between wild and cultivated rice accessions (Supplementary Table 4). Coding-region sequencing of 105 accessions (33 wild, 72 cultivated) identified 28 single-nucleotide polymorphisms (SNPs) and 2 indels, defining 3 haplotypes: HAP1 (31 wild), HAP2 (2 wild, 68 cultivated) and HAP3 (4 cultivated) (Supplementary Fig. 2). Amino acid alignment showed 18 variations between HAP1 and HAP2, but none between HAP2 and HAP3, indicating 2 main functional variants (Supplementary Fig. 3). Phylogenetic analysis based on OsRLK coding sequences from 105 rice accessions showed clear separation of HAP1 and HAP2, while HAP2 and HAP3 were closely related (Extended Data Fig. 7a). In addition, Acinetobacter (ASV1033682, ASV163285, ASV1936293) abundance was significantly higher in the HAP1 functional variant than in HAP2, with OsRLK in wild rice belonging to HAP1 and OsRLK in cultivated rice to HAP2 (Extended Data Fig. 7b,d,e and Supplementary Table 8).

This study shows that core rhizosphere bacteria from parental wild rice can be reinstated in RIL progeny under field conditions, supporting a genetic basis for microbiome assembly. This extends previous associations between plant loci and rhizosphere community composition in Zea mays¹⁸, Setaria italica¹⁹, tomato¹², barley¹¹ and Triticum aestivum¹³ by identifying the M-gene receptor-like kinase-encoding OsRLK and demonstrating its combinatorial beneficial effects with an Acinetobacter SynCom on NUE. Consistent with our previous findings that rice domestication reduced nitrogen-fixing functions and increased N₂O-related genes in rice rhizosphere¹⁵, and that indica-associated microbiota have greater nitrogen-metabolism potential than japonica²⁰, our findings indicate that wild-rice-based hybridization can partially restore beneficial microbial functions lost during domestication.

Beyond environmental filtering, our findings show that host genetic variation can influence the reassembly of functionally relevant microbial communities, suggesting that microbiome traits may be partially inherited and restored through breeding.

Although RILs are not cultivated accessions, they retain maternal cultivated genetic backgrounds while enabling recruitment of parental wild-rice core bacteria. They, therefore provide a genetic framework for restoring wild-derived microbial functions, including improved nitrogen acquisition and reduced greenhouse gas emissions, into modern rice breeding and potentially recovering host–microbe interactions disrupted during domestication.

Specific paternal wild loci on chromosomes 8, 11 and 12 were associated with multiple core bacterial taxa, highlighting the role of host genetics in microbiome assembly. However, loci for several core taxa could not be precisely mapped, potentially owing to QTL limitations²¹ or the absence of environmental stresses that favour their recruitment. For example, Gaiella, Paenarthrobacter and Nocardioides are enriched under salinity, heavy metal and cold stress^{22,23,24,25,26}. Thus, environmental cues required for gene positioning in plants^27,28 may also influence genetic control of rhizosphere microbiome assembly. Notably, QTLs for nitrogen metabolism-related core bacteria overlapped with reported rice NUE-associated genes, including OsAMT1;3(ref. ²⁹), which mediates ammonium uptake and associated nitrate–ammonium dynamics that alters root metabolites²¹, and shape rhizomicrobiome recruitment.

Comparing domestication loci with NUE-associated core-taxon QTLs suggests that microbiome assembly has been shaped by both domestication and natural selection. Among candidate genes, the receptor-like kinase-encoding gene OsRLK^30,31 emerged as a strong candidate based on its association with Acinetobacter, expression pattern and divergence between wild and cultivated rice. Recent transcriptomic analyses implicate OsRLK in immune regulatory networks³², while DUF26/CRK receptor-like kinases have been linked to symbiotic nitrogen fixation in legumes³³. Although its molecular mechanisms remain unsolved, we hypothesize that OsRLK regulates root exudation, thereby modifying the rhizosphere chemical environment and recruitment of beneficial microbiota. Future studies integrating root exudate metabolomics with host transcriptomics and microbial functional analyses will help identify the molecular signals and metabolic pathways underlying OsRLK-mediated microbiome recruitment.

The shattering gene Sh4, a major domestication locus, is located close to OsRLK (34.2 Mb), raising the possibility that favourable OsRLK alleles were indirectly counterselected through linkage with Sh4 (ref. ³⁴). Their genetic distance (~2 centimorgans (cM)) suggests that this linkage could be resolved through recombination, enabling introgression of OsRLK through fine mapping and marker-assisted selection while minimizing linkage drag. Genome editing could provide an alternative route to deploy beneficial alleles³⁵. However, practical application of these strategies requires further validation and optimization. This association may reflect retained interactions between plant loci and beneficial microbes rather than direct selection on nitrogen traits and may involve trade-offs, including immune regulation. Nevertheless, our findings provide an experimental basis for M-gene-based strategies to improve NUE by restoring beneficial microbiome functions, extending beyond rewilding towards microbiome-driven crop improvement.

Our study provides field-based evidence that beneficial microbiome functions lost during rice domestication can be recovered through inheritance of wild-rice genetic loci. Moving beyond correlation-based studies, we identify host loci that enable the re-domestication of beneficial microbiota and their NUE-related functions. We therefore propose ‘microbiome re-domestication’ as a genetic process whereby wild-rice alleles restore beneficial plant–microbiome interactions lost during domestication. Unlike microbiome rewilding, which primarily reintroduces wild microbial diversity, re-domestication uses host genetic variation to restore these interactions, providing a framework for microbiome-guided breeding of nutrient-efficient crops.

The experimental framework consists of three parts (Fig. 1): First, phenomenon identification, where we mapped genetic loci linked to wild-rice core rhizosphere microbiota and their overlap with domestication loci in RILs (Fig. 1a,b). Second, mechanism elucidation, using transcriptome analyses to identify plant genes regulating microbial interactions (Fig. 1c,d). Third, functional validation, where we engineered transgenic rice and tested a synthetic Acinetobacter community to assess their role in nitrogen-use efficiency (Fig. 1e,f). This integrative approach combines genetic mapping, transcriptomics and functional validation to reveal gene–microbiota interactions enhancing nitrogen utilization.

Identification of core bacterial groups in the rhizosphere of wild rice and linkage analysis

Plant material and growth conditions

Chinese Dongxiang common wild rice (Oryza rufipogon), representing the northernmost naturally distributed wild-rice germplasm, was used as paternal parent. This genotype, O. rufipogon Dongxiang, possesses cold resistance loci and nutrient deficiency resistance loci. Two F4 RIL populations, each consisting of 100 lines, were developed by crossing two O. sativa cultivars (female parents O. sativa L. ssp. japonica [Longgeng 31] and O. sativa L. ssp. indica [Huanghuazhan]) with male parent O. rufipogon Dongxiang.

Experimental plots for the parental lines and the two RIL populations were established at the demonstration base of the Northeast Institute of Geography and Agroecology, Chinese Academy of Sciences, located at the National Agricultural Station in Liaoheyuan (42° 54’ 15” N, 125° 26’ 10” E), Jilin Province, China. This site has temperate monsoonal climate and black soil, as classified by the United States Department of Agriculture (USDA) soil taxonomy.

Rice seeds from each RIL population and the parental lines were thoroughly mixed, and seeds were randomly selected. These seeds were soaked in tap water for 24 h, then disinfected with 70% ethanol followed by 2% sodium hypochlorite for 20 min. After disinfection, the seeds were germinated at 30 °C for 48 h. The germinated seeds were sown in 448-cell seedling trays (26 mm height, 16.5 mm diameter), with one seed per cell. One tray was used for each parental line and three trays for each RIL population.

One-month-old seedlings were transplanted to the field on 10 May 2022, with rows spaced 3–40 cm apart and one seedling per hill. The field was flooded to a depth of 10 cm until 10 days before full heading. Sampling occurred from 26–30 September 2022. To minimize border effects and variety interactions, plants were collected from the central areas of each plot. For microbiome and genome sequencing, 6 plants from the parental lines and 100 plants from each RIL population were randomly selected. Rhizosphere soil was collected for bacteria amplicon sequencing following established protocols¹⁵.

Rhizosphere microbial community

We assessed the microbial community composition in rice rhizosphere by 16S rRNA gene amplicon sequence. DNA was extracted from 0.5 g of rhizosphere soil using the Fast DNA SPIN kit (MP Biomedicals). The V3–V4 hypervariable regions of bacterial 16S rRNA gene were amplified using primers 338F (5’-ACTCCTACGGGAGGCAGCA-3’) and 806R (5’-GGACTACHVGGGTWTCTAAT-3’). The thermal cycling conditions included: initial denaturation at 98 °C for 8 min; 25 cycles of denaturation at 98 °C for 15 s, annealing at 55 °C for 30 s and extension at 72 °C for 30 s, followed by a final extension at 72 °C for 10 min. The PCR amplicons were purified using Agencourt AMPure Beads (Beckman Coulter), quantified and subjected to 250-bp paired-end sequencing on the Illumina NovaSeq platform. A total of 34,356,465 16S rRNA gene amplicon read pairs were generated from 424 samples. Raw reads were filtered using Trimmomatic v.0.33, and primer sequences were removed with cutadapt 1.9.1, resulting in 34,177,527 clean read pairs³⁶. These clean reads were then denoised using the Dada2 plugin in QIIME2, producing amplicon sequence variants (ASVs)³⁷. ASV taxonomic classification was performed on the basis of the Naive Bayes classifier in QIIME2 using the SILVA database (v.138.1)³⁸. ASVs identified as non-bacterial, chloroplast or mitochondrial were excluded. The FAPROTAX database was utilized to predict the potential functions of the core rhizobacteriome¹⁶. The gjam package in R (v.4.1.2) facilitated the identification of the core rhizobacteriome, consisting of microbes with regression coefficients that were positive and significantly different from zero across the populations. The model was run with a total of 9,000 iterations, of which the first 3,000 were discarded as burn-in to ensure convergence.

Plant physiological properties and nitrogen content

Leaves were used for chlorophyll concentration analysis, and the remaining plant material was used to measure fresh and dry weights of shoots and roots, as well as nitrogen content. Replicate counts varied by analysis: n = 6 for parental lines and bulk soil, n = 3 per RIL population sample. Rice parameters including dry weight, plant height, root length, main spike length, chlorophyll content and nitrogen content were determined. Plant height, main spike length and main root length were measured using standard rulers, and ground diameter was assessed using Vernier calipers. Fresh shoots and roots were weighed and then dried at 60 °C for 48 h to determine their dry weights. The dried shoots and roots were ground together and analysed as a composite sample to measure plant nitrogen content using a continuous-flow analyser (AutoAnalyzer3 Continuous-Flow Analyzer, Bran Luebbe). Chlorophyll was extracted from 0.1 g of fresh leaves using 95% ethanol, with three replicates per accession. Chlorophyll concentration (Chl) in the leaves was quantified using a spectrophotometer, adhering to a previously described method³⁹. The equations defining the relationship between chlorophyll concentration and optical density are as follows:

where D₆₆₅ and D₆₄₉ represent the optical densities at wavelengths of 665 nm and 649 nm, respectively; Ca, Cb and CT denote the concentrations of chlorophyll a, chlorophyll b and total chlorophyll, respectively. The coefficients 83.31 and 18.60 are the specific absorptions for chlorophyll a and chlorophyll b at 665 nm, respectively; 24.54 and 44.24 are the specific absorptions at 649 nm.

Plant DNA extraction and SNP calling for two RIL populations

Plant genomic DNA was extracted from fresh leaves of the two rice populations at rice maturity stage. Genomic DNA was extracted from young leaves using a cetyltrimethylam-monium bromide (CTAB)-based protocol⁴⁰. The concentration and quality of the total genomic DNA were assessed using a NanoDrop2000 spectrophotometer (Thermo Fisher). DNA libraries (350 bp) for Illumina/BGI sequencing were constructed for each accession according to manufacturer specifications. After DNA library construction, sequencing was performed on an Illumina HiSeq XTen/NovaSeq/BGI platform by a commercial service (Biomarker Technologies), with 150-bp read lengths. The raw reads were further processed with the bioinformatic pipeline tool BMKCloud (www.biocloud.net) online platform. Raw data (raw reads) of fastq were first processed through fastp software. In this step, clean data (clean reads) were obtained by removing reads containing adapter, reads containing ploy-N and low-quality reads (quality score <20) from raw data. At the same time, Q20, Q30, GC content and sequence duplication level of the clean data were calculated. All downstream analyses were based on clean data with high quality. Raw sequences were transformed into clean reads after data processing. Clean reads were then mapped to the Oryza sativa Nipponbare MSU 7.0 genome using Burrows–Wheeler Aligner software⁴¹, as this genome is widely recognized for its high-quality assembly and annotation, providing a reliable reference essential for accurate mapping and meaningful genetic insights. Duplicate reads were removed using Samtools⁴², followed by local realignment and base recalibration with Genome Analysis Toolkit (GATK)⁴³. Genotypes for 76 O. sativa L. ssp. japonica and 90 O. sativa L. ssp. indica accessions were determined using a combined GATK and Samtools SNP calling analysis under default settings.

SNPs exhibiting polymorphism between the two parentals were classified into eight segregation patterns, with only SNPs following the aa × bb pattern selected for further analysis. Recombination breakpoints were identified using a sliding-window approach to construct a physical bin map, with bins serving as markers for developing an ultra-high-density genetic map. These markers were grouped into linkage groups based on their physical locations on the Oryza sativa Nipponbare whole genome. To ensure linkage group integrity, modified logarithm of odds (MLOD) scores between markers were calculated. Markers with MLOD scores below 5 were excluded from ordering. The HighMap strategy⁴⁴ was used to order SNP markers and correct genotyping errors within linkage groups. Recombinant frequencies and LOD scores were determined through two-point analysis to infer linkage phases. Marker ordering was performed iteratively using simulated annealing algorithms^45,46. The SMOOTH strategy⁴⁷ was applied to correct errors based on parental genotype contributions, and a k-nearest neighbour algorithm was used to impute missing genotypes⁴⁸. Skewed markers were integrated into the map using a multipoint maximum likelihood approach, and map distances were estimated using the Kosambi mapping function⁴⁹.

QTL analysis

After analysing the core rhizosphere bacterial community of paternal wild-rice accessions and excluding the core microbiomes of bulk soil and cultivated parental lines, ASVs exclusively selected by the paternal wild-rice microbiome were mapped within the RIL populations. The QTL analysis linking selected core ASV features with plant loci was conducted using the R package R/qtl2 (ref. ⁵⁰), with genotypic and phenotypic data obtained from the same individual plants. Interval mapping using the expectation-maximization (EM) algorithm and a single-QTL model was implemented. A randomized permutation test with 1,000 permutations was performed to assess the distribution of LOD scores. The threshold for identifying significant LOD peaks was set at the 95th percentile using a Bayesian credible interval with a probability of P = 0.95.

Analysis of overlap between core bacterial taxa-related QTLs and domestication loci

Screening of rice domestication loci

To explore rice domestication loci, 117 wild rice (O. rufipogon) and 265 cultivated rice (O. sativa L. ssp. japonica and O. sativa L. ssp. indica) accessions were planted in the same field on 15 May 2023. The field was maintained with a water depth of 10 cm until 10 days before the full heading. Whole plants, including roots, were harvested between 20 and 30 September 2023. To minimize border effects and interference between accessions, plants of the same variety were collected from the centre of each planting area. Genomic sequences for 260 rice accessions were obtained from a published study⁵¹ (https://www.ncbi.nlm.nih.gov/sra/?term=PRJEB6180), while the remaining 122 accessions were sequenced in this study (Supplementary Table 4).

Plant genomic DNA was extracted from fresh leaves, and genomic sequences were processed as described in Experiment II. SNPs were identified within the 382 accessions using the HaplotypeCaller module in GATK (v.3.8)⁴³. These SNPs were filtered using the following parameters: QD < 2.0||MQ < 40.0||FS > 60.0||QUAL < 30.0||MQrankSum < −12.5||ReadPosRankSum < −8.0, -clusterSize 2, -clusterWindowSize 5. Only SNPs with a minor allele frequency greater than 5% and less than 80% missing data were retained. A total of 1,489,886 high-quality SNPs were obtained. Nucleotide diversity (π) was calculated for each SNP using vcftools⁵², employing a 100-kb sliding window with a 10-kb step size for each subpopulation. The genetic diversity ratio (π_wild/π_cultivar) between wild and cultivated rice was quantified across the genome, and a threshold value was determined using permutation tests. This ratio was used to identify domestication loci in the indica and japonica varieties population. Comparative analysis of the 20 core ASVs common to wild, RIL-J and RIL-I population-associated QTLs and plant domestication loci was then performed. Overlapping confidence intervals between core ASV-based QTLs and domesticated loci, where the genetic diversity ratio exceeded the threshold, indicate domestication selection to those loci.

Transcriptome analysis of the impact of wild-rice rhizosphere microbiome on plants

Pot experiment

To assess whether the rhizosphere microbiome of paternal wild-rice accessions enhances plant growth and development, Chinese Dongxiang common wild-rice accessions were planted in sterilized soil (W) and non-sterilized soil (WN). Rice seeds were surface sterilized with 2.5% NaClO and germinated in deionized water for 5 days at 26 °C in darkness. Seedlings were then transferred to sterilized and non-sterilized black soil, collected from the Northeast Institute of Geography and Agroecology experiment station and sterilized with gamma irradiation at 25 kGy. Rice plants were grown under controlled environmental conditions: 15 h of light at 35 °C and 9 h of darkness at 30 °C, with 70–80% relative humidity.

Transcriptome analyses of wild rice

Root samples were collected at 30 days post transplantation following previously described methods. Samples were washed with demi water and immediately frozen in liquid nitrogen for RNA extraction using a Promega RNA extraction kit (Promega, LS1040). The quality and concentration of the extracted RNA were assessed using an RNA Nano 6000 Assay kit on a Bioanalyzer 2100 system (Agilent). Total RNA was then used for complementary (c)DNA library preparation, involving mRNA isolation, fragmentation and reverse transcription. Double-stranded cDNA was synthesized, size-selected and amplified via PCR. The quality of the library was verified using the Agilent Bioanalyzer 2100 system. Index-coded samples were clustered and sequenced on an Illumina HiSeq XTen platform (Illumina) on a cBot Cluster Generation system following manufacturer instructions. Sequencing was performed on an Illumina Novaseq platform, producing 150-bp paired-end reads.

Raw fastq data were processed to remove adapters, poly-N sequences and low-quality reads using Perl scripts. Clean reads were aligned to the reference genome using Hisat2 v.2.0.5. The assembly of mapped reads for each sample was performed with StringTie (v.1.3.3b) using a network flow algorithm, which quantified full-length transcripts and splice variants for each gene locus⁵³. Differential expression analysis between the two soil treatments was conducted using edgeR. Genes exhibiting differential expression between sterilized and non-sterilized soil-planted wild rice were identified using corrected P values (q-values) < 0.05, controlled using the Benjamini–Hochberg method, with a significance threshold set at a q-value of 0.05 and an absolute fold change of 2.0. Gene Ontology (GO) enrichment analysis on differentially expressed genes was performed using the clusterProfiler R package, considering GO terms with q-values < 0.05 as significantly enriched. KEGG pathway annotations were facilitated using KAAS (KEGG Automatic Annotation Server: http://www.genome.jp/kegg/kaas/), and statistical enrichment of differentially expressed genes in KEGG pathways was tested using the same package.

Isolation of core bacterial taxa and construction of synthetic community

Following our findings that strong selective sweeps in domestication loci are associated with nitrogen fixation-related core microbes Acinetobacter ASV6549, this bacterium was notably linked to rice candidate gene segments related to nitrogen utilization efficiency. Given its high additive effect difference (ADD) value and phenotypic variance explained (PVE) value, Acinetobacter was selected for further validation as a core nitrogen-fixing bacterium in the rhizosphere of paternal wild-rice accessions, potentially influencing plant loci associated with nitrogen use efficiency.

Acinetobacter was isolated from rhizosphere soil samples collected from paternal wild rice. Approximately 0.5 g of rhizosphere soil samples were homogenized in distilled water, and the resultant solution was diluted and cultured on medium at 28 °C for 3 days. The medium includes potassium nitrate (1.0 g l⁻¹), asparagine (1.0 g ⁻¹), bromothymol blue (BTB) (0.05 g l⁻¹), sodium citrate (8.5 g l⁻¹), magnesium sulfate heptahydrate (1.0 g l⁻¹), ferric chloride hexahydrate (0.05 g l⁻¹), potassium dihydrogen phosphate (1.0 g l⁻¹) and calcium chloride dihydrate (0.2 g l⁻¹). After four rounds of re-plating on culture medium plates, individual colonies were selected and identified through 16S rRNA amplicon sequencing, following the methods described above. qPCR of the nifH gene was conducted to verify the nitrogen fixation capability of the isolated strains, employing primers nifHF (5’-AAAGGYGGWATCGGYAARTCCACCAC-3’) and nifHR (5’-TTGTTSGCSGCRTACATSGCCATCAT-3’). Plasmid DNA extracted from microbiomes harbouring nifH gene from rhizosphere samples of paternal wild-rice accessions was used to generate standard curves, which were then cloned into vectors. The standard curve involved 10-fold serial dilutions ranging from 10¹ to 10⁷. High reaction efficiencies (90–105%) were achieved for gene quantification, with R² values exceeding 0.99. qPCR analyses were performed in biological triplicates, each with 3 technical replicates. Concurrently, a minus-template qPCR using distilled water instead of rhizosphere DNA was conducted as a control. On the basis of their nitrogen fixation capabilities, 9 strains of Acinetobacter sp., including two A. baumannii strains, A. seifertii strain, A. halotolerans strain, A. sichuanensis strain, A. schindleri strain, A. geminorum strain, A. defluvii strain and A. lactucae strain were isolated from the rhizosphere soil of wild rice (Extended Data Fig. 5). PCR amplification was used to confirm the isolated strains. These strains were cultured overnight in tryptic soy broth medium at 37 °C. Cells from each strain were then washed and combined in equal concentrations to construct a synthetic community of Acinetobacter.

Transcriptomic analysis of the impact of Acinetobacter on plants

To identify candidate genes regulating the genus Acinetobacter, we planted Chinese Dongxiang common wild rice in two soil conditions: sterilized soil inoculated with Acinetobacter synthetic bacterial communities (WA) and sterilized soil (W) only. Acinetobacter synthetic bacterial communities (1.0 ml, with an optical density at 600 nm (OD₆₀₀) of 0.5) were added to 100 ml of sterilized water and combined with 100 g of black soil (resulting in a final concentration of 10⁶ colony-forming units (c.f.u.s) g⁻¹) for planting wild rice. The experimental design, including rice materials, seed sterilization and planting management, followed the protocol described above for the pot experiment assessing transcriptome analyses of wild rice. The transcriptomic analyses adhered to the methods described above.

Construction of genetically modified rice material

The selection of candidate genes was guided by three criteria: First, candidate genes must exhibit simultaneous upregulation in comparisons of both WN with W, and WA with W; second, they should demonstrate a log₂(fold change) > 2; third, the transcript abundance of candidate genes should exceed the average transcript abundance of all genes under the same conditions. On the basis of these criteria, Os04g0659300 (OsRLK) was selected as a candidate gene on chromosome 4 for regulating ASV6549 (Acinetobacter). To test the hypothesis that OsRLK contributes to nitrogen utilization efficiency and influences Acinetobacter, OsRLK overexpression lines and knockout mutants were generated using indica rice 1504 (wild-type) and Chinese Dongxiang wild rice (wild rice) as recipient materials, respectively. The pU1301 vector⁵⁴ containing a maize ubiquitin promoter was used for the overexpression studies. To construct the overexpression vector of OsRLK, the genomic sequence of OsRLK from Chinese Dongxiang wild rice was amplified and cloned into pU1301. The OsRLK mutants were generated using CRISPR/Cas9 technology as previously described⁵⁵. Two specific protospacer adjacent motif (PAM) sequences (PAM1: CACTTTGCTCGTTCTCCTCGTGG, PAM2: CGAGAACAACCTCCTGAACCTGG) of OsRLK were selected. Subsequently, the corresponding single-guide RNAs (sgRNAs) (sgRNA1: CACTTTGCTCGTTCTCCTCG, sgRNA2: CGAGAACAACCTCCTGAACC) were introduced into the pBWA (V) H2S-cas9pl (osu3) i vector backbone to construct the Cas9-OsRLK recombinant vector. The constructs were confirmed by sequencing, introduced into Agrobacterium tumefaciens strain EHA105 by electroporation, and transformed into the corresponding recipient materials by Agrobacterium-mediated transformation. The primers used for vector construction are listed in Supplementary Table 9. The PAM sequences of OsRLK mutant plants were sequenced to detect their mutation types, and homozygous mutants were selected for subsequent experiments (Extended Data Fig. 8).

For Southern blot hybridization to determine the copy number of the OsRLK^DY transgene, 10 µg of genomic DNA from both overexpressed and wild-type rice materials from 30-day-old seedlings was digested with 15 U of HindIII (NEB) for 16 h at 37 °C. The digested DNA was resolved on a 0.8% agarose gel at 30 V for 16 h, followed by capillary transfer onto a Hybond N+ nylon membrane (GE Healthcare) overnight. The blot was hybridized at 48 °C with a DNA probe complementary to the HPT gene, which encodes hygromycin B phosphotransferase (HPT) and serves as a selectable marker in genetic transformation, and then exposed to X-ray film at −80 °C for 4 h for imaging.

All transgenic plant materials were planted in Jiangxi Academy of Agricultural Sciences Biological Breeding Base at Yichun (28° 14’ 57” N, 115° 07’ 05” E) in 2024 for phenotype investigation. Each line was grown in a two-row plot with 10 plants per row, with a spacing of 16.7 cm between plants and 21.6 cm between rows. At the mature stage, effective panicle number per plant and panicle length were evaluated following the methods described previously⁵⁶. Spikelet number per panicle, seed-setting ratio, grain length, grain width, thousand-grain weight and grain weight per plant were measured using a yield traits scorer in a high-throughput phenotyping facility⁵⁷.

Validation of genetically modified rice materials in regulating Acinetobacter

Experiment design

The experiment utilized overexpressed rice genotype and synthetic microbial communities to confirm the regulatory relationship between gene OsRLK and the Acinetobacter synthetic community. This synthetic microbial community and a heat-killed synthetic microbial community (used as control) were inoculated into sterilized soil alongside both overexpressed and wild-type-rice material. Each pot received 1.0 ml of synthetic communities (with an OD₆₀₀ of 0.5) mixed with 100 g of black soil, resulting in a final concentration of 10⁶ c.f.u.s g⁻¹. Rice plants were grown under controlled conditions: 35 °C with 15 h of light per day and 30 °C during 9 h of darkness, maintaining 70–80% relative humidity. Rhizosphere soil samples were collected at 30 days post inoculation (dpi) to measure Acinetobacter abundance, nifH gene abundance and nitrogenase activity. Rhizosphere soil (0.5 g) was suspended in PBS at a 1:10 (w/v) ratio, vigorously vortexed for 40 s and then centrifuged at speed 4,000 revolutions per min for 15 min. The resulting supernatant was subjected to nitrogenase protein quantification using the Soil Nitrogenase (NITS) ELISA kit (ADANTI). At 30 dpi, rice plants were harvested for further analysis, with each treatment including 5 biological replicates.

Samples determination

DNA was extracted from 0.5 g of rhizosphere soil samples collected at 0, 5, 15 and 30 dpi using the Fast DNA SPIN kit (MP Biomedicals). qPCR of the Acinetobacter gene was conducted to verify Acinetobacter colonization, employing primers Ac436f (5’TTT AAG CGA GGA GGA GG 3’) and Ac676r (5’ATT CTA CCA TCC TCTCCC3’)⁵⁸. nifH primers and qPCR analysis was performed following the procedures outlined in Experiment Ⅴ. Leaf samples were also collected at 30 dpi to assess the expression of the OsRLK gene, using the RT–qPCR methods previously described⁵⁹. In addition, dry shoots and roots were combined and ground together to form a single sample for measuring plant nitrogen content. This analysis was conducted using a continuous-flow analyser (AutoAnalyzer3 Continuous-Flow Analyzer, Bran Luebbe).

Haplotype analysis and phylogenetic tree construction

The natural variations in the coding region of OsRLK in the world-wide rice core collection of 105 accessions (including 33 accessions of wild rice and 72 accessions of cultivated rice) were obtained from sequencing. Haplotype analysis and phylogenetic tree construction were performed on the basis of the sequencing data of 105 germplasm accessions. The EMBL software CLUSTAL 2.0.12 multiple sequence alignment tool was used to perform the coding of DNA sequences and protein alignment. Phylogenetic trees were constructed using the molecular evolutionary genetics analysis version 7.0 (MEGA7) software⁶⁰.

Statistical analysis

Differences in the rhizosphere bacterial community structure between the parental lines and the two RIL populations were assessed using non-metric multidimensional scaling (NMDS) based on the Bray–Curtis dissimilarity matrix in the vegan package in R (v.4.1.2). PERMANOVA was employed to evaluate the statistical significance of clustering patterns observed in the ordination plots. GJAM was conducted to determine the core bacteriome profiles and to identify core bacteriomes of different accessions recruited under similar soil conditions^61,62. The gjam package in R (v.4.1.2) facilitated the identification of the core rhizobacteriome, consisting of microbes with regression coefficients that were positive and significantly different from zero across the populations. Visualizations, including Sankey diagrams, were created using Origin 2023. On the basis of the taxonomic classification, a species tree was inferred using GraPhlAn, which displays taxonomic data generated by other metagenomic analysis tools without branch lengths^62,63. Correlations among plant physiological and biochemical properties were determined using Pearson’s correlation coefficients via the linKEt package in R (v.4.1.2). Mantel correlations between core microbial composition and plant physiological and biochemical data, with 9,999 permutations, were analysed using the vegan package in R (v.4.1.2)⁶⁴.

Reporting summary

Further information on research design is available in the Nature Portfolio Reporting Summary linked to this article.