Genome-wide association studies for a comprehensive understanding of the genetic architecture of culm strength and yield traits in rice

Lodging resistance in rice is a complex trait determined by culm morphological and culm physical strength traits, and these traits are a major determinant of yield. We made a detailed analysis of various component traits with the aim of deriving optimized parameters for measuring culm strength. Genotyping by sequencing (GBS)-based genome-wide association study (GWAS) was employed among 181 genotypes for dissecting the genetic control of culm strength traits. The VanRaden kinship algorithm using 6,822 filtered single-nucleotide polymorphisms (SNPs) revealed the presence of two sub-groups within the association panel with kinship values concentrated at<0.5 level, indicating greater diversity among the genotypes. A wide range of phenotypic variation and high heritability for culm strength and yield traits were observed over two seasons, as reflected in best linear unbiased prediction (BLUP) estimates. The multi-locus model for GWAS resulted in the identification of 15 highly significant associations (p< 0.0001) for culm strength traits. Two novel major effect marker–trait associations (MTAs) for section modulus and bending stress were identified on chromosomes 2 and 12 with a phenotypic variance of 21.87% and 10.14%, respectively. Other MTAs were also noted in the vicinity of previously reported putative candidate genes for lodging resistance, providing an opportunity for further research on the biochemical basis of culm strength. The quantitative trait locus (QTL) hotspot identified on chromosome 12 with the synergistic association for culm strength trait (section modulus, bending stress, and internode breaking weight) and grain number can be considered a novel genomic region that can serve a dual purpose of enhancing culm strength and grain yield. Elite donors in the indica background with beneficial alleles of the identified major QTLs could be a valuable resource with greater significance in practical plant breeding programs focusing on improving lodging resistance in rice.

1 Introduction

Lodging is a major constraint in rice, affecting both grain yield and quality, particularly in coastal areas prone to cyclonic weather. Lodged plants exhibit limited photosynthetic ability and dry matter production, as the stooping canopy causes losses in crop yield (Islam et al., 2007; Kashiwagi et al., 2010; Mariani and Ferrante, 2017). In severe cases, lodging can result in more than 50% yield loss or complete loss of the crop (Setter et al., 1997; Lang et al., 2012; Luo et al., 2022). It can also lead to pre-harvest sprouting of the panicles due to the high moisture content of the lodged plants. Lodging also poses difficulties in harvest operations, increases the cost of grain drying, and leads to rice mycotoxin contamination (Kashiwagi et al., 2005; Nakajima et al., 2008). The rice crop is prone to lodging from the mid-grain filling stage to the grain ripening, hardening, and maturity stages. The continuous translocation of organic matter to the spikelets from the culm sheath weakens the mechanical strength of the culm with an increase in mass in the top portions of the plant, i.e., increasing grain size grain filling, and grain weight, increasing the risk of lodging or making the culm vulnerable to lodging (Ishimaru et al., 2008).

A breakthrough “Green Revolution” in rice was successfully achieved with the introduction of high-fertilizer-responsive semi-dwarf varieties. In addition to doubling rice grain yield, the short stature conferred by the semidwarf1 (sd1) gene increased lodging resistance, making rice capable of supporting heavier panicles. Breeders were convinced of the suitability of semi-dwarf varieties for increased crop productivity, and many improved semi-dwarf varieties have been developed and cultivated worldwide (Peng et al., 2000; Ashikari et al., 2002; Asano et al., 2007; Hirano et al., 2017). The Sd1 gene, which encodes the gibberellin synthesis enzyme, and the loss-of-function mutant sd-1 allele reduce plant height and enhance lodging resistance (Sasaki et al., 2002). Thus, the semi-dwarf plant type has been for decades the main target in improving the lodging resistance and harvest index of rice (Keller et al., 1999; Peng and Khush, 2003). However, lodging remains a serious problem in widely adapted, high-yielding elite cultivars, despite the short stature conferred by thesd1 gene. The susceptibility of rice plants to lodging varies even among cultivars with short plant height, but with more biomass (Terashima et al., 1992). Reduced culm strength in semi-dwarf varieties due to decreased culm diameter, thickness, and plant biomass has been reported (Okuno et al., 2014; Ookawa et al., 2016) as sd1 is associated with short plant height; it also has negative pleiotropic effects on culm morphology in rice (Yano et al., 2015). Reducing the plant’s height also reduces its photosynthetic capacity and leads to a decrease in total biomass production, thus restricting the plant’s potential for further yield increases (Flintham et al., 1997; Murai et al., 2002; Okuno et al., 2014). These make it enigmatic to improve lodging resistance and, consequently, the harvest index by using the semi-dwarf trait alone. Further, attempts to increase the sink capacity, i.e., grain yield, will be successful only when the source efficiency is enhanced, which is the basis of the new plant type concept (Jyothi et al., 2018). Although an improved source–sink relationship is essential to maximize the harvest index, it also promotes biomass via increased stem and leaf elongation with an overall increase in plant height. Tall varieties are prone to lodging due to the breaking of the basal culm, while the plant type concept advocates semi-tall stature to achieve higher yields, thus making culm strength imperative to prevent discernible negative effects of lodging in high-yielding varieties.

Lodging resistance is known to be a complex trait influenced by many interacting agro-morphological, biochemical, and anatomical traits. Different techniques have been used to measure culm strength and lodging resistance in rice. The bending moment at breaking, which is a function of section modulus and bending stress, was used as a measure for the physical strength of the culm (Ookawa et al., 2010). Section modulus is directly influenced by the culm morphology (diameter and wall thickness), and bending stress is a function of culm cell wall components, such as cellulose and lignin content. Pushing resistance measured with a prostrate tester is the most common parameter used in assessing the culm strength (Terashima et al., 1992; Won et al., 1998; Berry et al., 2003; Kashiwagi and Ishimaru, 2004). However, because the tiller number influences the pushing resistance of the lower part of the culm, bending stress was measured considering both the tiller number and prostrate tester readings ranging from 0 to 40 (Hai et al., 2005; Yadav et al., 2017; Jyothi et al., 2018). Resistance to prostrating on a visually rated scale of 0 (prostrating) to 1 (no or little prostrating) was also used for measuring lodging resistance to typhoons (Kashiwagi and Ishimaru, 2004; Zhao et al., 2023). Recently, gravity center height, along with the breaking strength of the basal internode and culm morphology traits, was used in the measurement of lodging resistance traits (Yang et al., 2023). It is now evident that the culm strength parameter is described in different ways by researchers, and efforts are leading toward the development of an optimized parameter that truly reflects the strength of the culm.

To understand the physiological and morphological basis of culm strength and improve it further, prior knowledge of the genetic control of mechanisms that regulate culm strength is a prerequisite. Identification of genes can hasten precision breeding aimed at improving lodging resistance. A number of quantitative trait locus (QTL) mapping studies for culm length, culm strength, and culm thickness related to lodging resistance have been carried out using different segregating populations in rice (Kashiwagi and Ishimaru, 2004; Mu et al., 2004; Kashiwagi et al., 2008; Zhu et al., 2008; Ookawa et al., 2010; Yano et al., 2015; Yadav et al., 2017; Nomura et al., 2019; Yang et al., 2023; Zhao et al., 2023). To elucidate the molecular mechanisms of complex traits like lodging resistance and associated culm-related traits, a genome-wide association study (GWAS) is a promising strategy to exploit the abundant genetic diversity of rice captured among the diverse rice germplasm. GWAS uses numerous historic recombinations in a large natural population and offers the potential to localize the trait genetic determinants effectively to a narrower region (Bollinedi et al., 2020). Further, the application of single-nucleotide polymorphism (SNP) in a GWAS study provides dense coverage of markers in the entire genome, which helps to identify the functional variation governing the trait in a more precise manner. The results of GWAS can be used either for dissecting mechanisms of lodging resistance or for improving the prediction accuracy of genome-wide predictions in genomic selection programs (Anilkumar et al., 2023).

In this milieu, we hypothesized that the diverse germplasm panel tested over two seasons would have the potential to reveal the candidate genes that regulate culm strength and lodging resistance. We tested the hypothesis by adapting the GWAS approach with a multi-locus model for detecting genomic regions associated with the target traits. The present investigation aimed to identify (1) the genomic regions associated with culm strength, (2) the genetic associations between culm strength and yield traits, (3) the putative candidate genes within the identified genetic regions, and (4) a useful genetic resource for concurrently boosting culm strength and lodging resistance and provide a basis for further research on lodging resistance mechanisms in rice.

2 Materials and methods

2.1 Plant material and field experiments

The study consisted of an association mapping panel comprising 181 diverse genotypes, which included tropical japonica accessions (TrJ: 30), indica landraces (Ind-L: 7), indica cultivars (Ind-C: 26), indica/indica-derived lines (Ind-D: 55), and indica/tropical japonica-derived lines (ITrJ-D: 63). TrJ and Ind-L are from a global collection maintained at ICAR-IIRR, Hyderabad, India; Ind-C are released cultivars collected from all across the country; Ind-D are multiparent-derived, marker-assisted forward breeding introgression lines in diverse elite backgrounds; and ITrJ-D are from the crosses of multiparent elite indica/tropical japonica genotypes developed at ICAR-IIRR (Supplementary Table S1). The genetic purity of experimental material was maintained by continual self-pollination over the years. The complete set of genotypes was evaluated in two cropping seasons at ICAR-IIRR, Rajendranagar farm (17°19′N; 78°23′E, 542 m) during the dry seasons (DS) of 2022 and 2023. An augmented randomized complete block design (ARCBD) was adopted for field evaluation in both years. For the experimental design, 10 genotypes, including one TRJ accession (IRGC 39111), five released cultivars (i.e., DRR Dhan 54, DRR Dhan 48, Improved Samba Mahsuri, Samba Mahsuri, Swarna), and four Y-BILs (i.e., RMS 2077, RMS 2097, RMS 2495, and RMS 2509) were used as checks. The experiment was laid out in four blocks, replicating the checks block-wise and randomizing the test treatments across all blocks. The seeds were germinated on a raised seedbed to ensure uniform germination and 28-day-old seedlings were transplanted into a puddled rice field. In DS 2022, the experimental material was sown on 6 January 2022, and transplanted on 5 February 2022. In DS 2023, the experimental material was sown and transplanted on 26 December 2022, and 22 January 2023, respectively. Each genotype was planted in two rows at 14 hills per row with a 20-cm spacing between rows and a 15-cm spacing between plants. The recommended package of practices was adopted to ensure a good crop stand in both seasons. The lines were subjected to phenotypic evaluation of culm strength and yield-associated traits and genomic studies for the discovery of marker–trait associations and identification of putative candidate genes (schematically illustrated in Figure 1).

Figure 1

Figure 1 Schematic flow of the experimental strategy, including phenotypic evaluation, trait grouping for culm strength and yield traits, genome-wide studies for marker–trait associations, and identification of putative candidate genes. (A) Measurement of pushing resistance using a prostrate tester in the field. TR, tester reading with prostrate testing as a measure of pushing resistance. (B) Measurement of internode breaking weight in increments of 25 g on culm samples for basal internodes. (C) Measurement of culm morphology traits using digital Vernier calipers. ODMa, outer diameter in mm on the major axis; ODMi, outer diameter in mm on the minor axis; IDMa, inner diameter in mm on the major axis; IDMi, inner diameter in mm on the minor axis; CD, culm diameter in mm; CT, culm wall thickness in mm; PH, plant height in cm; IL, basal internode length in cm. (D) Color key and histogram depicting population structure. (E) A representative Manhattan plot for genomic regions associated with culm strength traits. (F) A representative candidate gene within the identified QTL region. GBS, genotyping by sequencing; MTAs, marker–trait associations. (G) A representative plant with a high yield (high grain number) and strong culm. PN, panicle number; GN, grain number; PW, panicle weight in g; GY, grain yield in g; CS, culm strength; QTL, quantitative trait locus.

2.2 Evaluation of morphological traits

A total of 12 quantitative traits were evaluated, viz., plant height, length of the basal internode, culm thickness, culm diameter, pushing resistance as a measure of prostrate tester reading, tiller number, bending stress, section modulus, breaking resistance as a measure of internode breaking weight, panicle number, grain number, and panicle weight. For culm parameters, two plants were uprooted 20 days after heading, and data were collected using destructive sampling in the lab as described below. For yield parameters, while data were collected at maturity in the field on plant height, tiller number, and panicle number, two panicles were collected randomly at harvest for post-harvest data. Plants on the edges were avoided, and plants from the middle of the rows were chosen for data collection to avoid the border effect. Plant height (in cm) was measured as the length between the plant base and the panicle tip, and tiller number was measured as the total number of tillers, panicle number as the number of panicle-bearing tillers, grain number as the total number of filled grains in a panicle, and panicle weight as the weight of a fully matured panicle.

2.2.1 Measurement of culm morphology traits

Destructive sampling was used for the measurement of culm traits 20 days after heading. Two plants were uprooted from the main plot, and the lower internode (second basal internode) of the main culms was sampled for recording the observations. The second basal internodes were dissected transversely at the mid-point, and the culm diameter was measured using digital Vernier calipers.

Culm thickness (in mm) was measured as the average difference between the outer and inner diameters on the major and minor axes.

$$\text{Culm thickness }(\text{CT})=\left[\frac{(a_1+b_1)}{2}\right]-\left[\frac{(a_2+b_2)}{2}\right].$$

Culm diameter (in mm) was measured as the average of the culm outer diameter on the minor and major axes using digital Vernier calipers (Figure 1C).

$$\text{Culm diameter }(\text{CD})= \frac{(a_1+b_1)}{2}.$$

Section modulus was determined as described by Ookawa et al. (2010); Ookawa and Ishihara (1992) using the following formula:

$$\text{Section modulus}(\text{SM})=\frac{\pi }{32}\times \frac{(a_1^{\mathrm{ 3}}{b}_1-a_2^{\mathrm{ 3}}{b}_2)}{a_1},$$

where a₁ is the outer diameter of the minor axis in an oval cross-section, b₁ is the outer diameter of the major axis in an oval cross-section, a₂ is the inner diameter of the minor axis in an ovalcross-section, and b₂ is the inner diameter of the major axis in an oval cross-section.

2.2.2 Measurement of culm physical strength traits

The prostrate tester (DIK-7400, Daiki Rika Kogyo Co. Ltd., Tokyo, Japan) was used to measure the pushing resistance of the culm (Figure 1A), and bending stress was calculated as described by Hai et al. (2005) using the following formula:

$$\text{Bending stress }(\text{BS})=\left(\frac{TR}{40}\right)\times \left(\frac{1,000}{TN}\right),$$

where TR is the prostrate tester reading value (a measure of pushing resistance) and TN is the tiller number.

Internode breaking weight (IBW) in increments of 25 g was used as a measure of breaking resistance. IBW was analyzed using a supporting structure composed of two burette stands and clumps spaced at the length of each internode analyzed (Figure 1B). The internode was arranged horizontally with the portions of the terminal nodes held on the clamps of each burette stand, with the addition of an empty container (previously weighed) hanging in the central part of the internode. Weights of the determined mass (25 g) were cumulatively added to the container at regular intervals of 1 second until the internode flexed and collapsed.

2.3 Genotyping by sequencing and marker data generation

Fresh leaf tissue samples from four to six plants per line were collected 20 days after transplanting. The next-generation sequencing (NGS) platform Illumina NovaSeq was used for genotyping by sequencing (GBS). Raw sequence data were processed to obtain high-quality clean reads using Trimmomatic, v. 0.38. The reads of the samples were aligned to the Oryza sativa indica reference genome, https://plants.ensembl.org/Oryza_indica/Info/Annotation/#assembly using BWA MEM (v. 0.7.17), with minimum seed length set to 32 and shorter split hits marked as secondary (parameters: -k 32 –M). The mpileup utility of Samtools (v. 0.1.18) was used to identify SNPs from the sorted BAM file of the sample. The SNPs were filtered based on a minimum read depth of 5 and a quality threshold of 25. Markers without clear physical position information on the chromosomes were removed. The genotyping data were filtered by removing markers with missing values >20% and minor allele frequency (MAF)<5%.

2.4 Genome-wide association analysis

The average measurements of all the traits recorded in each season were subjected to the estimation of best linear unbiased predictions (BLUP) that were used in further analysis. In the case of checks, the data of all the replications were averaged and used in BLUP estimates. Average values in each year were used to estimate, across years, the BLUP values for each trait. The across-year BLUP values for all the traits measured on experimental genotypes and SNP marker information were used as inputs for GWAS analysis. Genome-wide association analysis was performed following a multi-locus analysis method using the genome association and prediction integrated tool (GAPIT) v.3 software package in R (Lipka et al., 2012). The relatedness between individuals in the population was delineated by estimating a kinship similarity matrix (K) and population structure (Q) modeled using principal component analysis (PCA) as random and fixed effects, respectively (Yu et al., 2006), making the model more stringent and able to control false-positive marker–trait associations (MTAs). The multi-locus-based GWAS model of BLINK (Huang et al., 2019) was applied for MTA analysis. Population structure and the kinship similarity matrix were also incorporated according to the model type. The BLINK model is associated with an efficient mixed model approach (EMMA), which corrects for population structure and simultaneously identifies the marker–trait associations. The exploratory threshold with false discovery rate-corrected p ≤ 0.001 (−log₁₀p ≥ 3) was used to report significant marker–trait associations.

2.5 Identification of putative candidate genes

After significant MTAs were identified using data from multiple environments, the flanking sequence spanning 150 kb upstream and downstream of the significant SNP position was used as a query against the Nipponbare rice reference genome, IRSGP1.0 RefSeq, v.1.1. Subsequently, JBrowse was used to examine candidate genes residing in chromosomal regions harboring significant MTA. The genes within the search limits were filtered based on the reported functions, and only those genes that are known to have an impact on plant development were shortlisted as putative candidate genes that regulate culm strength.

3 Results

3.1 Variation for culm strength- and yield-related traits in the association panel

Broad-sense heritability was high (79 to 99) for all the studied traits except for culm thickness, which had moderate heritability (59) (Figure 2A). BLUP values were calculated to further remove the environmental effects from trait data, and highly significant differences (p< 0.001) were found among the genotypes in the association panel for all the traits (Table 1). For example, it varied from 61.62 cm to 142.5 cm for plant height, 6.48 cm to 20.87 cm for lower internode length, 2.68 mm to 7.29 mm for culm diameter, 0.01 mm to 2.21 mm for culm thickness, 5.11 to 29.63 for pushing resistance as a measure of prostrate tester reading, 3.22 to 31.72 for tiller number, 116 g to 1,506 g for breaking resistance as a measure of internode breaking weight, 2.38 to 31.43 for panicle number, 79.82 to 410 for grain number, and 0.97 to 5.51 for panicle weight. Third-degree statistics-skewness and fourth-degree statistics-kurtosis were employed to understand the distribution of phenotypes in the population. The skewness of the population for all the traits was positively significant except for lower internode length, culm diameter, and panicle weight. However, kurtosis for all the traits except lower internode length and panicle weight was greater than three, indicating a leptokurtic distribution of phenotypes in the population (Table 1). The mean phenotypic data on culm strength and yield traits among the 181 genotypes in the association panel are provided in Supplementary Table S2.

Figure 2

Figure 2 Genetic and phenotypic parameters for the components of culm strength and yield traits in the association panel. (A) Broad-sense heritability for seven culm strength (PH, plant height; IL, lower internode length; CD, culm diameter; CT, culm wall thickness; IBW, internode breaking weight [breaking resistance]; TR, prostrate tester reading value [pushing resistance]; TN, tiller number) traits and three yield (PN, panicle number; GN, grain number; PW, panicle weight) traits. (B) Correlation matrix showing correlation coefficients among seven culm strength traits and three yield traits. (C) Histograms showing the phenotypic distribution of section modulus (SM), bending stress (BS), and breaking resistance as a measure of internode breaking weight (IBW). *p < 0.05, **p < 0.01, ***p < 0.001.

Table 1

Table 1 Phenotypic variation and distribution pattern of culm strength and grain yield traits in the association mapping panel used in the study.

A correlation analysis among the seven culm strength traits and three yield traits was performed to understand the influence of variables on lodging resistance. A significant positive correlation was observed between culm morphology traits and breaking resistance (IBW), except for a negative correlation with internode length. Pushing resistance (TR) also correlated positively with culm morphology and breaking resistance, but it was not found to be statistically significant. Further, tiller number was negatively correlated with culm morphology traits and breaking resistance but positively correlated with pushing resistance (Figure 2B). These results suggest that culm morphology traits and breaking resistance have a direct influence on culm strength while pushing resistance is influenced by tiller number. Hence, both pushing resistance and tiller number should be considered in measuring the bending stress of the lower portion of the culm. Although breaking resistance showed a positive correlation with plant height, it showed a negative correlation with lower internode length, while pushing resistance showed a positive correlation with both height-related traits. Considering the indefinite relationship between culm strength and height traits in the present study, height traits were excluded from further analyses. While a positive correlation was observed for culm morphological and physical strength traits with grain number, correlations were negative with panicle number and had a non-significant association with panicle weight (Figure 1B). These findings suggest that the relationship of culm strength traits with grain number is synergistic, with panicle number being a tradeoff.

Based on the correlation analysis, the traits identified as key influencers of lodging resistance were: a) section modulus as a function of culm morphology traits; b) bending stress as a measure of pushing resistance and tiller number; and c) breaking resistance as a measure of internode breaking weight. The phenotypic distribution pattern in the form of histograms and density curves is depicted in Figure 2C. Section modulus varied from 1.64 mm³ to 30.6 mm³ with a mean of 11.18 mm³, bending stress ranged from 6.51 to 59.64 with a mean of 24.84, and internode breaking weight varied from 116 g to 1,506 g with a mean of 435 g. The normality of the trait distributions suggests that they may be controlled by several minor effect loci (Figure 2C). These three traits were further chosen for studying marker–trait associations.

3.2 Genome-wide distribution of SNPs

A total of 3,72,087 sites obtained by GBS from 181 genotypes in the association mapping panel were subjected to further filtration in order to retain the sites with SNPs. A total site count of 35,331 was achieved after all the filtration. After removing the SNPs with MAF< 5%, the final 6,822 SNPs were used in the construction of a linkage map depicting the genome-wide distribution of SNPs on the 12 rice chromosomes (Figure 3A). The SNPs spanned a total genome size of 369.72 Mb, with the highest number of 931 SNPs covering 38.04 Mb on chromosome 2 and the lowest number of 374 SNPs spanning 21.2 Mb on chromosome 9. The maximum genome coverage of 46.86 Mb with 876 SNPs was on chromosome 1, followed by 41.10 Mb with 471 SNPs on chromosome 3 (Figures 3A, B).

Figure 3

Figure 3 Genome-wide distribution of 6,822 SNPs across 12 chromosomes among 181 rice genotypes in the association panel. (A) Linkage map with 6,822 SNPs; the color of the vertical line indicates SNP density in a 1-Mb window. (B) Graphical representation of the chromosome-wide distribution of SNPs. SNPs, single-nucleotide polymorphisms.

3.3 Population structure and kinship analysis

The VanRaden kinship algorithm using the GAPIT was performed on the genome-wide 6,822 SNPs to ascertain the population structure and relatedness among genotypes in the association mapping panel. The relatedness in the association mapping (AM) panel was lower than the kinship index of 0.5, indicating less relatedness among the genotypes (Figure 4A). A lower level of relatedness among individuals in the GWAS panel reduces the false positives and increases the precision of the results. The kinship heatmap revealed the presence of two sub-groups within the AM panel and represented the distribution of genotypes within and between sub-populations (Figure 4A). There was no significant grouping of the genotypes, indicating an even distribution of the alleles in the panel, making it a perfect association mapping panel for GWAS. Principal component analysis provided insights into the SNP-based diversity in the association panel, which also detected the presence of two major population groups indicated by two significant components explaining the maximum variation in the association panel (Figure 4B). The scree plot (Figure 4C) showed that the first principal component (PC) explained the highest variation of 11%, followed by the second PC explaining 7%, and the third PC explaining 3.9% of the total variation. The variation explained by the first three PCs, which account for only 21.9% of the total variation, suggested the diverseness of the panel composition. The two major population groups in the panel were denoted as POP1 and POP2, and sub-populations within them as POP1A, POP1B, POP2A, and POP2B. It is interesting to note that there is no clear classification of indica and japonica in the sub-populations, and the genotypes belonging to both sub-species were interspersed in all the sub-population groups. POP2B was the largest and constituted 63.5% of the panel with 115 accessions, followed by POP1B, POP1A, and POP2A with 17.68, 11.60, and 7.18% of the panel, respectively. Sub-population-wise details of the genotypes of the AM panel in accordance with the kinship matrix are provided in Supplementary Table S1.

Figure 4

Figure 4 Population structure analysis of 181 genotypes in the association panel based on 6,822 SNPs. (A) Heatmap of the kinship matrix. The heatmap shows the level of relatedness among the population. The darker areas show a high level of relatedness between genotypes, and the dendrogram depicts the clustering of sub-populations. (B) 3D representation of principal component (PC) analysis showing no clear-cut grouping. (C) Scree plot depicting the number of significant PCs. There were three PCs that explained a cumulative variation of ∼22%, and nearly 10 PCs covered >95% of the variation in the population. SNPs, single-nucleotide polymorphisms.

3.4 Genome-wide association studies for culm strength traits

Genotyping data of 6,822 SNPs and two-season phenotypic data of three culm strength traits were subjected to association analysis using a multi-locus model of BLINK. The Q–Q plots depicted less deviation of the observed p-values from the expected p-values and were therefore chosen as the best fit (Figure 5). At p ≤ 0.0001 (−log₁₀p ≥ 4), 15 highly significant MTAs across seven out of 12 chromosomes were identified for the three culm strength traits (Table 2). Two MTAs (qSM2.1 and qBS12.1) showed a major effect with a phenotypic variance explained (PVE) of more than 10%, and the remaining 13 MTAs showed a minor effect with a PVE of less than 10%.

Figure 5

Figure 5 Quantile-quantile (Q–Q) plots and Manhattan plots for genomic regions associated with culm strength traits. Q–Q plots showing deviation of observed −log₁₀p values and expected −log₁₀p values, indicating the significant associations of the SNP with the trait. In Manhattan plots, chromosomes are represented on the x-axis explaining chromosome-wise SNP distribution, and significant associations at −log₁₀p values are represented on the y-axis. (A) Section modulus (SM), (B) bending stress (BS), and (C) breaking resistance as a measure of internode breaking weight (IBW). The threshold to select significant MTAs is represented as a black dotted horizontal line at LOD score 3. The putative candidate genes associated with lodging resistance are marked in the vicinity of significant MTAs in the Manhattan plots. SNPs, single-nucleotide polymorphisms; MTAs, marker–trait associations; LOD, logarithm of the odds.

Table 2

Table 2 Significant marker–trait associations for culm strength traits at p ≤ 0.0001 (−log₁₀p ≥ 4) were identified using the multi-locus BLINK model.

On further reducing the threshold to p ≤ 0.001 (−log₁₀p ≥ 3), a total of 198 MTAs were identified, corresponding to 78, 41, and 79 associations for section modulus, bending stress, and internode breaking weight, respectively. Several overlapping association signals with similar effects and sizes were observed. In earlier reports, a locus was defined as a chromosomal region at which adjacent pairs of associated SNPs were less than a certain physical distance without considering varying linkage disequilibrium levels (Yang et al., 2012). We considered adjacent SNPs spanning less than 300 kb as a single locus to reduce the redundancy of association signals similar to the previous reports (Huang et al., 2010; Chen et al., 2014; Xie et al., 2015; Guo et al., 2018; Guo et al., 2021) for identifying overlapping MTAs. After eliminating the redundant signals, 126 associations between 114 lead SNPs were identified (Supplementary Table S3), corresponding to 97 loci and 50, 31, and 45 MTAs for section modulus, bending stress, and internode breaking weight, respectively. MTAs were found on all the chromosomes, and chromosome 2 had the highest number of MTAs (23) followed by chromosome 5 (20 MTAs), chromosome 1 (15 MTAs), chromosome 12 (14 MTAs), chromosome 7 (12 MTAs), chromosomes 3 and 10 (11 MTAs each), chromosome 4 (10 MTAs), chromosomes 8 and 11 (7 MTAs each), and chromosome 6 (6 MTAs) (Supplementary Table S3). Manhattan plots depicting the significant association signals at the threshold p ≤ 0.001 (−log₁₀p ≥ 3) are provided in Figure 5.

A total of 50 MTAs were identified for section modulus, of which five were highly significant at p ≤ 0.0001 (−log₁₀p ≥ 4) (Figure 5A; Table 2), and the remaining 45 were significant at p ≤ 0.001 (−log₁₀p ≥ 3) (Supplementary Table S3). A single major effect QTL, qSM2.1 on chromosome 2, explained a maximum PVE of 21.87% significant at P_BLINK = 4.86 × 10⁻⁸ and an additive effect of 2.39 mm³ on the trait (Table 2). The highest additive effect of 4.23 mm³ for section modulus was identified with an association at snpM4010 (P_BLINK = 2.9 × 10⁻³) on chromosome 7 (Supplementary Table S3).

For bending stress, four MTAs highly significant at p ≤ 0.0001 (−log₁₀p ≥ 4) (Figure 5B and Table 2) and 27 MTAs significant at p ≤ 0.001 (−log₁₀p ≥ 3) were identified (Supplementary Table S3). qBS12.1 is a major effect of QTL at snpM6498 on chromosome 12, significant at P_BLINK = 1.32 × 10⁻⁸, with a PVE value of 10.14% and an additive effect of 4.61. It also colocalized the section modulus, which had an additive effect of 1.24 mm³ (Supplementary Table S3). Within the same locus, an association of internode breaking weight was found in close proximity of 126 bp at snpM6497 with an additive effect of 67.36 g and a PVE value of 4.02% (Table 2) (Supplementary Table S3). Further, the highest additive effect of 7.4 was identified for bending stress at snpM401 on chromosome 1 (Supplementary Table S3).

A total of 50 MTAs were identified for internode breaking weight, of which four were highly significant at p ≤ 0.0001 (−log₁₀p ≥ 4) (Figure 5C; Table 2), and the remaining 41 were significant at p ≤ 0.001 (−log₁₀p ≥ 3) (Supplementary Table S2). For internode breaking weight, the significant MTA, qIBW2.3 (P_BLINK = 3.29 × 10⁻⁵) on chromosome 2 at snpM1659, had a PVE value of 8.63% and an additive effect of 91.99 g (Table 2). It also colocalized the section modulus with an additive effect of 1.24 mm³ (Supplementary Table S3). The highest additive effect of 208.41 g was identified at snpM3981 on chromosome 7 for internode breaking weight, with a PVE value of 5.64% (Supplementary Table S3).

Tests for epistatic interactions were performed using a linear regression method implemented with PLINK v1.07. Manhattan plots depicting significantly associated genome-wide epistatic (SNP–SNP) interactions for CS and GY are presented in Figure 6. Epistatic tests identified three and one SNP–SNP interactions associated with SM and IBW, respectively (Supplementary Table S4). For SM, three significant SNPs (M715, M1454, and M6514) on chromosomes 1, 2, and 12 interacted with three SNPs (M3955, M5372, and M5616) on other chromosomes 7, 10, and 10, respectively. For IBW, the SNP (M5613) on chromosome 10 exhibited significant interactions with another SNP (M5485) on the same chromosome. However, there were no epistatic interactions among the 15 highly significant MTAs presented in Table 3. The epistatic interactions noted with the reduction of threshold value p ≤ 0.001 for MTAs of IBW and SM are likely due to numerous contributing loci for the complex culm strength trait, generally with small effect.

Figure 6

Figure 6 Manhattan plots for significantly associated genome-wide epistatic (SNP–SNP) interactions for section modulus (SM), bending stress (BS), breaking resistance as a measure of internode breaking weight (IBW), grain number (GN), and panicle weight (PW). SNP, single-nucleotide polymorphism.

Table 3

Table 3 Colocalization of MTAs with previously reported putative candidate genes related to lodging resistance.

3.5 Colocalization of MTA for culm strength and yield traits

Based on the significant correlations observed between culm strength traits and yield traits in the present investigation (Figure 2B), we explored the genetic bases of loci that synergistically enhance culm strength and yield. Culm strength traits correlated positively with grain number and panicle weight and negatively with panicle number (Figure 2B); thus, we considered grain number and panicle weight for GWAS along with culm strength traits. In total, we identified 16 grain yield loci associated with culm strength traits, including seven loci for grain number and nine loci for panicle weight (Table 4). Each locus harbored one-grain yield trait and one to three culm strength traits on seven out of 12 chromosomes. Out of 16 associations between grain yield and culm strength traits, chromosome 2 had three clusters; chromosomes 1, 4, 5, 7, 11, and 12 had two clusters each; chromosome 8 had a single cluster.

Table 4

Table 4 Genetic loci linking culm strength and yield traits at p ≤ 0.001 (−log₁₀p ≥ 3).

A QTL hotspot on chromosome 12 (snpM6498) colocalized grain number and all three culm strength traits. It could enhance grain yield and culm strength with an additive effect of 16.83, 1.24 mm³, 4.61, and 67.36 g for grain number, section modulus, bending stress, and internode breaking weight, respectively. Another six associations between grain number and culm strength traits showed a synergistic relationship on the same chromosome. The clusters for grain number and section modulus were identified on chromosomes 1, 2, 7, and 11. MTA clusters for grain number and bending stress were identified on chromosomes 4 and 12.

Among the nine genomic regions colocalizing panicle weight and culm strength traits, eight showed a synergistic relationship. An association between panicle weight and internode breaking weight at snpM839 on chromosome 1 was negatively correlated. Positive associations between panicle weight and internode-breaking weight were identified on chromosomes 2, 4, 11, and 12. Panicle weight colocalized with both bending stress and internode-breaking weight on chromosome 2. The two genomic regions that colocalized panicle weight and two culm strength traits, section modulus, and internode breaking weight, were identified on chromosome 5. The grain number locus associated with section modulus and bending stress was identified on chromosome 7.

3.6 Identification of candidate genes

The genomic regions within the same locus of the identified MTAs for culm strength were analyzed in silico for the presence of putative candidate genes previously reported for lodging resistance. A total of 481 putative genes were found, varying from three to 54 per QTL region. Many genes identified within the QTL region were responsible for the regulation of growth and development, resistance to abiotic and biotic stresses, amino acid transport, sugar transport, phytohormone synthesis, hypothetical and expressed proteins, glycoproteins, heat shock proteins, transcriptional factors, and precursors for various biochemical and metabolic pathways. It was interesting to note that six of the 15 MTAs identified in the present investigation were found to be in the vicinity of 11 previously reported putative candidate genes for lodging resistance (Table 3).

4 Discussion

The use of genes associated with lodging resistance in addition to the dwarfing gene is a promising approach to improving lodging resistance and further increasing productivity in rice (Ookawa et al., 2010). Owing to the complex nature of the lodging resistance in rice, it is relatively difficult to map multiple genes using a biparental-derived mapping population, wherein the genes have low individual effects and are sparsely distributed in the gene pool. Gene mapping in biparental populations is limited to the variation that exists only between the two parents used for the development of the mapping population. GWAS, therefore, offers the dual advantage of analyzing the extensive trait variation among the germplasm lines and identifying several genomic regions and alleles affecting the trait (Bollinedi et al., 2020). Globally, there are only a few GWAS reports for culm strength traits (Chigira et al., 2020; Meng et al., 2021; Nomura et al., 2021; Rashid et al., 2022). Furthermore, such studies are limited to only one report (Guo et al., 2021) for genetic associations between culm strength (CS) and grain yield (GY) traits. We report for the first time in the Indian sub-continent on genome-wide associations for culm strength traits and synergistic genetic associations between culm strength and grain yield traits.

4.1 Population structure

It is crucial to have control over population structure in GWAS to prevent spurious marker–trait associations (Nayak et al., 2022). The origin, selection pressure, and reproductive behavior of genotypes all have an impact on familial relatedness among individuals in an association panel (Atwell et al., 2010). In the present study, SNP-based diversity analysis provided insights into genetic relatedness among the individuals in the association panel. There was no clear grouping of indica and japonica in the population structure based on GBS data, and the two sub-populations were identified with admixtures of all the constituents of TrJ, Ind-L, Ind-C, Ind-D, and ITrJ-D (Supplementary Table S1). Similar to our findings, Nomura et al. (2021) also reported no clear and strong clusters in the GWAS panel for strong culm traits. The panel with a major focus on phenotypic variation was composed of culm strength traits. Population admixtures indicated no role for genetic grouping in the variation in culm strength. It is interesting to note the higher mean values for SM, BS, and IBW in ITrJ-D compared to TrJ, Ind-L, Ind-C, and Ind-D (Figure 7), which is similar to the findings of Nomura et al. (2019). Sub-populations arise due to allelic sharing between sub-populations, which is attributed to allelic accumulation due to spontaneous mutation over time (Agrama et al., 2007). PCA confirmed the presence of two sub-populations with a sparse distribution of alleles in the association panel. The kinship matrix generated by the VanRaden algorithm was plotted as a heatmap, showing relatedness values between 0 and +0.5, which indicated poor relationships existing between individuals in the association panel. These results assisted in understanding the population structure of the panel before proceeding to GWAS for the identification of putative genomic regions for culm strength traits. Based on the information about population structure, the multi-locus model BLINK with the EMMA approach has been selected for association analysis, which detects marker–trait associations while simultaneously addressing population structure to reduce the chances of false positives (Zhang et al., 2014; Wang et al., 2016).

Figure 7

Figure 7 Mean phenotypic performance of the culm strength traits in the constituent groups of the GWAS panel. TrJ, tropical japonica accessions; Ind-L, indica landraces; Ind-C, indica cultivars; Ind-D, indica/indica-derived lines; ITrJ: indica/tropical japonica-derived lines. GWAS, genome-wide association study.

4.2 Culm strength phenotype

High broad-sense heritability and a wide range of phenotypic variation in culm strength traits in the present study merit dissection of their genetic basis. The present investigation optimized the parameters for the measurement of culm strength in rice. Regarding the relevance of the traits for lodging resistance and considering the indefinite correlations of the lower internode length with culm-related traits for the set of genotypes in the association panel, we excluded the height traits from GWAS. With reference to the prostrate tester reading value as a measure of pushing resistance in the lower portion of the culm, it is evident that pushing resistance is more highly influenced by tiller number than the actual physical strength of the culm. For example, a strong culm genotype, IRGC 73031 (G123), identified in the present study has a pushing resistance value as low as 7, and its tiller number is only 4; however, it has the highest section modulus value (30.60 mm³) and highest culm diameter (7.23 mm). Hence, bending stress was measured, including both the tester reading value and tiller number, as described in Hai et al. (2005). A similar methodology was adopted in previous studies on the measurement of culm strength (Yadav et al., 2017; Jyothi et al., 2018). To further aid in assessing the culm’s physical strength, we also considered the internode breaking weight to assess the breaking resistance, an additional parameter for the culm’s physical strength. Breaking resistance (BR) as a measure of internode breaking weight is found to be a realistic parameter that depicts the true strength of the basal culm. There was a highly significant positive correlation of internode-breaking weight with other culm-related traits (Figure 2B), and IRGC 73031 (G123) was identified with the highest internode-breaking weight of 1,506 g (Supplementary Table S2).

4.3 Marker–trait associations and putative candidate genes

Epistasis is a genetic phenomenon of interaction that may enhance or reduce the expression (depending on degree and direction) of interacting loci underlying QTLs associated with the complex trait (Yadav et al., 2019). The role of epistasis phenomena in plant breeding has been discussed for a long time, but the extent of the expression level of quantitative traits due to this complex genetic phenomenon is poorly understood. For complex agronomic traits, estimating epistatic loci could help to clarify the complex genetic effects of GWAS loci and elucidate other types of interactions, such as genotype × environment effects (Assefa et al., 2019). Thus, there is a need to identify such epistatic gene interactions that would help to better understand the genotype–phenotype relationship of complex traits such as culm strength. The two major QTLs, qSM2.1 and qBS12.1, identified in the present study with no significant epistatic interactions suggest their direct use in MAB programs.

In the present investigation, the linkage disequilibrium (LD) ranged between 200 kb and 1 Mb in the data-driven analysis for the identified MTAs. For convenience of presentation and summary of GWAS results, a locus is usually defined as a chromosomal region at which adjacent pairs of associated SNPs are less than a certain physical distance without considering varying LD levels (Yang et al., 2012). To maintain the stringency of putative candidate gene identification in the present study, the threshold was fixed at 150 kb on either side of the identified MTAs, similar to Guo et al. (2021). The 167-kb and 123-kb LD decay rates in japonica and indica subpopulations have been reported in previous studies (Huang et al., 2010; Huang et al., 2012), and Rashid et al. (2022) followed a threshold criterion of at least three SNPs in a 170-kb LD block to define the QTLs and enhance the power of GWAS resolution. Notably, six MTAs were found in the vicinity of 11 previously reported putative candidate genes for lodging resistance in the present investigation. Among them, a putative candidate gene similar to the hydroxyproline-rich glycoprotein family protein was identified to be colocalized with qSM1.1 in close proximity to 54.79 kb. These proteins are important structural components of plant cell walls, as they can form a continuous glyco-network with non-cellulosic polysaccharides via covalent bonds or non-covalent interactions, thus strongly contributing to cell wall architecture (Hijazi et al., 2014). In the genomic region of qSM1.2, cinnamoyl-CoA reductase (OsCCR6) was identified at a distance of 117.45 kb. OsCCR6 is a key enzyme in lignin biosynthesis (Zhao et al., 2021). Lignin is the main structural component of vascular plants’ secondary cell wall, which is related not only to plant growth and development but also to mechanical strength. Further, qSM1.2 was found at a distance of 224 kb from the putative candidate gene regulating the synthesis of xylan O-acetyltransferase 1 (OsXOAT1). Acetylated xylans are the principal hemicelluloses in the cell walls of grass species (Qaseem and Wu, 2020). In the close vicinity of qIBW2.1, genes responsible for modulation of the gibberellin signaling pathway at a distance of 44.83 kb and a gene for auxin-responsive protein at 5.46 kb were identified, which have an effect on lodging resistance in rice. Auxin prompts the expression of genes related to lignin biosynthetic peroxidase (Prx) in Zinnia elegans and secondary growth/lignification (Gutiérrez et al., 2009). The effect of the plant hormone GA on lodging resistance in rice and increased total biomass signifies the positive impact of overexpression of GA on lodging resistance due to increased culm diameter and lignin deposition (Okuno et al., 2014). Two putative candidate genes associated with lodging resistance were identified in the qIBW5.1 locus. Genes regulating gibberellin and glucose transport (OsSWEET3a) at 17.92 kb and a gene for arabinogalactan protein 23 (OsAGP23) (Hijazi et al., 2014) at 116.4 kb were identified in the vicinity of qIBW5.1. In the locus harboring qIBW8.1, cinnamoyl-CoA reductase (OsCCR29) at 223 kb was identified. The MTA qBS10.1 was in close proximity (115 kb) to the known candidate gene cellulose synthase A7 (OsCesA7) on chromosome 10. OsCesA7 is mainly responsible for cellulose synthase (UDP-forming) activity, regulating the cellulose biosynthetic process and secondary cell wall formation, thus playing a major role in enhancing stem strength (Tanaka et al., 2003; Wang et al., 2010; Kotake et al., 2011; Noda et al., 2015; Fan et al., 2017). The same locus also harbored a gene for wall-associated kinases (OsWAK100) at a distance of 338.77 kb, which is known to regulate secondary cell wall thickening (Ma et al., 2022). The putative candidate genes identified in the present study further need expression-based validation and functional characterization. Such characterized candidate genes may be beneficial to enhance culm strength, leading to lodging resistance in rice.

Deciphering the relationship between lodging resistance and yield traits is crucial to increasing crop productivity. In crop breeding, relationships among traits should be considered to enhance desirable correlated traits and simultaneously reduce undesirable trade-offs (Chen and Lubberstedt, 2010). To our knowledge, few studies have focused on genetic associations between culm strength and yield traits (Ookawa et al., 2010; Wang et al., 2016; Guo et al., 2021). Among the genetic associations between culm strength and grain yield traits identified in the present study, except for one negative association between internode breaking weight and panicle weight on chromosome 1, the remaining 15 associations between culm strength and panicle weight or grain number were synergistic. Some of them were colocalized with the same SNP, and some of them were within the close vicinity of 46 bp to 169.28 kb. The same QTL hotspot on chromosome 12 harboring all the culm strength traits (section modulus, bending stress, and internode breaking weight) at snpM6498 (pos 8010763) also harbors grain number. Also, it is a novel MTA. Internode breaking weight and grain number colocalized on chromosome 11 at a close proximity of 46 bp. Another colocalization of culm strength traits, section modulus, and internode breaking weight with panicle weight was identified on chromosome 5 at a proximal distance of 0.94 kb. Previously, significant phenotypic correlations between culm strength and yield traits were observed, and 63 loci linking them were detected (Guo et al., 2021). Interestingly, we also identified a putative candidate gene, rice floricaula/aberrant panicle organization 2 (APO2/RFL), linked to the MTA colocalizing BS and GN on chromosome 4 at a close proximity of 373.47 kb. APO2/RFL is a probable transcription factor controlling inflorescence, flower development, and short and solid culm, thus enhancing both grain yield and lodging resistance (Wang et al., 2016). QTL SCM2, a weaker allele of APO1, increases spikelet number without reducing panicle number (Ookawa et al., 2010). These colocalized MTAs can be targeted for the simultaneous improvement of culm strength and grain yield traits in indica rice, and favorable alleles of the major effect loci associated with only culm strength can be introgressed in the elite high-yielding indica cultivar background.

4.4 Elite donors

Culm morphological and physical strength differences in indica and temperate japonica were reported earlier; the section modulus of indica was greater than that of japonica, and the bending stress of indica was inferior to that of japonica cultivars (Ookawa et al., 2010). Superior alleles of STRONG CULM 1 and 2 (SCM1 and SCM2) QTLs, which increase culm strength, were detected in an indica variety Habataki (Ookawa et al., 2010); similarly, SCM3 and SCM4 QTLs were detected in a tropical japonica variety Chugoku 117 (Yano et al., 2015; Yang et al., 2023). Although we observed genetic and phenotypic variations for culm strength traits in the indica, tropical japonica, and breeding lines, we also identified genotypes with superior trait values for both culm morphology and physical strength traits in both indica and tropical japonica, including their derived lines. The present study identified three genotypes, IRGC 73031 (G123), JBB 4547 (G61), and JBB 4514 (G14), with high values for all the culm strength traits (section modulus, bending stress, and internode breaking weight) and panicle weight (Supplementary Table S2). Additionally, 11 genotypes were identified with high values for one to four traits. Among the genotypes with superior value for all the culm strength traits, G123 is a tropical japonica accession, while G14 is an indica line, suggesting the presence of superior alleles for wider culms in japonica and the presence of superior alleles for physical strength in indica genotypes. G61 is a derived line possibly combining the superior trait values for culm morphology and physical strength from indica and japonica. Nomura et al. (2019) reported very high culm strength in Takanari, a derived line from the cross of indica and temperate japonica. Further, out of the aforementioned 14 strong culm genotypes, JBB 4514 (G14), JBB 6436 (G165), and JBB 1120 (G20) were identified with a high panicle weight of more than 5 g and G36 with a high grain number (Supplementary Table S2). The novel donors identified in the present study with beneficial alleles for enhanced culm strength and high yield can be utilized in the simultaneous improvement of yield and lodging resistance in rice.

5 Conclusion

We report for the first time GWAS for lodging resistance traits and synergistic associations for culm strength and grain yield from the Indian subcontinent. The two novel QTLs, qSM2.1 for section modulus and qBS12.1 for bending stress, with a major effect identified in the present study, can be further validated and employed in marker-assisted breeding for improvement of lodging resistance in elite high-yielding lodging-prone cultivars. We propose that the validation and introgression of a novel genomic region on chromosome 12, which is a QTL hotspot for all culm strength traits (section modulus, bending stress, and internode breaking weight) and grain number, can serve a dual purpose of enhancing culm strength and grain yield. MTAs found in the vicinity of previously reported putative candidate genes for lodging resistance would provide scope for further research on the biochemical basis of culm strength.

Data availability statement

The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found below: https://www.ncbi.nlm.nih.gov/sra/PRJNA1026065.

Author contributions

JB: Conceptualization, Data curation, Funding acquisition, Investigation, Methodology, Project administration, Resources, Supervision, Validation, Visualization, Writing – original draft, Writing – review & editing. RP: Data curation, Formal analysis, Investigation, Software, Writing – review & editing. AC: Data curation, Formal analysis, Software, Writing – review & editing. AM: Writing – review & editing. SI: Investigation, Writing – review & editing. SA: Resources, Writing – review & editing. RS: Funding acquisition, Project administration, Resources, Writing – review & editing.

Funding

The author(s) declare financial support was received for the research, authorship, and/or publication of this article. The research was conducted as part of ICAR-Indian Institute of Rice Research (ICAR-IIRR)- Institute/IRC Project-GEY/CI/BR/25: “Broadening the genetic base of indica rice and modify plant type by introgressing traits from Tropical japonica” and funding support from ICAR-IIRR and DBT project BT/PR32853/AGIII/103/1159/2019- “Development of Haplotype based Near Isogenic Lines (Haplo-NILs) for Enhanced Genetic Gain in Rice”.

Acknowledgments

The authors greatly acknowledge ICAR-IIRR and DBT for their funding support.

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.

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.2023.1298083/full#supplementary-material

References

Agrama, H. A., Eizenga, G. C., Yan, W. (2007). Association mapping of yield and its components in rice cultivars. Mol. Breed. 19 (4), 341–356. doi: 10.1007/s11032-006-9066-6

CrossRef Full Text | Google Scholar

Anilkumar, C., Muhammed Azharudheen, T. P., Sah, R. P., Sunitha, N. C., Devanna, B. N., Marndi, B. C., et al. (2023). Gene based markers improve precision of genome-wide association studies and accuracy of genomic predictions in rice breeding. Heredity 130 (5), 335–345. doi: 10.1038/s41437-023-00599-5

PubMed Abstract | CrossRef Full Text | Google Scholar

Asano, K., Takashi, T., Miura, K., Qian, Q., Kitano, H., Matsuoka, M., et al. (2007). Genetic and molecular analysis of utility of sd1 alleles in rice breeding. Breed. Sci. 57 (1), 53–58. doi: 10.1270/jsbbs.57.53

CrossRef Full Text | Google Scholar

Ashikari, M., Sasaki, A., Ueguchi-Tanaka, M., Itoh, H., Nishimura, A., Datta, S., et al. (2002). Loss-of-function of a rice gibberellin biosynthetic gene, GA20 oxidase (GA20ox-2), led to the rice ‘green revolution’. Breed. Sci. 52 (2), 143–150. doi: 10.1270/jsbbs.52.143

CrossRef Full Text | Google Scholar

Assefa, T., Otyama, P. I., Brown, A. V., Kalberer, S. R., Kulkarni, R. S., Cannon, S. B. (2019). Genome-wide associations and epistatic interactions for internode number, plant height, seed weight and seed yield in soybean. BMC Genomics 20, 527. doi: 10.1186/s12864-019-5907-7

PubMed Abstract | CrossRef Full Text | Google Scholar

Atwell, S., Huang, Y. S., Vilhjálmsson, B. J., Willems, G., Horton, M., Li, Y., et al. (2010). Genome-wide association study of 107 phenotypes in Arabidopsis thaliana inbred lines. Nature 465, 627–631. doi: 10.1038/nature08800

PubMed Abstract | CrossRef Full Text | Google Scholar

Berry, P. M., Spink, J. H., Gay, A. P., Craigon, J. A. (2003). Comparison of root and stem lodging risks among winter wheat cultivars. J. Agric. Sci. 141, 191–202. doi: 10.1017/S002185960300354

CrossRef Full Text | Google Scholar

Bollinedi, H., Yadav, A. K., Vinod, K. K., Gopala Krishnan, S., Bhowmick, P. K., Nagarajan, M., et al. (2020). Genome-wide association study reveals novel marker-trait associations (MTAs) governing the localization of Fe and Zn in the rice grain. Front. Genet. 11. doi: 10.3389/fgene.2020.00213

CrossRef Full Text | Google Scholar

Chen, W., Gao, Y., Xie, W., Gong, L., Lu, K., Wang, W. (2014). Genome-wide association analyses provide genetic and biochemical insights into natural variation in rice metabolism. Nat. Genet. 46 (7), 714–721. doi: 10.1038/ng.3007

PubMed Abstract | CrossRef Full Text | Google Scholar

Chen, Y., Lubberstedt, T. (2010). Molecular basis of trait correlations. Trends Plant Sci. 15, 454–461. doi: 10.1016/j.tplants.2010.05.004

PubMed Abstract | CrossRef Full Text | Google Scholar

Chigira, K., Kojima, N., Yamasaki, M., Yano, K., Adachi, S., Nomura, T., et al. (2020). Landraces of temperate japonica rice have superior alleles for improving culm strength associated with lodging resistance. Sci. Rep. 10 (1), 19855. doi: 10.1038/s41598-020-76949-8

PubMed Abstract | CrossRef Full Text | Google Scholar

Fan, C., Feng, S., Huang, J., Wang, Y., Wu, L., Li, X. (2017). At CesA8-driven OsSUS3 expression leads to largely enhanced biomass saccharification and lodging resistance by distinctively altering lignocellulose features in rice. Biotechnol. Biofuels. 10 (1), 1–12. doi: 10.1186/s13068-017-0911-0

PubMed Abstract | CrossRef Full Text | Google Scholar

Flintham, J. E., Borner, A., Worland, A. J., Gale, M. D. (1997). Optimizing wheat grain yield: effects of Rht (gibberellin-insensitive) dwarfing genes. J. Agric. Sci. 128, 11–25. doi: 10.1017/S0021859696003942

CrossRef Full Text | Google Scholar

Guo, Z., Liu, X., Zhang, B., Yuan, X., Xing, Y., Liu, H., et al. (2021). Genetic analyses of lodging resistance and yield provide insights into post-Green-Revolution breeding in rice. Plant Biotechnol. J. 19 (4), 814–829. doi: 10.1111/pbi.13509

PubMed Abstract | CrossRef Full Text | Google Scholar

Guo, Z., Yang, W., Chang, Y., Ma, X., Tu, H., Xiong, F., et al. (2018). Genome-wide association studies of image traits reveal genetic architecture of drought resistance in rice. Mol. Plant 11, 789–805. doi: 10.1016/j.molp.2018.03.018

PubMed Abstract | CrossRef Full Text | Google Scholar

Gutiérrez, J., Núñez-Flores, M. J. L., Gómez-Ros, L. V., Uzal, E. N., Carrasco, A. E., Díaz, J., et al. (2009). Hormonal regulation of the basic peroxidase isoenzyme from Zinnia elegans. Planta 230, 767–778. doi: 10.1007/s00425-009-0982-7

PubMed Abstract | CrossRef Full Text | Google Scholar

Hai, L., Guo, H., Xiao, S., Jiang, G., Zhang, X., Yan, C., et al. (2005). Quantitative trait loci (QTL) of stem strength and related traits in a doubled-haploid population of wheat (Triticum aestivum L.). Euphytica 141, 1–9. doi: 10.1007/s10681-005-4713-2

CrossRef Full Text | Google Scholar

Hijazi, M., Velasquez., S. M., Jamet, E., Estevez, J. M., Albenne, C. (2014). An update on post translational modifications of hrdoxyproline-rich glycoproteins: toward a model highlighting their contribution to plant cell wall architecture. Front. Plant Sci. 5. doi: 10.3389/fpls.2014.00395

PubMed Abstract | CrossRef Full Text | Google Scholar

Hirano, K., Ordonio, R. L., Matsuoka, M. (2017). Engineering the lodging resistance mechanism of post-Green Revolution rice to meet future demands. Phys. Biol. Sci. 93, 220–233. doi: 10.2183/pjab.93.014

CrossRef Full Text | Google Scholar

Huang, M., Liu, X., Zhou, Y., Summers, R. M., Zhang, Z. (2019). BLINK: a package for the next level of genome-wide association studies with both individuals and markers in the millions. Gigascience 8 (2), 154. doi: 10.1093/gigascience/giy154

CrossRef Full Text | Google Scholar

Huang, X., Wei, X., Sang, T., Zhao, Q., Feng, Q., Zhao, Y. (2010). Genome-wide association studies of 14 agronomic traits in rice landraces. Nat. Genet. 42 (11), 961–967. doi: 10.1038/ng.695

PubMed Abstract | CrossRef Full Text | Google Scholar

Huang, X., Zhao, Y., Li, C., Wang, A., Zhao, Q., Li, W., et al. (2012). Genome-wide association study of flowering time and grain yield traits in a worldwide collection of rice germplasm. Nat. Genet. 44 (1), 32–39. doi: 10.1038/ng.1018

CrossRef Full Text | Google Scholar

Ishimaru, K., Togawa, E., Ookawa, T., Kashiwagi, T., Madoka, Y., Hirotsu, N. (2008). New target for rice lodging resistance and its effect in a typhoon. Planta 277, 601–609. doi: 10.1007/s00425-007-0642-8

CrossRef Full Text | Google Scholar

Islam, M. S., Peng, S. B., Visperas, R. M., Ereful, N., Bhuiya, M. S. U., Julfiquar, A. W. (2007). Lodging-related morphological traits of hybrid rice in a tropical irrigated ecosystem. Field Crops Res. 101 (2), 240–248. doi: 10.1016/j.fcr.2006.12.002

CrossRef Full Text | Google Scholar

Jyothi, B., Divya, B., Subba Rao, L. V., Laxmi Bhavani, P., Revathi, P., Raghuveer Rao, P., et al. (2018). New plant type trait characterization and development of core set among indica and tropical japonica genotypes of rice. Plant Genet. Resour. 16 (6), 504–512. doi: 10.1017/S1479262118000084

CrossRef Full Text | Google Scholar

Kashiwagi, T., Hirotsu, N., Ujiie, K., Ishimaru, K. (2010). Lodging resistance locus prl5 improves physical strength of the lower plant part under different conditions of fertilization in rice (Oryza sativa L.). Field Crops Res. 115 (1), 107–115. doi: 10.1016/j.fcr.2009.10.011

CrossRef Full Text | Google Scholar

Kashiwagi, T., Ishimaru, K. (2004). Identification and functional analysis of a locus for improvement of lodging resistance in rice. Plant Physiol. 134, 676–683. doi: 10.1104/pp.103.029355

PubMed Abstract | CrossRef Full Text | Google Scholar

Kashiwagi, T., Sasaki, H., Ishimaru, K. (2005). Factors responsible for decreasing sturdiness of the lower part in lodging of rice (Oryza sativa L.). Plant Production. Sci. 8 (2), 166–172. doi: 10.1626/pps.8.166

CrossRef Full Text | Google Scholar

Kashiwagi, T., Togawa, E., Hirotsu, N., Ishimaru, K. (2008). Improvement of lodging resistance with QTLs for stem diameter in rice (Oryza sativa L.). Theor. Appl. Genet. 117, 749–757. doi: 10.1007/s00122-008-0816-1

PubMed Abstract | CrossRef Full Text | Google Scholar

Keller, M., Karutz, Ch., Schmid, J. E., Stamp, P., Winzeler, M., Keller, B., et al. (1999). Quantitative trait loci for lodging resistance in a segregating wheat × spelt population. Theor. Appl. Genet. 98, 1171–1182. doi: 10.1007/s001220051182

CrossRef Full Text | Google Scholar

Kotake, T., Aohara, T., Hirano, K., Sato, A., Kaneko, Y., Tsumuraya, Y., et al. (2011). Rice Brittle culm 6 encodes a dominant-negative form of CesA protein that perturbs cellulose synthesis in secondary cell walls. J. Exp. Bot. 62, 2053–2062. doi: 10.1093/jxb/erq395

PubMed Abstract | CrossRef Full Text | Google Scholar

Lang, Y., Yang, X., Wang, M., Zhu, Q. (2012). Effects of lodging at different filling stages on rice yield and grain quality. Rice Sci. 19 (4), 315–319. doi: 10.1016/S1672-6308(12)60056-0

CrossRef Full Text | Google Scholar

Lipka, A. E., Tian, F., Wang, Q., Peiffer, J., Li, M., Bradbury, P. J., et al. (2012). GAPIT: genome association and prediction integrated tool. Bioinform 28 (18), 2397–2399. doi: 10.1093/bioinformatics/bts444

CrossRef Full Text | Google Scholar

Luo, X., Wu, Z., Fu, L., Dan, Z., Yuan, Z., Liang, T., et al. (2022). Evaluation of lodging resistance in rice based on an optimized parameter from lodging index. Crop Sci. 62 (3), 1318–1332. doi: 10.1002/csc2.20712

CrossRef Full Text | Google Scholar

Ma, Y., Stafford, L., Ratcliffe, J., Bacic, A., Johnson, K. L. (2022). WAKL8 regulates arabidopsis stem secondary wall development. Plants 11, 2297. doi: 10.3390/11172297

PubMed Abstract | CrossRef Full Text | Google Scholar

Mariani, L., Ferrante, A. (2017). Agronomic management for enhancing plant tolerance to abiotic stresses—drought, salinity, hypoxia, and lodging. Horticulturae 3, 52. doi: 10.3390/horticulturae3040052

CrossRef Full Text | Google Scholar

Meng, B., Wang, T., Luo, Y., Xu, D., Li, L., Diao, Y., et al. (2021). Genome-wide association study tdentified novel candidate loci/genes affecting lodging resistance in rice. Genes 12 (5), 718. doi: 10.3390/genes12050718

PubMed Abstract | CrossRef Full Text | Google Scholar

Mu, P., Li, Z., Li, C., Wang, X. (2004). QTL analysis for lodging resistance in rice using a DH population under lowland and upland ecosystems. Acta.Gen. Sin. 31, 717–723.

Google Scholar

Murai, M., Takamura, I., Sato, S., Tokutome, T., Sato, Y. (2002). Effects of the dwarfing gene originating from Dee–geo–woo–gen on yield and its related traits in rice. Breed. Sci. 52, 95–100. doi: 10.1270/jsbbs.52.95

CrossRef Full Text | Google Scholar

Nakajima, T., Yoshida, M., Tomimura, K. (2008). Effect of lodging on the level of mycotoxins in wheat, barley, and rice infected with the Fusarium graminearum species complex. J. Gen. Plant Pathol. 74 (4), 289–295. doi: 10.1007/s10327-008-0103

CrossRef Full Text | Google Scholar

Nayak, A. K., Anilkumar, C., Behera, S., Sah, R. P., Lavanya, G. R., Kumar, A. (2022). Genetic dissection of grain size traits through genome-wide association study based on genic markers in rice. Rice Sci. 29 (5), 462–472. doi: 10.1016/j.rsci.2022.07.006

CrossRef Full Text | Google Scholar

Noda, S., Koshiba, T., Hattori, T., Yamaguchi, M., Suzuki, S., Umezawa, T. (2015). The expression of a rice secondary wall-specific cellulose synthase gene, OsCesA7, is directly regulated by a rice transcription factor, OsMYB58/63. Planta 242, 589–600. doi: 10.1007/s00425-015-2343-z

PubMed Abstract | CrossRef Full Text | Google Scholar

Nomura, T., Arakawa, N., Yamamoto, T., Ueda, T., Adachi, S., Yonemaru, J. I., et al. (2019). Next generation long-culm rice with superior lodging resistance and high grain yield, Monster Rice 1. PloS One 14 (8), 0221424. doi: 10.1371/journal.pone.0221424

CrossRef Full Text | Google Scholar

Nomura, T., Seki, Y., Matsuoka, M., Yano, K., Chigira, K., Adachi, S., et al. (2021). Potential of rice landraces with strong culms as genetic resources for improving lodging resistance against super typhoons. Sci. Rep. 11 (1), 15780. doi: 10.1038/s41598-021-95268-0

PubMed Abstract | CrossRef Full Text | Google Scholar

Okuno, A., Hirano, K., Asano, K., Takase, W., Masuda, R., Morinaka, Y., et al. (2014). New approach to increasing rice lodging resistance and biomass yield through the use of high gibberellin producing varieties. PloS One 9 (2), 86870. doi: 10.1371/journal.pone.0086870

CrossRef Full Text | Google Scholar

Ookawa, T., Aoba, R., Yamamoto, T., Ueda, T., Takai, T., Fukuoka, et al. (2016). Precise estimation of genomic regions controlling lodging resistance using a set of reciprocal chromosome segment substitution lines in rice. Sci. Rep. 6 (1), 30572. doi: 10.1038/srep30572

PubMed Abstract | CrossRef Full Text | Google Scholar

Ookawa, T., Hobo, T., Yano, M., Murata, K., Ando, T., Miura, H., et al. (2010). New approach for rice improvement using a pleiotropic QTL gene for lodging resistance and yield. Nat. Commun. 1, 132. doi: 10.1038/ncomms1132

PubMed Abstract | CrossRef Full Text | Google Scholar

Ookawa, T., Ishihara, K. (1992). Varietal difference of physical characteristics of the culm related to lodging resistance in paddy rice. Japanese. J. Crop Sci. 61 (3), 419–425. doi: 10.1626/jcs.61.419

CrossRef Full Text | Google Scholar

Peng, S., Khush, G. S. (2003). Four decades of breeding for varietal improvement of irrigated lowland rice in the International Rice Research Institute. Plant Prod. Sci. 6, 157–164. doi: 10.1626/6.157

CrossRef Full Text | Google Scholar

Peng, S., Laza, R. C., Visperas., R. M., Sanico, A. L., Cassmam, K. G., Khush, G. S. (2000). Grain yield of rice cultivars and lines developed in the Philippines since 1966. Crop Sci. 40, 307–314. doi: 10.2135/cropsci2000.402307x

CrossRef Full Text | Google Scholar

Qaseem, M. F., Wu, A. M. (2020). Balanced xylan acetylation is the key regulator of plant growth and development, and cell wall structure and for industrial utilization. Int. J. Mol. Sci. 21 (21), 7875. doi: 10.3390/21217875

PubMed Abstract | CrossRef Full Text | Google Scholar

Rashid, M. A. R., Zhao, Y., Azeem, F., Zhao, Y., Ahmed, H. G. M., Atif, R. M., et al. (2022). Unveiling the genetic architecture for lodging resistance in rice (Oryza sativa. L) by genome-wide association analyses. Front. Genet. 13. doi: 10.3389/2022.960007

CrossRef Full Text | Google Scholar

Sasaki, A., Ashikari, M., Ueguchi-Tanaka, M., Itoh, H., Nishimura, A., Swapan, D., et al. (2002). Green revolution: a mutant gibberellin-synthesis gene in rice. Nature 416, 701–702. doi: 10.1038/416701a

PubMed Abstract | CrossRef Full Text | Google Scholar

Setter, T. L., Laureles, E. V., Mazaredo, A. M. (1997). Lodging reduces yield of rice by self-shading and reductions in canopy photosynthesis. Field Crops Res. 49 (2-3), 95–106. doi: 10.1016/S0378-4290(96)01058-1

CrossRef Full Text | Google Scholar

Tanaka, K., Murata, K., Yamazaki, M., Onosato, K., Miyao, A., Hirochika, H. (2003). Three distinct rice cellulose synthase catalytic subunit genes required for cellulose synthesis in the secondary wall. Plant Physiol. 133 (1), 73–83. doi: 10.1104/pp.103.022442

PubMed Abstract | CrossRef Full Text | Google Scholar

Terashima, K., Ogata, T., Akita, S. (1992). Eco–physiological characteristics related with lodging tolerance of rice in direct sowing cultivation in comparison of root lodging tolerance among cultivars by the measurement of pushing resistance. Jpn. J. Crop Sci. 61, 380–387. doi: 10.1626/jcs.61.380

CrossRef Full Text | Google Scholar

Wang, L., Guo, K., Li, Y., Tu, Y., Hu, H., Wang, B., et al. (2010). Expression profiling and integrative analysis of the CESA/CSL superfamily in rice. BMC Plant Biol. 10, 282. 175. doi: 10.1186/1471-2229-10-282

CrossRef Full Text | Google Scholar

Wang, Y. H., Zheng, Y. M., Cai, Q. H., Liao, C. J., Mao, X. H., Xie, H. G., et al. (2016). Population structure and association analysis of yield and grain quality traits in hybrid rice primal parental lines. Euphytica 212 (2), 261–273. doi: 10.1007/s10681-016-1766-3

CrossRef Full Text | Google Scholar

Won, J. G., Hirahara, Y., Yoshida, T., Imabayashi, S. (1998). Selection of rice lines using SPGP seedling method for direct seeding. Plant Prod. Sci. 1, 280–285. doi: 10.1626/pps.1.280

CrossRef Full Text | Google Scholar

Xie, W., Wang, G., Yuan, M., Yao, W., Lyu, K., Zhao, H. (2015). Breeding signatures of rice improvement revealed by a genomic variation map from a large germplasm collection. Proc. Natl. Acad. Sci. 112 (39), E5411–E5419. doi: 10.1073/pnas.1515919112

CrossRef Full Text | Google Scholar

Yadav, S., Sandhu, N., Majumder, R. R., Dixit, S., Kumar, S., Singh, S. P., et al. (2019). Epistatic interactions of major effect drought QTLs with genetic background loci determine grain yield of rice under drought stress. Sci. Rep. 9 (1), 2616. doi: 10.1038/s41598-019-39084-7

PubMed Abstract | CrossRef Full Text | Google Scholar

Yadav, S., Singh, U. M., Naik, S. M., Venkateshwarlu, C., Ramayya, P. J., Raman, A., et al. (2017). Molecular mapping of QTLs associated with lodging resistance in dry direct seeded rice (Oryza sativa L.). Front. Plant Sci. 8. doi: 10.3389/01431

CrossRef Full Text | Google Scholar

Yang, J., Ferreira, T., Morris, A. P., Medland, S. E., Madden, P. A. F., Heath, A. C., et al. (2012). Conditional and joint multiple-SNP analysis of GWAS summary statistics identifies additional variants influencing complex traits. Nat. Genet. 44, 369–375. doi: 10.1038/ng.2213

PubMed Abstract | CrossRef Full Text | Google Scholar

Yang, X., Lai, Y., Wang, L., Zhao, M., Wang, J., Li, M., et al. (2023). Isolation of a novel QTL, qSCM4, associated with strong culm affects lodging resistance and panicle branch number in rice. International Journal of Molecular Sciences 24 (1), 812. doi: 10.3390/ijms24010812

PubMed Abstract | CrossRef Full Text | Google Scholar

Yano, K., Ookawa, T., Aya, K., Ochiai, Y., Hirasawa, T., Ebitani, T., et al. (2015). Isolation of a novel lodging resistance QTL gene involved in strigolactone signaling and its pyramiding with a QTL gene involved in another mechanism. Mol. Plant 8, 303–314. doi: 10.1016/j.molp.10.009

PubMed Abstract | CrossRef Full Text | Google Scholar

Yu, J., Pressoir, G., Briggs, W. H., Vroh Bi, I., Yamasaki, M., Doebley, J. F., et al. (2006). A unified mixed-model method for association mapping that accounts for multiple levels of relatedness. Nat. Genet. 38 (2), 203–208. doi: 10.1038/ng1702

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhang, P., Liu, X. D., Tong, H. H., Lu, Y. G., Li, J. Q. (2014). Association mapping for important agronomic traits in core collection of rice (Oryza sativa L.) with SSR markers. PloS One 9 (10), 111508. doi: 10.1371/journal.pone.0111508

CrossRef Full Text | Google Scholar

Zhao, D. D., Jang, Y. H., Kim, E. G., Park, J. R., Jan, R., Lubna, et al. (2023). Identification of a major locus for lodging resistance to typhoons using QTL Analysis in rice. Plants 12 (3), 449. doi: 10.3390/plants12030449

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhao, D. D., Son, J. H., Lee, G. S., Kim, K. M. (2021). Screening for a novel gene, OsPSLSq6, using QTL analysis for lodging resistance in rice. Agronomy 11, 334. doi: 10.3390/agronomy11020334

CrossRef Full Text | Google Scholar

Zhu, L. H., Zhong, D. B., Xu, J. L., Yu, S. B. (2008). Differential expression of lodging resistance related QTL in rice (Oryza sativa L.). Plant Sci. 175, 898–905. doi: 10.1016/j.plantsci.09.001

CrossRef Full Text | Google Scholar