Front. Plant Sci. Frontiers in Plant Science Front. Plant Sci. 1664-462X Frontiers Media S.A. 10.3389/fpls.2025.1616424 Plant Science Original Research Genetic dissection of root traits in a rice ‘global MAGIC’ population for candidate traits to breed for reduced methane emission Roy Ripon Kumar 1 2 3 Misra Gopal 1 Sharma Shaina 1 4 Pahi Bandana 1 4 Khatibi Seyed Mahdi Hosseiniyan 1 Trijatmiko Kurniawan Rudi 1 Kim Sung-Ryul 1 Hernandez Jose E. 3 Henry Amelia 1 Sreenivasulu Nese 1 Diaz Maria Genaleen Q. 3 Ocampo Eureka Teresa M. 3 Sinha Pallavi 1 4 Kohli Ajay 1 * 1 International Rice Research Institute (IRRI), Los Baños, Philippines 2 Bangladesh Rice Research Institute (BRRI), Gazipur, Bangladesh 3 University of the Philippines Los Baños, College, Laguna, Philippines 4 International Rice Research Institute, South Asia Hub, Hyderabad, Telangana, India

Edited by: Dayun Tao, Yunnan Academy of Agricultural Sciences, China

Reviewed by: Harun Bektaş, Siirt University, Türkiye

Subroto Das Jyoti, Texas A&M University College Station, United States

*Correspondence: Ajay Kohli, ajakoy@yahoo.com

09 07 2025 08 09 2025 2025 16 1616424 22 04 2025 04 06 2025 Copyright © 2025 Roy, Misra, Sharma, Pahi, Khatibi, Trijatmiko, Kim, Hernandez, Henry, Sreenivasulu, Diaz, Ocampo, Sinha and Kohli. 2025 Roy, Misra, Sharma, Pahi, Khatibi, Trijatmiko, Kim, Hernandez, Henry, Sreenivasulu, Diaz, Ocampo, Sinha and Kohli

This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

Rice cultivation is critical for global food security. The largely practiced method of rice cultivation by transplantation under flooded fields contributes significantly to methane (CH4) emissions, posing challenges to climate-smart agriculture. This study uses a multi-parent advanced generation inter-cross (MAGIC) population of 250 rice genotypes to understand the genetic basis of root traits that may govern CH4 mitigation. Phenotyping under controlled greenhouse conditions revealed significant variation in root diameter (0.122–0.481 mm) and porosity (5.344–56.793%), and strong correlations between root diameter and porosity traits (r = 0.40–0.50, p < 0.001). Association studies revealed key candidate genes including Os05g0411200 (thermosensitive chloroplast development), Os10g0177300 (chalcone synthase), and Os04g0405300 (alcohol dehydrogenase), which regulate aerenchyma formation and auxin homeostasis. Protein-protein interaction networks linked these genes to flavonoid biosynthesis (KEGG map00941) and N-glycan pathways, earlier identified as critical for root architecture. Haplotype-phenotype analysis revealed 8 superior haplotypes across 7 genes for average root porosity, base root porosity, root diameter, and tip root porosity. These findings provide the foundation for breeding high-yielding rice varieties with reduced methane emissions, addressing the challenges of food security and climate change.

 root diameter root porosity genome-wide association analysis superior haplotype protein-protein interaction methane emission section-in-acceptance Plant Breeding Introduction

Rice is the staple food for more than 50% of the world’s population (Abeysekara and Rathnayake, 2024). Growing rice comes with the present understanding that it accounts for 12% of the global anthropogenic methane (CH4) emissions. Methane is a greenhouse gas (GHG) that is nearly 28 times more potent than carbon dioxide as a GHG. With the rising global population, especially in regions where rice is a staple food, meeting the growing demand for rice while reducing methane emissions from paddy fields will be an increasingly challenging task (Hertel, 2011; Sapkota et al., 2019).

Methane emission in flooded rice soils is a result of bacterial processes. These include the production and consumption of CH4 by methanogenic and methanotrophic bacteria, respectively, in anaerobic/aerobic microenvironments, and the oxidation of CH4 in aerobic microenvironments of the roots and the rhizosphere (Davamani et al., 2020; Luo et al., 2022). The rice plant plays a regulatory role in all these processes. Methane is produced by methanogenesis under reducing conditions (Kusin et al., 2025). Carbon sources in flooded rice fields provide substrates for soil methanogens, promoting a decline in soil redox potential (Eh) and creating suitable conditions for methanogen proliferation (Rago et al., 2015). Just after flooding, the soil transitions from an oxidized to a reduced state. Within a few hours to a day, electron acceptors such as O2, NO3 -, SO4 -2, Mn4+, and Fe3+ are depleted by various soil processes (Sahrawat, 2005).

Oxygen (O2) is utilized by aerobes and facultative anaerobes for aerobic respiration, nitrate (NO3 -) by denitrifiers for denitrification, iron-reducing bacteria reduce Fe3+ to Fe2+ through iron reduction, sulfate (SO4 -2) is reduced to sulfides (H2S, S2 –, and HS−) by sulfate-reducing bacteria, manganese-reducing bacteria reduce Mn4+ to Mn2+, and methanogenic bacteria produce CH4 using CO2 (Ponnamperuma, 1972; Sahrawat, 2005). The correlation between pH and Eh parameters at each observation stage has a negative correlation value, indicating an inverse relationship, which enhances the methanogenic bacterial population as well as CH4 production (Husson, 2013). On the other hand, the reduced forms can be re-oxidized by root-released O2 in the rhizosphere (Yang et al., 2012). Oxygen, being the most potent oxidizing agent, ensures the availability of electron acceptors in the soil (Husson, 2013). This availability limits the food intake of methanogens, which is toxic to them (Xu et al., 2003). Consequently, the methanogenic population and its activities are inhibited, thereby hampering CH4 production.

Rice plants play a role in both methane production and oxidation by diffusing O2 to the rhizosphere (Colmer, 2003). In rice fields, methanotrophic bacteria oxidize CH4 as part of their energy-gaining process (Fazli et al., 2013). However, they cannot utilize CH4 directly in the absence of O2. Methane first reacts with O2 to produce formaldehyde, which methanotrophic bacteria then use through the ribulose monophosphate (RuMP) pathway and the serine pathway to obtain energy. This process increases the methanotrophic population, facilitating more CH4 oxidation. Sharp counter-gradients of oxidized and reduced species characterize the soil surface in rice fields. Where these gradients overlap, methanotrophic bacteria oxidize ≧̸90% of potentially emitted methane before it is released into the atmosphere (Reim et al., 2012).

High-yielding rice varieties with higher root porosity have greater oxidation potential as they release more oxygen into the soil (Zheng et al., 2005; Jiang et al., 2013, 2017). Atmospheric oxygen is transported from the shoot through the well-developed rice aerenchyma to the roots and finally diffuses into the rhizosphere (Mei et al., 2009; Win et al., 2011). Aerenchyma, a cortical airspace (porosity), provides a low-resistance internal pathway for the movement of O2 from the shoot to the roots (Armstrong, 1971, 1980). The amount of radial oxygen loss (ROL) is determined by the oxygen concentration gradient, the physical resistance to radial oxygen diffusion between the aerenchyma and the soil, and the consumption of oxygen by cells along the radial diffusion path (Ejiri et al., 2021).

The ability for radial oxygen loss (ROL) is stronger when aerenchyma is well developed, as higher root porosity enhances the oxygen diffusion to the soil (Colmer, 2003; Mei et al., 2009; Jiang et al., 2017). Root porosity, defined as the ratio of the root space volume to the mass volume of the root, is a primary factor controlling root internal O2 concentration and ROL (Luxmoore et al., 1970). Porosity in plant tissues results from intercellular gas-filled spaces formed during development and can be further enhanced by the formation of aerenchyma (Armstrong, 1980; Raven, 1996). Rice varieties with higher porosity tend to have higher rates of ROL (Mei et al., 2009). Root porosity is significantly associated with root diameter, with increased porosity depending on increased root diameter (Striker et al., 2007). It has been stated that root porosity depends on root diameter because a larger aerenchyma area develops in roots with greater diameter, ensuring more O2 secretion to the rhizosphere (Kirk, 2003; Kim et al., 2018).

From the above discussion, high-yielding rice varieties with higher radial oxygen loss capacity might be a viable option to mitigate CH4 emissions from rice fields. Developing such varieties requires genetic information about root diameter and porosity. Only a few genes, including LESION SIMULATING DISEASE 1 (LSD1), ENHANCED DISEASE SUSCEPTIBILITY 1 (EDS1), and PHYTOALEXIN DEFICIENT 4 (PAD4), have been reported to regulate lysigenous aerenchyma formation in Arabidopsis in response to hypoxia through H2O2 and ethylene signaling (Mühlenbock et al., 2007). The QTLs Qaer1.02-3, Qaer1.07, Qaer5.09, and Qaer8.06–7 have been reported for maize root aerenchyma formation under non-flooding conditions (Mano et al., 2007). However, there is a lack of studies on the genetic control (QTLs, QTNs, or genes) of rice root diameter, root porosity, and root aerenchyma.

Genome-wide association studies (GWAS) are becoming a popular method for linking genotypic variation with corresponding trait differences in several crops (Ingvarsson and Street, 2011; Xiao et al., 2017). This method facilitates the identification of major allelic variants and haplotypes in candidate genes (Sinha et al., 2020). Therefore, in the current study, we will perform GWAS to reveal the quantitative trait nucleotides (QTNs) and candidate genes responsible for higher root diameter and root porosity, aiding in marker-assisted breeding to develop low methane-emitting rice varieties.

 Materials and methods Plant materials

A multi-parent advanced generation inter-cross known as “global MAGIC” population of 250 inbred lines, including 234 lines and 16 parents (eight indica and eight japonica), was developed by the International Rice Research Institute (IRRI) (Bandillo et al., 2013). The elite parents used in the study were originated in the different geographical regions of the world (Colombia, China, IRRI-Philippines, India, USA, Latin America, Korea, WARDA- Ivory Coast, Uruguay) having desirable agronomic traits (high yield, good grain quality), biotic (blast and bacterial blight diseases) and abiotic (drought, salinity, submergence) stress tolerance. Experiments were carried out in the greenhouse of the International Rice Research Institute (IRRI), Los Banos, Laguna, Philippines, for 44 days (23 December 2020 to 6 February 2021 - 1st Experiment) and 40 days (26 February to 7 April 2021 - 2nd Experiment). The temperature was 26.78°C and 28.79°C, and the relative humidity was 76.88% and 71.55% in the first and second experiments, respectively, measured by HOBO Pro v2 U23-001A - Temperature/Relative Humidity Data Logger (ONSET, 470 MacArthur Blvd., Bourne, MA 02532, USA). The experiments were laid out in a Complete Randomized Design (CRD) with three replicates. Three germinated seeds were sown in a pot to ensure a single plant per pot.

 Evaluation of the MAGIC population for the targeted traits

To investigate specific root traits that may impact methane emissions the study focused on traits such as average root porosity, base root porosity, tip root porosity, root diameter, tiller number, root number, root per tiller, shoot dry weight, root dry weight, root shoot ratio, root length, fine lateral root, thick lateral root, lateral root and nodal root on a diverse subset of 250 genotypes from the MAGIC panel ( Supplementary Table S1 ). Three roots from each plant were scanned to get an 8-bit grayscale image in an Epson Perfection 7000 scanner at 600 dots per inch resolution. The images were analyzed using the WinRHIZO root system analyzer software (WinRHIZO 2004, Regent Instruments, Inc., Quebec, Canada) to generate data on average root diameter (RDia, mm). The same roots together were used to measure the root porosity traits, such as base root porosity (BP), middle root porosity (MP), and tip root porosity (TP), to minimize the error. Three cm each were taken from the base, the middle, and from the tip of the root, and the porosity was measured by following the microbalance method (Visser and Bögemann, 2003). In the case of tip porosity measurement, 3 cm was taken from 0.5 cm away from the root tip to minimize the handling error. Three cm roots were chopped into 4 equal segments to fit into a 0.5 ml microfuge tube and easily evacuate air from the root in the vacuum concentrator. We have used a 0.5 ml microfuge tube (Eppendorf) instead of a gelatin capsule for measuring the sample. Additionally root length (RL), lateral root (LR) were counted and total root dry weight (RDW) were taken.

 Genomic data quality control

A total of 560 genotypes, including sixteen parents of the global magic population, were sequenced at an average 5x depth of coverage using the genotyping by sequencing (GBS) technique. Whole-genome sequence data of 16 parents were also present in the 3K rice project with a higher depth of coverage. Individual samples’ read data were demultiplexed into individual paired-end fastq files. Raw reads mapping and variant calling were done using a combination of BWA (Jung and Han, 2022) and GATK (McKenna et al., 2010) for every line of the population to the reference Japonica genome (version 7 of the Rice Genome Annotation Project, http://rice.plantbiology.msu.edu/).

Individual VCF files were generated for complete base calls across the entire genome for individual lines of the population, and the individual VCF files for the population were merged into a single genotype file. SNPs were filtered from this genotype file with a 30% missing rate. Single genotype files were also created for the sixteen parents, and SNPs were retained when a base pair position was present in all sixteen parents. This resulted in a set of 5 million high-quality SNPs for all sixteen parents. An SNP set with a 30% missing rate was input into the BEAGLE (Browning et al., 2018) pipeline to impute the low depth of coverage data by using parent SNP data as a reference set. The imputation process resulted in 4,40,350 SNPs across the population. Finally, 250 genotypes were selected based on genomic diversity, with 3,53,641 SNPs retained based on a 1% missing rate and a 5% minor allele frequency. The filtering method was carried out using PLINK (Purcell et al., 2007).

 Multi-locus association mapping

Multi-locus association mapping was performed on 250 genotypes with 3,53,641 SNPs. The analysis was conducted using five different Multilocus GWAS models, including mrMLM (Wang et al., 2016), FASTmrMLM, FASTmrEMMA, ISIS EM-BLASSO, and pLARmEB (Zhang et al., 2017). The final set of significant SNPs was obtained using a threshold criterion of a LOD of 3 and above.

 Candidate genes and functional annotation

To identify candidate genes underlying the MTAs, we considered all annotated genes located within a ±100kb length around each significant SNP, based on the estimated LD decay in the MAGIC population. Functional annotation of these genes was used as the first pass towards candidate genes, and the candidature was further refined on the basis of network analysis. The reference genome used was Oryza sativa cv. Nipponbare (IRGSP-1.0; RAP database: http://rapdb.dna.affrc.go.jp/download/irgsp1.html). The functional annotation of the candidate genes was ascertained using RAPDP (https://rapdb.dna.affrc.go.jp/), Rice Genome Annotation Project (http://rice.uga.edu/), KEGG (https://www.genome.jp/kegg/), NCBI (https://www.ncbi.nlm.nih.gov/gene), KnetMiner (https://knetminer.com/), and by reviewing the literature.

 Haplo-pheno analysis

The total candidate genes were taken into the haplo-pheno analysis for the identification of superior haplotypes. Haplotype analysis was carried out using in-house scripts in R programming. To identify the robust donors having superior haplotypes of the key genes, best linear unbiased prediction (BLUP) analysis was done. The BLUP values were utilized for the haplo-pheno analysis. Haplo-pheno analysis was performed to associate the identified haplotypes of the selected genes with superior phenotypes. Haplotypes present in fewer than three accessions and heterozygous SNPs were removed from the analysis. The genotypes were then categorized based on haplotype groups, and superior haplotypes were identified using the phenotypic data of the individuals in each haplotype group. Then the analysis included One-way ANOVA with haplotype as a fixed factor. Subsequently, Duncan’s multiple range test (DMRT) and ANOVA were used to test the statistical significance among the means of haplotype groups using the Agricolae package in R (de Mendiburu, 2019) Groups with different letters in the graphs indicate significant differences among the groups at a p < 0.05 level of significance.

 Protein-protein interaction network and tissue-specific expression profiling

To understand the protein interactions, String 12.0 (https://string-db.org/) (Szklarczyk et al., 2023), a database of known and predicted protein-protein interactions, was used. The input gene set was the candidate genes identified by our present study. The parameter of the PPI score was set as 0.4 (indicating medium confidence). Therefore, the PPI from String was collected for the construction of a differential protein interaction network among the candidate genes. The subset of potential genes found in the PPI network was further used for identifying Tissue-specific expression profiles using RiceXPro database, which provides high-resolution expression data across root and shoot tissues in rice (Sato et al., 2011).

 Results

Phenotypic values of fifteen root traits from 250 MAGIC lines were analyzed across two experiments to assess phenotypic variation. A broad range of variation was observed among the tested lines, as detailed in Supplementary Table S1 . In the first experiment, the traits average root porosity (ARP), base root porosity (BRP), and tip root porosity (TRP) ranged from 9.457 to 47.098%, 11.41 to 53.93%, and 5.344 to 43.911%, respectively ( Figure 1A ). Root diameter (Rdia) varied from 0.122 to 0.481 mm, while root number (RN) ranged from 39 to 250. Shoot dry weight (SDW) and root dry weight (RDW) varied from 1.13 to 6.49 g and 0.18 to 1.88 g, respectively. The root-to-shoot ratio (RSRatio) ranged from 0.061 to 0.605. Root length (RL) varied from 445.265 to 2054.341 cm, and lateral root number (LR) ranged from 64.206 to 87.4.

In other experiments, ARP, BRP, and TRP ranged from 11.19 to 47.77%, 10.738 to 53.859%, and 7.168 to 56.793%, respectively. Root diameter (RDia) varied from 0.122 to 0.449 mm, while root number (RN) ranged from 79 to 332. Shoot dry weight (SDW) and root dry weight (RDW) varied from 2.69 to 8.48 g and 0.37 to 3.11 g, respectively. The root-to-shoot ratio (RSRatio) ranged from 0.11 to 0.58. Root length (RL) varied from 343.421 to 1496.119 cm, and lateral root number (LR) ranged from 64.448 to 86.745.

The standard deviation ranges between the first and second experiments demonstrated minimal variation for traits such as TN (1.893-1.635), RPT (6.320-7.025), SDW (0.840-0.919), RDW (0.273-0.355), RS Ratio (0.071-0.055), RDia (0.051-0.052), BRP (6.656-6.788), TRP (6.236-6.141), ARP (5.460-5.284), FLR (6.092-5.415), TLR (4.651-4.275), LR (3.471-3.421), and NR (3.471-3.421). In contrast, traits such as RN (32.303-44.721) and RL (168.080-158.594) exhibited greater standard deviation ranges, indicating a higher degree of diversity within the Global MAGIC population ( Supplementary Table S1 ). The frequency distributions of phenotypes revealed continuous distributions, characteristic of quantitative traits, in both experiments. These results demonstrate significant phenotypic variation in root traits and methane emissions, providing a foundation for further genetic and environmental impact studies.

Correlation among the traits

Correlation trends among the studied traits were consistent across both experiments ( Supplementary Figures S1 , S2 ). RN exhibited strong correlations with RPT, RS ratio, TN, and RDW. Root porosity traits (BRP, TRP, and ARP) showed significant correlations with RDia. The correlation values of RL and FLR with NR were -0.38 and -0.61, respectively, indicating negative correlations. The correlation values between BRP and RDia, TRP and RDia, and ARP and RDia were approximately 0.45, 0.40, and 0.50, respectively. Additionally, BRP, TRP, and ARP were positively correlated, with ARP being the average of BRP and TRP. The correlation values between TRP and BRP, TRP and ARP, and BRP and ARP were around 0.45, 0.85, and 0.90, respectively. These results suggest that candidate genes associated with any of the studied traits could effectively enhance root porosity.

 Unveiling significant MTAs through multi-locus GWAS

GWAS analysis was done using 3,49,594 SNPs on 250 genotypes to unveil meaningful MTAs. Employing five multilocus GWAS methods mrMLM, FASTmrMLM, FASTmrEMMA, ISIS EM-BLASSO, and pLARmEB, significant associations were identified by considering a LOD value ≥3 along with the phenotypic variance (PVE) as a threshold for significance. In total, GWAS successfully identified 193 and 195 significant QTNs in the 1st and 2nd experiments, respectively, associated with four root traits (RDia, BRP, TRP, and ARP) across the 12 chromosomes in all four multi-locus GWAS methods.

For Average root porosity (ARP), a total of 48 and 41 significant MTAs were identified from the GWAS analysis in both experiments, respectively. These MTAs were distributed across all chromosomes in both experiments ( Supplementary Table S2 ). Among these significant MTAs, 30 and 25 were found within the genic regions in both experiments, respectively ( Supplementary Table S3 ). The phenotypic variance explained (PVE) by the MTAs ranged from 0.34% to 12.59% in 1st experiment and 0.2% to 12.55% in 2nd experiment. Furthermore, we examined the consistency of MTAs between the two experiments and identified 12 common MTAs. The soundness of the model is evident from the well-fitted Manhattan and Q-Q plots ( Supplementary Figure S3 ).

For Base root porosity, a total of 100 significant MTAs were identified from the GWAS analysis in both experiments. These MTAs were distributed across all chromosomes in both experiments except chromosome 4 in 2nd experiment ( Figure 1B , Supplementary Table S4 ). Among these significant MTAs, 29 were found within the genic regions in both experiments ( Supplementary Table S5 ). The phenotypic variance explained (PVE) by the MTAs ranged from 0.34% to 10.71% in 1st experiment and 0.23% to 9.67% in 2nd experiment. Furthermore, we examined the consistency of MTAs between the two experiments and identified 20 common MTAs. The soundness of the model is evident from the well-fitted Manhattan and Q-Q plots ( Supplementary Figure S4 ).

For root diameter, a total of 45 and 60 significant MTAs were identified from the analysis in the 1st and 2nd experiments, respectively ( Figure 1B ). These MTAs were distributed across all chromosomes in both experiments except chromosome 11 in 1st experiment ( Supplementary Table S6 ). Among these significant MTAs, 20 and 30 were found within the genic regions in both experiments respectively ( Supplementary Table S7 ). The phenotypic variance explained (PVE) by the MTAs ranged from 0.74% to 10.2% in 1st experiment and 0.00% to 9.57% in 2nd experiment. Furthermore, we examined the consistency of MTAs between the two experiments and identified 11 common MTAs. The soundness of the model is evident from the well-fitted Manhattan and Q-Q plots ( Supplementary Figure S5 ).

For tip root porosity, a total of 50 and 44 significant MTAs were identified from the analysis in both experiments, respectively. These MTAs were distributed across all chromosomes except chromosome 5 in both experiments ( Supplementary Table S8 ). Among these significant MTAs, 29 and 27 were found within the genic regions in both experiments respectively ( Supplementary Table S9 ). The phenotypic variance explained (PVE) by the MTAs ranged from 0.22% to 10.65% in 1st experiment and 1.01% to 16.24% in 2nd experiment. Furthermore, we examined the consistency of MTAs between the two experiments and identified 13 common MTAs. The soundness of the model is evident from the well-fitted Manhattan and Q-Q plots ( Supplementary Figure S6 ).

For Tiller number (TN), a total of 44 and 48 significant MTAs were identified from the analysis in the first and second experiments respectively, spanning all chromosomes ( Supplementary Table S10 ). Among these significant MTAs, 25 and 30 were located within genic regions in both experiments respectively ( Supplementary Table S11 ). The phenotypic variance explained (PVE) by the MTAs ranged from 0.9% to 9.1% in the first experiment and 0.67% to 9.4% in the second experiment. The robustness of the model is evidenced by the well-fitted Manhattan and Q-Q plots ( Supplementary Figure S7 ). Similarly, for root number (RN), a total of 55 and 52 significant MTAs were identified from the analysis in the first and second experiments respectively, spread across all chromosomes (Supplementary Table S12). Among these significant MTAs, 24 were present within genic regions in both the experiments (Supplementary Table S13). The phenotypic variance explained (PVE) by the MTAs ranged from 0.48% to 11.6% in the first experiment and 0.38% to 8.57% in the second experiment. The reliability of the model is apparent from the well-fitted Manhattan and Q-Q plots ( Supplementary Figure S8 ).

For root per tiller (RPT), a total of 39 and 56 significant MTAs identified from the analysis in 1st and 2nd experiments respectively. These MTAs were distributed across all chromosomes in both experiments except 3 and 9 chromosomes in 1st experiment ( Supplementary Table S14 ). Among these significant MTAs, 26 and 24 were found within the genic regions in both experiments respectively ( Supplementary Table S15 ). The phenotypic variance explained (PVE) by the MTAs ranged from 0.47% to 10.63% in 1st experiment and 0.34% to 8.25% in 2nd experiment. The soundness of the model is evident from the well-fitted Manhattan and Q-Q plots ( Supplementary Figure S9 ). For shoot dry weight (SDW), a total of 38 and 61 significant MTAs identified from the analysis in 1st and 2nd experiments respectively. These MTAs were distributed across all chromosomes in both experiments except chromosome 1 in 1st experiment ( Supplementary Table S16 ). Among these significant MTAs, 20 and 37 were found within the genic regions in both experiments respectively ( Supplementary Table S17 ). The phenotypic variance explained (PVE) by the MTAs ranged from 0.26% to 12.01% in 1st experiment and 3.8% to 11.55% in 2nd experiment. The soundness of the model is evident from the well-fitted Manhattan and Q-Q plots ( Supplementary Figure S10 ).

For root dry weight (RDW), a total of 41 and 49 significant MTAs identified from the analysis in 1st and 2nd experiments respectively. Both experiments distributed These MTAs across all chromosomes ( Supplementary Table S18 ). Among these significant MTAs, 21 and 32 were found within the genic regions in both experiments respectively ( Supplementary Table S19 ). The phenotypic variance explained (PVE) by the MTAs ranged from 1.33% to 9.23% in 1st experiment and 0.62% to 10.2% in 2nd experiment. Furthermore, we examined the consistency of MTAs between the two experiments and identified 2 common MTAs. The soundness of the model is evident from the well-fitted Manhattan and Q-Q plots ( Supplementary Figure S11 ). For root shoot ratio (RS Ratio), a total of 51 and 58 significant MTAs identified from the analysis in 1st and 2nd experiments respectively. Both experiments distributed These MTAs across all chromosomes ( Supplementary Table S20 ). Among these significant MTAs, 32 and 31 were found within the genic regions in both experiments respectively ( Supplementary Table S21 ). The phenotypic variance explained (PVE) by the MTAs ranged from 0% to 10.97% in 1st experiment and 0.04% to 10.13% in 2nd experiment. Furthermore, we examined the consistency of MTAs between the two experiments and identified 1 common MTAs. The soundness of the model is evident from the well-fitted Manhattan and Q-Q plots ( Supplementary Figure S12 ).

For root length (RL), a total of 65 and 60 significant MTAs identified from the analysis in 1st and 2nd experiments respectively. These MTAs were distributed across all chromosomes in both experiments ( Supplementary Table S22 ). Among these significant MTAs, 30 and 31 were found within the genic regions in both experiments respectively ( Supplementary Table S23 ). The phenotypic variance explained (PVE) by the MTAs ranged from 0.5% to 9.2% in 1st experiment and 0.91% to 10.4% in 2nd experiment. Furthermore, we examined the consistency of MTAs between the two experiments and identified 9 common MTAs. The soundness of the model is evident from the well-fitted Manhattan and Q-Q plots ( Supplementary Figure S13 ). For fine lateral root (FLR), a total of 39 and 51 significant MTAs identified from the analysis in 1st and 2nd experiments respectively. These MTAs were distributed across all chromosomes in both experiments ( Supplementary Table S24 ). Among these significant MTAs, 19 and 23 were found within the genic regions in both experiments respectively ( Supplementary Table S25 ). The phenotypic variance explained (PVE) by the MTAs ranged from 0.9% to 12.08% in 1st experiment and 0.6% to 9.02% in 2nd experiment. Furthermore, we examined the consistency of MTAs between the two experiments and identified 8 common MTAs. The soundness of the model is evident from the well-fitted Manhattan and Q-Q plots ( Supplementary Figure S14 ).

For thick lateral root (TLR), a total of 49 and 35 significant MTAs identified from the analysis in 1st and 2nd experiments respectively. These MTAs were distributed across all chromosomes in both experiments except chromosome 9 in 2nd experiment ( Supplementary Table S26 ). Among these significant MTAs, 31 and 20 were found within the genic regions in both experiments respectively ( Supplementary Table S27 ). The phenotypic variance explained (PVE) by the MTAs ranged from 0.96% to 9.8% in 1st experiment and 0.2% to 10.4% in 2nd experiment. Furthermore, we examined the consistency of MTAs between the two experiments and identified 8 common MTAs. The soundness of the model is evident from the well-fitted Manhattan and Q-Q plots ( Supplementary Figure S15 ).

For lateral root (LR), a total of 53 and 50 significant MTAs identified from the analysis in 1st and 2nd experiments respectively. These MTAs were distributed across all chromosomes in both experiments ( Supplementary Table S28 ). Among these significant MTAs, 31 and 32 were found within the genic regions in both experiments respectively ( Supplementary Table S29 ). The phenotypic variance explained (PVE) by the MTAs ranged from 0.75% to 9.6% in 1st experiment and 0.73% to 12.83% in 2nd experiment. Furthermore, we examined the consistency of MTAs between the two experiments and identified 6 common MTAs. The soundness of the model is evident from the well-fitted Manhattan and Q-Q plots ( Supplementary Figure S16 ).

For Nodal root (NR), a total of 48 and 55 significant MTAs identified from the analysis in 1st and 2nd experiments respectively. These MTAs were distributed across all chromosomes in both experiments except in chromosome 6 in 1st experiment ( Supplementary Table S30 ). Among these significant MTAs, 30 and 33 were found within the genic regions in both experiments respectively ( Supplementary Table S31 ). The phenotypic variance explained (PVE) by the MTAs ranged from 1.32% to 9.8% in 1st experiment and 1.2% to 9.9% in 2nd experiment. Furthermore, we examined the consistency of MTAs between the two experiments and identified 6 common MTAs. The soundness of the model is evident from the well-fitted Manhattan and Q-Q plots ( Supplementary Figure S17 ).

 Prediction and annotation of potential candidate genes

The analysis of 15 targeted traits across two experiments, examining the number of MTAs, genic MTAs, and expected PVE ranges is presented in Table 1 . ARP yielded 197 candidate genes from 48 significant MTAs in 1st experiment and 148 candidate genes from 41 significant MTAs in 2nd experiment. For BRP, 257 candidate genes emerged from 50 significant MTAs in the initial experiment, with a similar count of 255 candidate genes from 50 significant MTAs in the subsequent experiment. In the case of RDia, the first experiment identified 189 candidate genes from 45 significant MTAs, while the second experiment revealed 323 candidate genes from 60 significant MTAs. TRP showcased 203 candidate genes from 50 significant MTAs in the first experiment, contrasting with 155 candidate genes from 44 significant MTAs in the second experiment. In the first experiment, TN yielded 233 candidate genes from 44 significant MTAs, while the second experiment identified 189 candidate genes from 48 significant MTAs.

List of identified genic MTAs, candidate genes and PVE % for all traits.

Trait Experiment No of MTAs No of Candidate genes PVE (%)
Average Root Porosity Exp1 48 30 0.34 - 12.59
Average Root Porosity Exp2 41 25 0.20 - 12.55
Base Root Porosity Exp1 50 29 0.34 - 10.71
Base Root Porosity Exp2 50 29 0.23 - 09.67
Root Diameter Exp1 45 20 0.74 - 10.20
Root Diameter Exp2 60 30 0.00 - 09.57
Tip Root Porosity Exp1 50 29 0.22 - 10.65
Tip Root Porosity Exp2 44 27 1.01 - 16.24
Tiller Number Exp1 44 25 0.90 - 09.10
Tiller Number Exp2 48 30 0.67 - 09.40
Root Number Exp1 55 24 0.48 - 11.60
Root Number Exp2 52 24 0.38 - 08.57
Root Per Tiller Exp1 39 26 0.47 - 10.63
Root Per Tiller Exp2 56 54 0.34 - 08.25
Shoot Dry Weight Exp1 38 20 0.26 - 12.01
Shoot Dry Weight Exp2 61 37 3.82 - 11.55
Root Dry Weight Exp1 41 21 1.33 - 09.23
Root Dry Weight Exp2 49 32 0.62 - 10.20
Root Shoot Ratio Exp1 51 32 0.00 - 10.97
Root Shoot Ratio Exp2 58 31 0.04 - 10.13
Root Length Exp1 65 30 0.50 - 09.20
Root Length Exp2 60 31 0.91 - 10.40
Fine Lateral Root Exp1 39 19 0.90 - 12.08
Fine Lateral Root Exp2 51 23 0.60 - 09.02
Thick Lateral Exp1 49 31 0.96 - 09.80
Thick Lateral Exp2 35 20 0.20 - 10.40
Lateral Root Exp1 53 31 0.75 - 09.60
Lateral Root Exp2 50 32 0.73 - 12.86
Nodal Root Exp1 48 30 1.32 - 09.80
Nodal Root Exp2 55 33 1.20 - 09.90

Similarly, RN reveled 256 candidate genes from 55 significant MTAs in the first experiment and 236 candidate genes from 52 significant MTAs in the second experiment. RPT resulted in 193 candidate genes from 39 significant MTAs initially and 275 candidate genes from 56 significant MTAs subsequently. In the case of SDW, 172 candidate genes from 38 significant MTAs and 344 candidate genes from 61 significant MTAs were identified from the first and second experiments respectively. The RDW trait showed 181 candidate genes from 41 significant MTAs in the first experiment and 276 candidate genes from 49 significant MTAs in the second experiment. Further, RS ratio unveiled 268 candidate genes from 51 significant MTAs in the first experiment and 271 candidate genes from 58 significant MTAs in the second experiment. Lastly, RL exhibited 239 candidate genes from 65 significant MTAs in the first experiment and 258 candidate genes from 60 significant MTAs in the second experiment.

For FLR, the first experiment identified 187 candidate genes from 39 significant MTAs, while the second experiment reveled 253 candidate genes from 51 significant MTAs. In the case of TLR, the initial experiment uncovered 297 candidate genes from 49 significant MTAs, whereas the subsequent experiment yielded 165 candidate genes from 35 significant MTAs. LR analysis showcased 233 candidate genes from 53 significant MTAs and 205 candidate genes from 50 significant MTAs from the first and second experiments, respectively. Similarly, a thorough investigation revealed 203 candidate genes from 48 significant MTAs and 273 candidate genes from 55 significant MTAs associated with NR in the first and second experiments respectively. All candidate genes were annotated using the SnpEff software.

 Identification of superior haplotypes

To further explore the functional implications of these genetic loci, we investigated haplotype-level variation and gene expression dynamics among the identified candidate genes. A total of 219 candidate genes, 55,58,50, and 56 genes showed MTAs for the traits ARP, BRP, RDia, and TRP, respectively, in both experiments ( Table 1 ). Among these, 7 of these genes exhibited superior haplotypes while the remaining genes showed no variation across the 250 MAGIC lines. The distribution of haplotypes associated with target traits reveals a notable variation across different loci. A total of 8 superior haplotypes were identified for seven target genes associated with these traits. The distribution of lines for superior haplotype across the identified genes ranges from 38 to 76 in ARP, 6 to 76 in BRP,18 to 143 in RDia and 1 to 76 in TRP. Meanwhile, the range of frequency distribution of the targeted traits was as low as 0.03 in LOC_Os10g02080-H2 for BRP ( Figures 1C, D ), while as high as 0.77 in LOC_Os01g29250-H1 for RDia. Frequency distribution for the observed traits varied, with ARP ranging from 0.29 (LOC_Os01g39740-H5) to 0.44 (LOC_Os07g02630-H1), BRP from 0.44 (LOC_Os07g02630-H1) to 0.03 (LOC_Os10g02080-H2), RDia from 0.11 (LOC_Os11g40860-H1) to 0.77 (LOC_Os01g29250-H1) and for TRP 0.44 (LOC_Os07g02630-H1) ( Supplementary Table S32 ). The diversity analysis of haplotypes across identified MTA can provide valuable insight into the genetic variation of the population.

(A) Phenotypic diversity of base root porosity across experiments. (B) Manhattan plot for root diameter. (C) Haplotype diversity of gene LOC_Os10g02080 associated with base root porosity. (D) The box plot depicts the base root porosity variation among the haplotypes of gene LOC_Os10g02080.

Panel A shows violin plots comparing TRP values for two experiments (Exp1 and Exp2), alongside histograms. Panel B presents Manhattan plots with P-values mapped across chromosomes. Panel C displays a pie chart illustrating haplotype distribution, with H1 to H3 represented in various colors. Panel D features a box plot of BRP_E2 values across haplotypes for the gene LOC_Os10g02080. Network analysis of candidate genes associated with rice roots

Expression of identified candidate genes using RiceXPro (https://ricexpro.dna.affrc.go.jp/) (Sato et al., 2013) in various tissues ( Figure 2 ) showed 16 differentially expressed genes in root tissues. This gene set consists of many genes that are involved in plant growth and development, such as U-box containing protein, basic Leucine Zipper Domain (bZIP) containing transcription factor, F-box proteins, and cytochrome P450 genes (Abd-Hamid et al., 2020; Mao et al., 2022; Han et al., 2023), and the details of the identified genes are given in Supplementary Table S33 . We have performed a protein-protein interaction network analysis of these genes, and Figure 3 shows the interaction network as determined by the STRING database (Szklarczyk et al., 2023) (version 11.0; https://string-db.org/). We observed 3 main clusters with three auxiliary networks that consist of 9 genes from our queries. The main cluster has 4 of the query genes and though most of the genes in this cluster are involved in the porphyrin and chlorophyll metabolism some genes belong to the Chalcone/stilbene synthases family genes (A0A0P0XSC7: Os10g0177300, Q6L4H9_ORYSJ: Os05g0212500, and A0A0N7KKC2: Os05g0212600) that are known for their involvement in flavonoid biosynthesis and in the regulation of auxin transport and in the root gravitropism (Hu et al., 2021). Two other clusters (2 & 4) mainly consist of putative genes, and they are found to be interacting with one of the peptide transferase proteins (Os05g0410900) we identified. The second cluster has genes that are all involved in the N-glycan biosynthesis pathways, which have an effect on the cellulose content (Strasser, 2014). The genes in this cluster are interacting with one of the candidate genes, Os05g0411200 (OsTCD5). It is a thermosensitive chloroplast development gene. The details of the genes identified in this network and the interaction score for this network are given in Supplementary Tables S34 , S35 , respectively.

Tissue-specific expression profiles of candidate genes associated with root traits using RiceXPro. Expression patterns of 16 candidate genes identified through GWAS were examined across various rice tissues. The heatmap illustrates gene expression intensity based on median-centered values generated from RiceXPro, with a focus on root tissue.

Heatmap showing gene expression levels with hierarchical clustering. Rows represent gene IDs, columns represent different plant tissues and stages. Colors range from blue (low expression) to red (high expression), with a color key from -5 to 5.

Protein-protein interaction (PPI) network of 16 candidate genes involved in root development. The network consists of three primary clusters and three auxiliary subnetworks.

Network diagram illustrating protein interactions. Four clusters are shown, with Cluster 1 linked to porphyrin and chlorophyll metabolism and Clusters 2 and 3 associated with N-glycan biosynthesis. Interaction types include known interactions (curated databases, experimentally determined) and predicted interactions (gene neighborhood, fusions, co-occurrence). Other connections involve text mining, co-expression, and protein homology. Each node represents a protein.

To determine how the identified candidate genes interact with reported root trait-related genes, we retrieved 26 proteins already known to be involved in root development and functionality from the literature (Wimalagunasekara et al., 2023) and the oryzabase database (https://shigen.nig.ac.jp/rice/oryzabase/). The details of these reported genes are given in Supplementary Table S36 . The combined input of the 26 reported genes and the 16 candidate genes, resulted in a fully connected network with 3 distinct clusters ( Figure 4 ). Four of the candidate genes were found to be interacting with 11 of the reported genes. All four of the candidate genes are present in the first cluster, which contains genes involved in three major pathways (Peroxisome, Pantothenate and CoA biosynthesis, and Cutin, suberine and wax biosynthesis). The second cluster has two reported genes which are involved in the N-glycan biosynthesis pathway. The third cluster contains 4 reported genes. Based on the genes present in the closed cluster, we can interpret that the 4 candidate genes might have a role in root branching and/or the number of adventitious and lateral roots. The details of the genes identified in this network and the interaction score for this network are given in Supplementary Tables S37 , S38 , respectively. Interestingly, one gene (Os04g0405300) among the 4 candidate genes showed a significant difference in auxin and Jasmonic acid production with respect to the different periods of time ( Figure 5 ). This gene regulates alcohol dehydrogenase/reductase, which is crucial for maintaining homeostasis of active auxin levels within the plants.

Combined protein-protein interaction network of candidate genes and reported root development genes.

A network diagram shows three clusters of genes involved in various biosynthetic pathways. Cluster 1 focuses on peroxisome, pantothenate and CoA biosynthesis, and cutin, suberine, and wax biosynthesis. Cluster 2 is related to N-glycan biosynthesis. Cluster 3 contains additional gene interactions. Nodes are color-coded to indicate known and predicted interactions, with specific genes highlighted in red or black boxes to denote reported and candidate genes. A legend explains interaction types such as curated databases, experimental determination, gene neighborhood, fusions, and co-occurrence.

Temporal expression pattern of Os04g0405300 in response to hormonal cues (Sato et al., 2013).

Line graphs display the Cy5/Cy3 ratio (log2) over time for different hormones in root and shoot tissues. Root graphs show data at 15 minutes to 6 hours, while shoot graphs range from 1 to 12 hours. Hormones include abscisic acid, gibberellin, auxin, brassinosteroid, cytokinin, and jasmonic acid. All graphs show relatively stable trends with some variations. Discussion

O2 plays a regulatory role in CH4 production, oxidation, and emission from flooded rice fields (Colmer, 2003; Yang et al., 2012; Fazli et al., 2013). Atmospheric O2 reaches the rhizospheric soil through the well-developed rice aerenchyma system (Armstrong, 1971, 1980). In most cases, O2 diffusion into the soil depends on the root’s radial oxygen loss capacity and porosity (Luxmoore et al., 1970). Root porosity generally depends on aerenchyma formation and root diameter (Armstrong, 1980; Jeffree and Fry, 1986). These factors are crucial because the root acts as a conduit, holding and transporting various gases (O2, CH4, CO2, N2O) between the atmosphere and the soil. Recent studies continue to emphasize the importance of these traits in managing greenhouse gas emissions from rice paddies (Mei et al., 2009; Jiang et al., 2017; Zheng et al., 2018).

Our study enabled the identification of essential QTNs and genes located within these QTNs. We evaluated the Global Multi-parent Advanced Generation Inter-Cross (MAGIC) population, focusing on four traits related to root porosity and diameter: RDia, BP, TP, and ARP. Notably, significant variations were observed among these traits and genotypes. This diversity could be attributed to the use of diverse genotypes derived from multiple elite parents across the globe. The Global MAGIC population, generated from eight indica and eight japonica patents, integrates multiple traits from both gene pools, which have adapted to diverse ecotypes worldwide. Recombination among the 16 genotypes resulted in novel genetic variation. Key advantages of the global MAGIC population include ample genetic diversity, a higher allelic balance frequency compared to biparental populations, and minimal impact on population structure (Valdar et al., 2006; Mackay and Powell, 2007; Kover et al., 2009). Additionally, it enhances mapping accuracy through historical and synthetic recombination (Meng et al., 2016).

The evaluation of root traits across two experiments using the Global MAGIC population revealed a consistent pattern of phenotypic variation, reinforcing the reliability of the observed trait distribution. However, the minor variation detected, particularly in root number and root length, reflects the plastic nature of root architecture, which is likely due to environmental influence. In our study, we observed significant correlations among the traits related to root porosity: basal root porosity, total root porosity, and apical root porosity. These porosity traits were strongly associated with root diameter (RDia). Previous research by Striker et al. (2007) reported that increased root porosity is directly dependent on larger root diameter. Additionally, we found that BRP, TRP, and ARP were closely linked. Root porosity is intricately tied to the formation of root aerenchyma. In rice, lysigenous aerenchyma forms in the cortex through cell death and subsequent lysis (Armstrong, 1980; Jackson et al., 1985). The process begins in mid-cortex cells and spreads radially to surrounding cortical cells (Kawai et al., 1998). Aerenchyma formation initiates at the apical parts of rice roots and gradually extends to the basal regions (Ranathunge et al., 2003). Ultimately, fully developed aerenchyma is observed at the basal part of the root (Armstrong and Armstrong, 1994; Kozela and Regan, 2003; Ranathunge et al., 2003). Despite the remarkable difference in BRP compared to TRP, these traits remain closely associated. These factors had a direct or indirect positive effect on the amount of oxygen diffused into the rhizosphere, promoting methane oxidation.

We utilized 250 global MAGIC population lines, each characterized by 3,53,641 Single SNPs. Our goal was to identify QTNs and candidate genes associated with root diameter and root porosity. To achieve this, multi-locus GWAS models outperform single-locus models statistically, with a lower False Positive Rate (FPR) (Segura et al., 2012; Wang et al., 2016; Zhang et al., 2019). These models closely resemble true genetic architectures in plants and animals. By incorporating polygenic effects and accounting for population structure, multi-locus models mitigate bias in effect estimation (Zhang et al., 2005, 2010; Yu et al., 2006). The five multi-locus GWAS methods in this study is an integration of FASTmrEMMA (Wen et al., 2018), ISIS EM-BLASSO (Tamba et al., 2017) mrMLM (Wang et al., 2016), pLARmEB (Zhang et al., 2017) and pKWmEB (Ren et al., 2018).

The present study has identified a substantial number of MTAs across key root traits in rice, with candidate genes ranging from 19 to 54 for each trait ( Table 1 ). PVE values were generally low to moderate (0.2% to 16%), indicating that these traits are influenced by many small-effect loci rather than a few significant genes. We filtered down from 274 candidate genes for Rdia, ARP, BRP and TRP to 16 candidate genes based on the expression of the identified genes in root tissues. The protein-protein interaction of these genes revealed three primary clusters ( Figure 3 ). The primary cluster includes key genes regulating chlorophyll metabolism, flavonoid biosynthesis, and affecting auxin transport. While another cluster, 2 and 4, has putative genes interacting with peptide transferase. Combining GWAS with protein-protein analysis enhances the validation of identified MTAs and deepens the understanding of the molecular framework regulating traits of interest.

 Conclusion

This study provides critical insights into the genetic architecture of rice root traits that can mitigate methane emissions in paddy ecosystems. By utilizing a genetically diverse Global MAGIC population and GWAS, we identified several MTAs and candidate genes associated with root diameter and porosity, two traits directly influencing root-mediated oxygen transport and CH4 oxidation. Our haplotype-phenotype analysis identified eight superior haplotypes across seven key genes, with wide variation in frequency, confirming substantial allelic diversity within the population. Furthermore, network analysis of 16 root-expressed candidate genes revealed their functional integration with known root developmental pathways, including auxin regulation, flavonoid biosynthesis, and cell wall modification, reinforcing their biological significance. Together, these results establish a comprehensive framework for selecting and developing rice varieties with improved root traits and lower greenhouse gas emissions, paving the way for sustainable and climate-smart rice cultivation.

 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 in the article/ Supplementary Material .

 Author contributions

RR: Conceptualization, Methodology, Data curation, Writing – review & editing, Writing – original draft, Formal analysis. GM: Formal analysis, Data curation, Writing – review & editing. SS: Data curation, Formal analysis, Writing – review & editing, Visualization, Methodology. BP: Methodology, Formal analysis, Writing – review & editing, Data curation, Resources. SH: Methodology, Writing – review & editing, Resources, Visualization. KT: Resources, Writing – review & editing, Methodology. SK: Writing – review & editing, Resources, Methodology. JH: Writing – review & editing, Methodology, Resources. AH: Resources, Writing – review & editing. NS: Writing – review & editing, Resources. MD: Methodology, Investigation, Writing – review & editing, Resources. EO: Resources, Methodology, Writing – review & editing. PS: Resources, Formal analysis, Writing – original draft, Investigation, Conceptualization, Validation, Methodology, Supervision. AK: Resources, Project administration, Conceptualization, Validation, Visualization, Methodology, Funding acquisition, Supervision, Writing – review & editing, Formal analysis, Writing – original draft, Investigation, Data curation.

 Funding

The author(s) declare that financial support was received for the research and/or publication of this article. The author acknowledges the financial support of the Indian Council of Agricultural Research (ICAR) under the IRRI-ICAR Work Plan and National Agricultural Technology Program-Phase II Project (NATP-2) BARC, Bangladesh.

 Acknowledgments

The authors express their sincere gratitude to the IRRI authorities for their valuable support in providing greenhouse facilities for the experiments. We also appreciate the help from the Plant Molecular Biology Laboratory and the Plant Physiology Laboratory at IRRI for their contributions to the experimental work and data processing. We would also like to acknowledge the NATP-2 project for providing a scholarship that supported the author during the course of the experiment.

 Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

 Correction note

A correction has been made to this article. Details can be found at: 10.3389/fpls.2025.1692969.

 Generative AI statement

The author(s) declare that no Generative AI was used in the creation of this manuscript.

 Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

 Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fpls.2025.1616424/full#supplementary-material

Phenotypic variation and correlation among traits in experiment 1. The correlation plot compares Average Root Porosity (ARP), Base Root Porosity (BRP), Root Diameter (RDia), Tip Root Porosity (TRP), Tiller Number (TN), Root Number (RN), Root Per Tiller (RPT), Shoot Dry Weight (SDW), Root Dry Weight (RDW), Ratio Between Root Shoot (RS ratio), Root Length (RL), Fine Lateral Root (FLR), Thick Lateral Root (TLR), Lateral Root (LR) and Nodal Root (NR). Bar charts show the frequency distributions of each trait on the diagonal. Pearson correlation values are above the diagonal with their corresponding P-values, and scatter plots are below the diagonal. Data correspond to mean values of 250 MAGIC population in two different experiments (a, b). Symbols indicate the significance levels: * **P < 0.01 and ***P < 0.005.

Phenotypic variation and correlation among traits in experiment 2. The correlation plot compares Average Root Porosity (ARP), Base Root Porosity (BRP), Root Diameter (RDia), Tip Root Porosity (TRP), Tiller Number (TN), Root Number (RN), Root Per Tiller (RPT), Shoot Dry Weight (SDW), Root Dry Weight (RDW), Ratio Between Root Shoot (RS ratio), Root Length (RL), Fine Lateral Root (FLR), Thick Lateral Root (TLR), Lateral Root (LR) and Nodal Root (NR). Bar charts show the frequency distributions of each trait on the diagonal. Pearson correlation values are above the diagonal with their corresponding P-values, and scatter plots are below the diagonal. Data correspond to mean values of 250 MAGIC population in two different experiments (a, b). Symbols indicate the significance levels: * **P < 0.01 and ***P < 0.005.

(A) Phenotypic diversity of population structure and average root porosity of 250 accessions. (B) Manhattan plot for root length (C) Q-Q plot for average root porosity.

(A) Phenotypic diversity of population structure and base root porosity of 250 accessions. (B) Manhattan plot for base root porosity (C) Q-Q plot for base root porosity.

(A) Phenotypic diversity of population structure and root diameter of 250 accessions. (B) Manhattan plot for root diameter (C) Q-Q plot for root diameter.

(A) Phenotypic diversity of population structure and tip root porosity of 250 accessions. (B) Manhattan plot for tip root porosity (C) Q-Q plot for tip root porosity.

(A) Phenotypic diversity of population structure and tiller number of 250 accessions. (B) Manhattan plot for tiller number (C) Q-Q plot for tiller number.

(A) Phenotypic diversity of population structure and root number of 250 accessions. (B) Manhattan plot for root number (C) Q-Q plot for root number.

(A) Phenotypic diversity of population structure and root per tiller of 250 accessions. (B) Manhattan plot for root per tiller (C) Q-Q plot for root per tiller.

(A) Phenotypic diversity of population structure and shoot dry weight of 250 accessions. (B) Manhattan plot for shoot dry weight (C) Q-Q plot for shoot dry weight.

(A) Phenotypic diversity of population structure and root dry weight of 250 accessions. (B) Manhattan plot for root dry weight (C) Q-Q plot for root dry weight.

(A) Phenotypic diversity of population structure and root shoot ratio of 250 accessions. (B) Manhattan plot for root shoot ratio (C) Q-Q plot for root shoot ratio.

(A) Phenotypic diversity of population structure and root length of 250 accessions. (B) Manhattan plot for root length (C) Q-Q plot for root length.

(A) Phenotypic diversity of population structure and fine lateral root of 250 accessions. (B) Manhattan plot for fine lateral root (C) Q-Q plot for fine lateral root.

(A) Phenotypic diversity of population structure and thick lateral root of 250 accessions. (B) Manhattan plot for thick lateral root (C) Q-Q plot for thick lateral root.

(A) Phenotypic diversity of population structure and lateral root of 250 accessions. (B) Manhattan plot for lateral root (C) Q-Q plot for lateral root.

(A) Phenotypic diversity of population structure and nodal root of 250 accessions. (B) Manhattan plot for nodal root (C) Q-Q plot for nodal root.

Co-expression analysis of genes associated with root traits.

 References Abd-Hamid N.-A. Ahmad-Fauzi M.-I. Zainal Z. Ismail I. (2020). Diverse and dynamic roles of F-box proteins in plant biology. Planta 251, 68. doi: 10.1007/s00425-020-03356-8, PMID: 32072251 Abeysekara I. Rathnayake I. (2024). Global trends in rice production, consumption and trade. doi: 10.2139/ssrn.4948477 Armstrong J. Armstrong W. (1994). Chlorophyll development in mature lysigenous and schizogenous root aerenchymas provides evidence of continuing cortical cell viability. New Phytol. 126, 493–497. doi: 10.1111/j.1469-8137.1994.tb04246.x, PMID: 33874471 Armstrong W. (1971). Radial oxygen losses from intact rice roots as affected by distance from the apex, respiration and waterlogging. Physiologia Plantarum 25, 192–197. doi: 10.1111/j.1399-3054.1971.tb01427.x Armstrong W. (1980). “Aeration in Higher Plants,” in Advances in Botanical Research. Ed. Woolhouse H. W. (The University of Hull, Hull, England: Academic Press), 225–332. doi: 10.1016/S0065-2296(08)60089-0 Bandillo N. Raghavan C. Muyco P. A. Sevilla M. A. L. Lobina I. T. Dilla-Ermita C. J. . (2013). Multi-parent advanced generation inter-cross (MAGIC) populations in rice: progress and potential for genetics research and breeding. Rice (N Y) 6, 11. doi: 10.1186/1939-8433-6-11, PMID: 24280183 Browning B. L. Zhou Y. Browning S. R. (2018). A one-penny imputed genome from next-generation reference panels. Am. J. Hum. Genet. 103, 338–348. doi: 10.1016/j.ajhg.2018.07.015, PMID: 30100085 Colmer T. (2003). Long-distance transport of gases in plants: aperspective on internal aeration and radial oxygen loss from roots. Plant Cell Environ. 26, 17–36. doi: 10.1046/j.1365-3040.2003.00846.x Davamani V. Parameswari E. Arulmani S. (2020). Mitigation of methane gas emissions in flooded paddy soil through the utilization of methanotrophs. Sci. Total Environ. 726, 138570. doi: 10.1016/j.scitotenv.2020.138570, PMID: 32305766 de Mendiburu F. (2019). Agricolae: Statistical Procedures for Agricultural Research. R Package version 1.3-1. Available online at: https://CRAN.R-project.org/package=agricolae (Accessed on June 10, 2024). Ejiri M. Fukao T. Miyashita T. Shiono K. (2021). A barrier to radial oxygen loss helps the root system cope with waterlogging-induced hypoxia. Breed Sci. 71, 40–50. doi: 10.1270/jsbbs.20110, PMID: 33762875 Fazli P. Man H. C. Shah U. K. Idris A. (2013). Characteristics of methanogens and methanotrophs in rice fields: A review. Asia-Pacific Journal of Molecular Biology and Biotechnology. 21, 3–17. Han H. Wang C. Yang X. Wang L. Ye J. Xu F. . (2023). Role of bZIP transcription factors in the regulation of plant secondary metabolism. Planta 258, 13. doi: 10.1007/s00425-023-04174-4, PMID: 37300575 Hertel T. W. (2011). The global supply and demand for agricultural land in 2050: A perfect storm in the making? Am. J. Agric. Economics 93, 259–275. doi: 10.1093/ajae/aaq189 Hu Y. Omary M. Hu Y. H. Doron O. Hoermayer L. Chen Q. . (2021). Cell kinetics of auxin transport and activity in Arabidopsis root growth and skewing. Nat. Commun. 12, 1657. doi: 10.1038/s41467-021-21802-3, PMID: 33712581 Husson O. (2013). Redox potential (Eh) and pH as drivers of soil/plant/microorganism systems: a transdisciplinary overview pointing to integrative opportunities for agronomy. Plant Soil 362, 389–417. doi: 10.1007/s11104-012-1429-7 Ingvarsson P. K. Street N. R. (2011). Association genetics of complex traits in plants. New Phytol. 189, 909–922. doi: 10.1111/j.1469-8137.2010.03593.x, PMID: 21182529 Jackson M. B. Fenning T. M. Drew M. C. Saker L. R. (1985). Stimulation of ethylene production and gas-space (aerenchyma) formation in adventitious roots of Zea mays L. by small partial pressures of oxygen. Planta 165, 486–492. doi: 10.1007/BF00398093, PMID: 24241221 Jeffree C. E. Fry S. C. (1986). The genesis of intercellular spaces in developing leaves of Phaseolus vulgaris L. Protoplasma 132, 90–98. doi: 10.1007/BF01275795 Jiang Y. van Groenigen K. J. Huang S. Hungate B. A. van Kessel C. Hu S. . (2017). Higher yields and lower methane emissions with new rice cultivars. Glob Chang Biol. 23, 4728–4738. doi: 10.1111/gcb.13737, PMID: 28464384 Jiang Y. Wang L. Yan X. Tian Y. Deng A. Zhang W. (2013). Super rice cropping will enhance rice yield and reduce CH4 emission: A case study in Nanjing, China. Rice Sci. 20, 427–433. doi: 10.1016/S1672-6308(13)60157-2 Jung Y. Han D. (2022). BWA-MEME: BWA-MEM emulated with a machine learning approach. Bioinformatics 38, 2404–2413. doi: 10.1093/bioinformatics/btac137, PMID: 35253835 Kawai M. Samarajeewa P. K. Barrero R. A. Nishiguchi M. Uchimiya H. (1998). Cellular dissection of the degradation pattern of cortical cell death during aerenchyma formation of rice roots. Planta 204, 277–287. doi: 10.1007/s004250050257 Kim Y. X. Ranathunge K. Lee S. Lee Y. Lee D. Sung J. (2018). Composite transport model and water and solute transport across plant roots: an update. Front. Plant Sci. 9. doi: 10.3389/fpls.2018.00193, PMID: 29503659 Kirk G. J. D. (2003). Rice root properties for internal aeration and efficient nutrient acquisition in submerged soil. New Phytol. 159, 185–194. doi: 10.1046/j.1469-8137.2003.00793.x, PMID: 33873689 Kover P. X. Valdar W. Trakalo J. Scarcelli N. Ehrenreich I. M. Purugganan M. D. . (2009). A Multiparent Advanced Generation Inter-Cross to fine-map quantitative traits in Arabidopsis thaliana. PloS Genet. 5, e1000551. doi: 10.1371/journal.pgen.1000551, PMID: 19593375 Kozela C. Regan S. (2003). How plants make tubes. Trends Plant Sci. 8, 159–164. doi: 10.1016/S1360-1385(03)00050-5, PMID: 12711227 Kusin K. Ogawa K. Doi H. Tokida T. Hirano T. Jaya A. . (2025). Increasing the pH of tropical peat can enhance methane production and methanogenic growth under anoxic conditions. CATENA 250, 108791. doi: 10.1016/j.catena.2025.108791 Luo D. Li Y. Yao H. Chapman S. J. (2022). Effects of different carbon sources on methane production and the methanogenic communities in iron rich flooded paddy soil. Sci. Total Environ. 823, 153636. doi: 10.1016/j.scitotenv.2022.153636, PMID: 35124061 Luxmoore R. J. Stolzy L. H. Letey J. (1970). Oxygen diffusion in the soil-plant system IV. Oxygen concentration profiles, respiration rates, and radial oxygen losses predicted for rice roots. Agron. J. 62, 329–332. doi: 10.2134/agronj1970.00021962006200030006x Mackay I. Powell W. (2007). Methods for linkage disequilibrium mapping in crops. Trends Plant Sci. 12, 57–63. doi: 10.1016/j.tplants.2006.12.001, PMID: 17224302 Mano Y. Omori F. Takamizo T. Kindiger B. Bird R.M. Loaisiga C. H. . (2007). QTL mapping of root aerenchyma formation in seedlings of a maize × rare teosinte “Zea nicaraguensis“ cross. Plant Soil 295, 103–113. doi: 10.1007/s11104-007-9266-9 Mao X. Yu C. Li L. Wang M. Yang L. Zhang Y. . (2022). How many faces does the plant U-box E3 ligase have? Int. J. Mol. Sci. 23, 2285. doi: 10.3390/ijms23042285, PMID: 35216399 McKenna A. Hanna M. Banks E. Sivachenko A. Cibulskis K. Kernytsky A. . (2010). The Genome Analysis Toolkit: a MapReduce framework for analyzing next-generation DNA sequencing data. Genome Res. 20, 1297–1303. doi: 10.1101/gr.107524.110, PMID: 20644199 Mei X. Q. Ye Z. H. Wong M. H. (2009). The relationship of root porosity and radial oxygen loss on arsenic tolerance and uptake in rice grains and straw. Environ. pollut. 157, 2550–2557. doi: 10.1016/j.envpol.2009.02.037, PMID: 19329236 Meng L. Guo L. Ponce K. Zhao X. Ye G. (2016). Characterization of three rice multiparent advanced generation intercross (MAGIC) populations for quantitative trait loci identification. Plant Genome 9. doi: 10.3835/plantgenome2015.10.0109, PMID: 27898805 Mühlenbock P. Plaszczyca M. Plaszczyca M. Mellerowicz E. Karpinski S. (2007). Lysigenous aerenchyma formation in Arabidopsis is controlled by LESION SIMULATING DISEASE1. Plant Cell 19, 3819–3830. doi: 10.1105/tpc.106.048843, PMID: 18055613 Ponnamperuma F. N. (1972). “The Chemistry of Submerged Soils,” in Advances in Agronomy. Ed. Brady N. C. (Los Baños, Laguna, Philippines: Academic Press), 29–96. doi: 10.1016/S0065-2113(08)60633-1 Purcell S. Neale B. Todd-Brown K. Thomas L. Ferreira M. A. Bender D. . (2007). PLINK: a tool set for Whole-Genome association and Population-Based linkage analyses. Am. J. Hum. Genet. 81, 559–575. doi: 10.1086/519795, PMID: 17701901 Rago L. Ruiz Y. Baeza J. A. Guisasola A. Cortés P. (2015). Microbial community analysis in a long-term membrane-less microbial electrolysis cell with hydrogen and methane production. Bioelectrochemistry 106, 359–368. doi: 10.1016/j.bioelechem.2015.06.003, PMID: 26138343 Ranathunge K. Steudle E. Lafitte R. (2003). Control of water uptake by rice (Oryza sativa L.): role of the outer part of the root. Planta 217, 193–205. doi: 10.1007/s00425-003-0984-9, PMID: 12783327 Raven J. A. (1996). Into the voids: the distribution, function, development and maintenance of gas spaces in plants. Ann. Bot. 78, 137–142. doi: 10.1006/anbo.1996.0105 Reim A. Lüke C. Krause S. Pratscher J. Frenzel P. (2012). One millimetre makes the difference: high-resolution analysis of methane-oxidizing bacteria and their specific activity at the oxic–anoxic interface in a flooded paddy soil. ISME J. 6, 2128–2139. doi: 10.1038/ismej.2012.57, PMID: 22695859 Ren W.-L. Wen Y.-J. Dunwell J. M. Zhang Y.-M. (2018). pKWmEB: integration of Kruskal-Wallis test with empirical Bayes under polygenic background control for multi-locus genome-wide association study. Heredity (Edinb) 120, 208–218. doi: 10.1038/s41437-017-0007-4, PMID: 29234158 Sahrawat K. L. (2005). Iron toxicity in wetland rice and the role of other nutrients. J. Plant Nutr. 27, 1471–1504. doi: 10.1081/PLN-200025869 Sapkota T. B. Vetter S. H. Jat M. L. Sirohi S. Shirsath P. B. Singh R. . (2019). Cost-effective opportunities for climate change mitigation in Indian agriculture. Sci. Total Environ. 655, 1342–1354. doi: 10.1016/j.scitotenv.2018.11.225, PMID: 30577126 Sato Y. Antonio A. B. Namiki N. Takehisa H. Minami H. Kamatsuki K. . (2011). RiceXPro: a platform for monitoring gene expression in japonica rice grown under natural field conditions. Nucleic Acids Res. 39, D1141–D1148. doi: 10.1093/nar/gkq1085, PMID: 21045061 Sato Y. Takehisa H. Kamatsuki K. Minami H. Namiki N. Ikawa H. . (2013). RiceXPro Version 3.0: expanding the informatics resource for rice transcriptome. Nucleic Acids Res. 41, D1206–D1213. doi: 10.1093/nar/gks1125, PMID: 23180765 Segura V. Vilhjálmsson B. J. Platt A. Korte A. Seren Ü. Long Q. . (2012). An efficient multi-locus mixed-model approach for genome-wide association studies in structured populations. Nat. Genet. 44, 825–830. doi: 10.1038/ng.2314, PMID: 22706313 Sinha P. Singh V. K. Saxena R. K. Khan A. W. Abbai R. Chitikineni A. . (2020). Superior haplotypes for haplotype-based breeding for drought tolerance in pigeonpea (Cajanus cajan L.). Plant Biotechnol. J. 18, 2482–2490. doi: 10.1111/pbi.13422, PMID: 32455481 Strasser R. (2014). Biological significance of complex N-glycans in plants and their impact on plant physiology. Front. Plant Sci. 5. doi: 10.3389/fpls.2014.00363, PMID: 25101107 Striker G. G. Insausti P. Grimoldi A. A. Vega A. S. (2007). Trade-off between root porosity and mechanical strength in species with different types of aerenchyma. Plant Cell Environ. 30, 580–589. doi: 10.1111/j.1365-3040.2007.01639.x, PMID: 17407536 Szklarczyk D. Kirsch R. Koutrouli M. Nastou K. Mehryary F. Hachilif R. . (2023). The STRING database in 2023: protein-protein association networks and functional enrichment analyses for any sequenced genome of interest. Nucleic Acids Res. 51, D638–D646. doi: 10.1093/nar/gkac1000, PMID: 36370105 Tamba C. L. Ni Y.-L. Zhang Y.-M. (2017). Iterative sure independence screening EM-Bayesian LASSO algorithm for multi-locus genome-wide association studies. PloS Comput. Biol. 13, e1005357. doi: 10.1371/journal.pcbi.1005357, PMID: 28141824 Valdar W. Flint J. Mott R. (2006). Simulating the collaborative cross: power of quantitative trait loci detection and mapping resolution in large sets of recombinant inbred strains of mice. Genetics 172, 1783–1797. doi: 10.1534/genetics.104.039313, PMID: 16361245 Visser E. J. W. Bögemann G. M. (2003). Measurement of porosity in very small samples of plant tissue. Plant Soil 253, 81–90. doi: 10.1023/A:1024560322835 Wang S.-B. Feng J.-Y. Ren W.-L. Huang B. Zhou L. Wen Y.-J. . (2016). Improving power and accuracy of genome-wide association studies via a multi-locus mixed linear model methodology. Sci. Rep. 6, 19444. doi: 10.1038/srep19444, PMID: 26787347 Wen Y.-J. Zhang H. Ni Y.-L. Huang B. Zhang J. Feng J.-Y. . (2018). Methodological implementation of mixed linear models in multi-locus genome-wide association studies. Brief Bioinform. 19, 700–712. doi: 10.1093/bib/bbw145, PMID: 28158525 Wimalagunasekara S. S. Weeraman J. W. J. K. Tirimanne S. Fernando P. C. (2023). Protein-protein interaction (PPI) network analysis reveals important hub proteins and sub-network modules for root development in rice (Oryza sativa). J. Genet. Eng. Biotechnol. 21, 69. doi: 10.1186/s43141-023-00515-8, PMID: 37246172 Win K. T. Nonaka R. Win A. T. Sasada Y. Toyota K. Motobayashi T. . (2011). Comparison of methanotrophic bacteria, methane oxidation activity, and methane emission in rice fields fertilized with anaerobically digested slurry between a fodder rice and a normal rice variety. Paddy Water Environ. 10, 281–289. doi: 10.1007/s10333-011-0279-x Xiao Y. Liu H. Wu L. Warburton M. Yan J. (2017). Genome-wide association studies in maize: praise and stargaze. Mol. Plant 10, 359–374. doi: 10.1016/j.molp.2016.12.008, PMID: 28039028 Xu H. Cai Z. C. Tsuruta H. (2003). Soil moisture between rice-growing seasons affects methane emission, production, and oxidation. Soil Sci. Soc. America J. 67, 1147–1157. doi: 10.2136/sssaj2003.1147 Yang J.-X. Liu Y. Ye Z.-H. (2012). Root-induced changes of pH, eh, fe(II) and fractions of pb and zn in rhizosphere soils of four wetland plants with different radial oxygen losses. Pedosphere 22, 518–527. doi: 10.1016/S1002-0160(12)60036-8 Yu J. Pressoir G. Briggs W. H. Vroh Bi I. Yamasaki M. Doebley J. F. . (2006). A unified mixed-model method for association mapping that accounts for multiple levels of relatedness. Nat. Genet. 38, 203–208. doi: 10.1038/ng1702, PMID: 16380716 Zhang Z. Ersoz E. Lai C.-Q. Todhunter R. J. Tiwari H. K. Gore M. A. . (2010). Mixed linear model approach adapted for genome-wide association studies. Nat. Genet. 42, 355–360. doi: 10.1038/ng.546, PMID: 20208535 Zhang J. Feng J.-Y. Ni Y.-L. Wen Y.-J. Niu Y. Tamba C. L. . (2017). pLARmEB: integration of least angle regression with empirical Bayes for multilocus genome-wide association studies. Heredity 118, 517–524. doi: 10.1038/hdy.2017.8, PMID: 28295030 Zhang Y.-M. Jia Z. Dunwell J. M. (2019). Editorial: the applications of new multi-locus GWAS methodologies in the genetic dissection of complex traits. Front. Plant Sci. 10. doi: 10.3389/fpls.2019.00100, PMID: 30804969 Zhang Y.-M. Mao Y. Xie C. Smith H. Luo L. Xu S. (2005). Mapping quantitative trait loci using naturally occurring genetic variance among commercial inbred lines of maize (Zea mays L.). Genetics 169, 2267–2275. doi: 10.1534/genetics.104.033217, PMID: 15716509 Zheng G. Freidlin B. Li Z. Gastwirth J. L. (2005). Genomic control for association studies under various genetic models. Biometrics 61, 186–192. doi: 10.1111/j.0006-341X.2005.t01-1-.x, PMID: 15737092 Zheng H. Fu Z. Zhong J. Long W. (2018). Low methane emission in rice cultivars with high radial oxygen loss. EBSCOhost. 431, 119. doi: 10.1007/s11104-018-3747-x