Novel Genes and Genetic Loci Associated With Root Morphological Traits, Phosphorus-Acquisition Efficiency and Phosphorus-Use Efficiency in Chickpea

Chickpea—the second most important grain legume worldwide—is cultivated mainly on marginal soils. Phosphorus (P) deficiency often restricts chickpea yields. Understanding the genetics of traits encoding P-acquisition efficiency and P-use efficiency will help develop strategies to reduce P-fertilizer application. A genome-wide association mapping approach was used to determine loci and genes associated with root architecture, root traits associated with P-acquisition efficiency and P-use efficiency, and any associated proxy traits. Using three statistical models—a generalized linear model (GLM), a mixed linear model (MLM), and a fixed and random model circulating probability unification (FarmCPU) —10, 51, and 40 marker-trait associations (MTAs), respectively were identified. A single nucleotide polymorphism (SNP) locus (Ca1_12310101) on Ca1 associated with three traits, i.e., physiological P-use efficiency, shoot dry weight, and shoot P content was identified. Genes related to shoot P concentration (NAD kinase 2, dynamin-related protein 1C), physiological P-use efficiency (fasciclin-like arabinogalactan protein), specific root length (4-coumarate–CoA ligase 1) and manganese concentration in mature leaves (ABC1 family protein) were identified. The MTAs and novel genes identified in this study can be used to improve P-use efficiency in chickpea.

Introduction

Phosphorus (P) is an essential nutrient for crop production. Using spatially explicit global maps for more than 100 crops, MacDonald et al. (2011) reported that 29% of the global cropland area is P deficient, while 71% has surplus P. For example, 42% of field soil in India is low in P, 38% is medium and 20% is high (Motsara, 2002). Excessive P fertilizer and manure application in industrialized countries have led to low P-use efficiency (PUE), with surplus P retained in soil as residual P (Syers et al., 2008) or lost to the environment where it is causing significant water quality problems. Globally 51–86% more P input will be needed for sustainable crop production by 2050 (Mogollón et al., 2018), unless we work toward more P-efficient crops and cropping systems (Cong et al., 2020).

Chickpea (Cicer arietinum L.), is one of the most important grain legumes cultivated by smallholder farmers in more than 50 countries. Advances in chickpea genomics in the last decade have made large-scale genomic resources available to the chickpea research community including molecular markers (Thudi et al., 2011; Hiremath et al., 2012), genetic maps (Nayak et al., 2010), genome sequences (Varshney et al., 2013b), and resequencing of several germplasm lines (Thudi et al., 2016a, b; Varshney et al., 2019). These resources have improved our understanding of both abiotic (Varshney et al., 2014) and biotic stress tolerance in chickpea (Sabbavarapu et al., 2013) and enabled fine mapping of traits (Jaganathan et al., 2015; Kale et al., 2015). Furthermore, the resources have been successfully used to develop new varieties with enhanced tolerance or resistance (see Varshney et al., 2013a; Mannur et al., 2019; Bharadwaj et al., 2020; Roorkiwal et al., 2020). Twelve chickpea genotypes with well-known responses under drought and irrigation were evaluated for profuse root length density (RLD) in surface soil and root dry weight (RDW) and root:shoot ratio (RSR) in deeper soil layers (Purushothaman et al., 2017). This study revealed that drought stress increased RLD below 300 mm soil depth, deep RDW, and RSR, but decreased root diameter. Moisture-conservation practices and optimum P levels to enhance PUE were reported recently in chickpea grown in vertisols in central India (Chaudhary et al., 2018). In soybean, root length is positively correlated with P accumulation in well-watered and water-stressed conditions (He et al., 2017). Phenotypic plasticity and genetic variability in root architectural traits of chickpea and their role in drought tolerance using a novel semi-hydroponic system have been reported (Chen et al., 2017). A recent study revealed that root system plasticity affects P acquisition efficiency, PUE, and photosynthetic PUE in 266 chickpea genotypes (Pang et al., 2018c). Root exudates like carboxylate enable P acquisition from the soils which are low in available P. Manganese concentration in mature chickpea leaves is positively correlated with the amount of rhizosheath carboxylates, offering an easily measurable proxy for assessing rhizosheath carboxylates in 100 chickpea genotypes grown under low P availability (Pang et al., 2018a). In the context of global climate change, Pang et al. (2018b) summarized the factors affecting PUE and enhancing P-acquisition efficiency in legumes, and key areas for future research. Understanding the genetics of these traits, identification of genomic regions, molecular markers and or marker trait associations for PUE efficiency related traits will help improving these traits through marker-assisted selection or genomics-assisted breeding.

Various type of molecular markers have been used to establish marker-trait associations (MTA) in chickpea. For instance, 1,072 Diversity Arrays Technology (DArTs), 651 single nucleotide polymorphisms (SNPs), 113 gene-based SNPs, and 36 simple sequence repeats (SSRs) were used to establish 312 MTAs for drought- and heat-tolerance related traits (Thudi et al., 2014). SNPs were used extensively for association studies in chickpea due to their abundance and amenability for high-throughput genotyping (Diapari et al., 2014; Jadhav et al., 2015; Upadhyaya et al., 2016; Varshney et al., 2019; Sab et al., 2020).

In legumes, genome-wide association studies (GWAS) are gaining momentum. A recent study reported two SNP markers tightly linked to seed iron (Fe) and one to seed zinc (Zn) concentration in lentils (Lens culinaris Medik.) (Khazaei et al., 2017). Another study reported 159 quantitative trait nucleotides (QTNs) and 52 candidate genes associated with the photosynthetic response to low-P stress in soybean [Glycine max (L.) Merr.] (Lü et al., 2018). In chickpea, association mapping between SNP markers and seed copper (Cu), P, and potassium (K) concentrations identified eight SNPs associated with variation in three nutrients in more than two environments (Ozkuru et al., 2018). Similarly, seed mineral concentration in pea (Pisum sativum L.; Gali et al., 2019), Fe chlorosis in soybean; Mamidi et al., 2014; Assefa et al., 2020), and Fe bioavailability in cooked dry beans (Phaseolus vulgaris L.; Katuuramu et al., 2018) have been mapped. Novel genes involved in the accumulation of P in Lotus japonicas have been reported using GWAS analysis (Giovannetti et al., 2019).

In view of above, this study was conducted to undertake GWAS analysis to identify MTAs for (i) root architectural traits evaluated in a high-throughput semi-hydroponic root phenotyping platform, and (ii) root morphological and physiological traits related to P-acquisition efficiency and P−use efficiency under low P supply. This is the first study that reports genomic regions associated with above mentioned traits by using three different models, namely fixed and random model circulating probability unification (FarmCPU), mixed linear model (MLM), and generalized linear model (GLM) in the GAPIT-R package. The reference genome has been used to identify the candidate genes in the identified MTAs associated with above traits.

Materials and Methods

Germplasm Lines and Phenotyping

The chickpea reference set (Upadhyaya et al., 2008) comprising of 300 diverse accessions (267 landraces, 13 advanced lines and cultivars, 7 wild Cicer accessions, and 13 accessions with unknown biological status) was used for phenotyping datasets as following: (i) 270 genotypes of the reference set were evaluated for 30 root architectural traits in a semi-hydroponic phenotyping system (Chen et al., 2017); (ii) Two hundred and sixty-six genotypes (including 255 from the chickpea reference set along with 11 Australian chickpea cultivars (Ambar, Almaz, Neelam, Genesis 079, Genesis 090, Genesis 509, Genesis 836, Genesis Kalkee, PBA Boundary, PBA Slasher and PBA Striker) were evaluated for P-acquisition efficiency and P−use efficiency with P supplied as insoluble FePO₄ (Pang et al., 2018c); and (iii) a selected subset of 100 chickpea genotypes of the reference set—showing visual differences in plant size or leaf symptoms of P deficiency—phenotyped for shoot/root morphological and physiological traits to understand the relative roles of root morphology and physiology in P-acquisition efficiency (Pang et al., 2018a, c).

Genotyping and Determining Population Structure

The SNP dataset, based on whole-genome resequencing data on the reference dataset was filtered for missing values (≥20%) and minor-allele frequency <5% using vcftools and imputed by BEAGLE-5.0 (Browning and Browning, 2016; Varshney et al., 2019). As the number of genotypes were different in data sets, Population structure was separately determined for the chickpea reference set and 91 genotypes (from the 100 selected subset of 199 genotypes) using ADMIXTURE v1.3.0 (Zhou et al., 2011).

Genome-Wide Association Analysis

GWAS for Root Traits and Seven P-Acquisition Efficiency Traits

To identify significant MTAs and avoid spurious associations of the 270 genotypes mentioned above, 233 genotypes (with phenotyping data and genotyping data) were considered for GWAS analysis. GWAS analysis was performed using 698,183 SNPs (SNP calls obtained based on aligning 233 genotypes to reference genome CDC Frontier) and phenotyping data for 37 traits (i.e., 30 root traits and seven P-acquisition efficiency traits). MTAs were determined using three models, namely fixed and random model circulating probability unification (FarmCPU), mixed linear model (MLM), and generalized linear model (GLM) in the GAPIT-R package (Lipka et al., 2012).

GWAS for Biochemical Traits/Proxy Traits for PUE

A total of 706,865 SNPs and phenotyping data generated on 91 chickpea germplasm lines was for analysis in GAPIT-R using FarmCPU and MLM models to determine MTAs for biochemical traits.

The Bonferroni correction threshold of 7.07E-08 was used to avoid spurious associations. The genes involving significant SNP markers were aligned against the NCBI non-redundant (nr) protein database taxon Viridiplantae using BLASTX, to obtain functional annotations. GO and KEGG (Kyoto Encyclopedia of Genes and Genomes) pathway identification were conducted on these sequences in the KEGG pathways in-built in BLAST2GO. The SNPEff- 4.3T open source program was used for variant annotation and prediction of significant SNP effects.

Results

Population Structure and Genome-Wide Association Study (GWAS)

Using genome-wide SNP data for 233 chickpea genotypes (with both genotypic and phenotypic data), three subpopulations were identified (Supplementary Figures 1A,B) using ADMIXTURE. Similarly, in the smaller subset of 91 genotypes, three subpopulations were identified (Supplementary Figures 1C,D). Three statistical models—GLM, MLM, and FarmCPU—enabled identification of 10, 51, and 40 MTAs, respectively (Table 1 and Supplementary Table 1) after applying Bonferroni and FDR corrections. Spurious MTAs were excluded by examining Q/Q plots. Forty-two of the MTAs identified by the GLM and MLM models were robust (>10% of the phenotypic variation explained); the phenotypic variation explained by the MTAs identified in the FarmCPU model was not computed. The GLM model identified MTAs for four traits (root growth rate, root mass density, specific root length, and shoot P concentration); this low number of traits could be due to the many spurious MTAs (Supplementary Figure 2). The FarmCPU model identified MTAs for the most traits, followed by the MLM model.

TABLE 1

Table 1. Summary of marker trait associations identified using three statistical models namely GLM, MLM and FarmCPU.

GWAS Signals for Root Architectural Traits

Fifty-seven MTAs were identified for 19 root architectural traits phenotyped in our earlier studies (Chen et al., 2017; Pang et al., 2018c; Table 1). The FarmCPU model identified one significant MTA each in root depth zone 1 (TRL_z1; cm) and root depth zone 2 (TRL_z1; cm) on Ca3 and Ca6 for taproot length (i.e., root depth, z1 + z2), one significant MTA (Ca3_26114159) for root mass (RM, mg), two significant MTAs (Ca4_29694614 and Ca6_57970784) for branch density (BD; cm^–1 taproot), and one significant MTA (Ca2_6340118) for branch intensity (BI; cm^–1 root). No significant MTAs were detected for the above root traits in the GLM and MLM methods. Nevertheless, the MLM model identified six significant MTAs and FarmCPU identified one MTA for root mass ratio (RMR), and the MLM and FarmCPU models identified two significant MTAs (one each) for root length ratio (RLR_top/sub) and three significant MTAs (two and one, respectively) for average branch length (ABL; cm per branch). The FarmCPU model identified one significant MTA each for root diameter s2 (RD_s2; mm), subsoil root diameter (RD_sub; mm), and topsoil root diameter (RD_top; mm) on Ca6, Ca1, and Ca2, respectively. The FarmCPU model also identified two significant MTAs for water-use efficiency (WUE; A/g_s) on Ca6 and Ca8. Based on the physical position of the SNP loci (MTAs) associated with different traits on Ca4, none of the MTAs was mapped in “QTL-hotspot_a” or “QTL-hotspot_b” that harbors several drought-tolerance-related root traits (Kale et al., 2015). The MTAs revealed key root traits for efficiently acquiring soil resources and adapting to drought and other abiotic stresses. Of the 13 significant MTAs (12 from MLM and one from GLM) identified for RGR (cm d^–1), six (46.2%) were identified on Ca3 and three (23.1%) on Ca6. Six MTAs identified on Ca6 explained 3.6–10.0% of the phenotypic variation associated with five genes (Ca_08259, Ca_01151, Ca_01152, Ca_01156, and Ca_16553). Both statistical models (GLM and MLM) detected one significant MTA (Ca1_4716136), which explained higher phenotypic variation (11.8–13.1%) than the other MTAs (Supplementary Table 1). Eight significant MTAs were identified for specific root length (SRL, cm) in all three models; of these, one MTA on Ca7 (Ca7_3606123), explaining about 17.8% of the phenotypic variation, was consistent in all three models used for GWAS analysis. The MTA was associated with the Ca_03107 gene that encodes pectinesterase/pectinesterase inhibitor PPE8B. Six MTAs for total rhizosheath dry soil (g plant^–1) were identified on Ca8 (3), Ca5 (1), Ca2 (1), and Ca1 (1).

GWAS Signals for P-Acquisition Efficiency and PUE-Related Traits

A total of 10 significant MTAs (three based on GLM and seven based on MLM) for shoot P content—six on Ca4, two on Ca6 and one each on Ca7 and Ca2 were identified (Figure 1A). Two significant MTAs on Ca4 (Ca4_38518152 and Ca4_8269508) were identified with the GLM and MLM models (Supplementary Table 1). The 15 MTAs identified—14 from MLM and one from FarmCPU—explained 5–21% of the phenotypic variation. Five of the MTAs were identified on Ca6 followed by Ca7 (4), Ca4 (2), Ca2 (2), Ca1 (1), and C3 (1) (Figure 1B). A SNP locus, Ca7_33808891, associated with physiological P-Use efficiency present in gene Ca_16189 on Ca7 explained 20.23% PVE (Figure 1C; Supplementary Table 1). A SNP locus, Ca4_38518152, associated with Shoot phosphorus content present in gene Ca_13110 on Ca4 explained 31.5% PVE (Figure 1D; Supplementary Table 1). In the case of physiological PUE one significant MTA (Ca1_12310101) on Ca1 explaining 12.0% of the phenotypic variation was identified. A single nucleotide polymorphism locus (Ca1_12310101) on Ca1 associated with three traits, i.e., physiological P-use efficiency, shoot dry weight, and shoot P content was identified (Figure 2).

FIGURE 1

Figure 1. GWAS signals for phosphorus (P) -acquisition and P-use efficiency-related traits. (A) Shoot P concentration, (B) physiological P-use efficiency, (C) P-utilization efficiency, and (D) total shoot P content. Three different statistical models GLM, MLM, and FarmCPU were used to identify the MTAS. The significant MTAs were determined using Bonferroni correction (Table 1).

FIGURE 2

Figure 2. An SNP locus (Ca1_12310101) on Ca1 showing association with three traits (A) physiological phosphorus (P)-use efficiency, (B) shoot dry weight, and (C) shoot P content. Three different statistical models GLM, MLM, and FarmCPU were used to identify the MTAS. The significant MTAs were determined using Bonferroni correction (Table 1).

GWAS Signals for Proxy Traits

One significant MTA for P concentration in mature leaves (P_ML; mg g^–1) was identified on Ca4. Seven significant MTAs for Mn concentration in mature leaves were identified on Ca2 (3), Ca4 (2), and Ca7 (Figure 3). The SNP loci on Ca7 (Ca7_32383349) and Ca4 (Ca4_1791932) were associated with manganese concentration in mature leaves with two different p-values (Figure 3), while three SNPs associated with the Mn_ML on Ca2 (Ca2_7561143, Ca2_866639, and Ca2_359984). One MTA each for citric (μmol plant^–1), Ci, carboxylate conc (μmol g^–1 root dw), malonic (μmol plant^–1), Pn_area (μmol m^–2s^–1), Pn_mass (μmol g^–1 s^–1), specific rhizosheath weight (g g^–1root DW) were identified using FarmCPU. No significant MTAs were identified for these traits using GLM and MLM models.

FIGURE 3

Figure 3. GWAS signal for manganese concentration in mature leaves. Of seven MTAs identified, three on Ca2 were unique and two each on Ca4 and Ca7 were the same SNP loci associated with the trait at different significance levels. (A) Manhattan plot showing significant MTAs and (B) Q-Q plots for the trait. Three different statistical models GLM, MLM, and FarmCPU were used to identify the MTAS. The significant MTAs were determined using Bonferroni correction (Table 1).

Discussion

Understanding the genetics of traits associated with enhanced PUE, including P-acquisition efficiency and P-use efficiency, is essential for its manipulation. In our earlier studies on chickpea, we gained insight into root system architecture, shoot/root traits associated with P-acquisition efficiency and P-use efficiency, and associated proxy traits. Many SNPs associated with root-hair length in P-sufficient and P-deficient conditions have been reported recently (Kohli et al., 2020); here we report SNPs associated with PUE and P-acquisition efficiency related root traits as well as proxy traits, which can be deployed for breeding.

As population structure is important for avoiding spurious associations, we identified three subpopulations using ADMIXTURE in the reference set (233 genotypes) and a smaller subset (91 genotypes). Similarly, earlier studies reported three subpopulations using different marker systems (Thudi et al., 2014) and genome-wide SNP markers (Varshney et al., 2019). In the smaller subset of 91 genotypes, we also identified three subpopulations (Supplementary Figures 1A–D). Four subpopulations were recently reported in a diverse set of 186 genotypes (including 20 Iranian landraces and 166 Kabuli advanced breeding lines from ICRISAT and ICARDA) using DArTseq markers (Farahani et al., 2019). Two subpopulations were reported in a set of 92 (77 landraces and five elite cultivars) chickpea germplasm lines that represent arid, semi-arid, and tropical climates using 8,113 genotyping-by-sequencing based SNPs (Sani et al., 2017).

Three statistical models (GLM, MLM, and FarmCPU) were used to identify genome-wide association signals for root architectural traits, P-acquisition efficiency and PUE-related traits. Similarly, a recent study used different statistical models to compare selected traits with different heritabilities in soybean and maize—the FarmCPU model provided a closer number of QTL than those in the literature and known genomic regions (Kaler et al., 2020). In this study, identification of significant MTAs and determining their usefulness for chickpea improvement was a major focus, rather than comparative utilities of different statistical models. In the present study, significant MTAs for 17 of the 30 root traits analyzed was reported. Earlier, in chickpea, MTAs for drought-tolerance-related root traits were based on screening a reference set in polyvinyl chloride pipes or semi-automated root screening facility at ICRISAT (Thudi et al., 2014; Varshney et al., 2019). In this study, MTAs for root traits studied under semi-hydroponic conditions and the role of root and shoot traits in P-acquisition efficiency and PUE are reported.

The MTA (Ca1_4716136) identified for RGR, is present in a gene (Ca_00555) that encodes a receptor-like cytosolic serine/threonine-protein kinase RBK2 involved in protein phosphorylation in Arabidopsis thaliana¹. Rac-like GTP-binding protein (ARAC5), with 87.8% homology, is localized on the plasma membrane of Arabidopsis thaliana root tips. The Ca_01156 gene encodes Acyl-CoA-binding domain-containing protein 4; in Arabidopsis thaliana; overexpression of Acyl-CoA -binding protein 3 (ACBP3) results in leaf senescence (Xiao et al., 2010). In chickpea, the pectinesterase gene was downregulated in the roots of salt-tolerant genotypes (Kaashyap et al., 2018). Pectinestarases also play a key role in root-hair initiation and elongation (Cosgrove, 2016). Similarly, another MTA on Ca1 (Ca1_8712480) explained about 20% of the phenotypic variation detected in the GLM and MLM models. The MTA was present in the Ca_02941 gene that encodes 4-coumarate–CoA ligase 1. In tobacco, root length increased by 64% compared with the wildtype, on overexpression of Fm4CL-like 1 [4-coumarate:CoA ligase 4 (4CL-like 1) from Fraxinus mandshurica] under mannitol-simulated drought stress (Chen et al., 2019).

Two MTAs (Ca4_38518152 and Ca4_8269508) for shoot P content were present in two genes, Ca_13110 and Ca_08315, which encode NAD kinase 2, chloroplastic-like isoform X2 and piezo-type mechanosensitive ion channel homolog isoform X2, respectively. NAD kinase 2 is involved in phosphorylation. In general, NADK genes show tissue specificity in expression. In Arabidopsis, NADK2 is expressed in leaves, while TaNADK2 is highly expressed in wheat pistils, caryopses, and endosperm during the reproductive stage (Li et al., 2018). All MTAs for shoot P content were robust and explained 27.7–41.4% of the phenotypic variation. QTL for shoot P content and PUE were located on chromosomes 3 and 4, respectively (Hammond et al., 2009). Of the five MTAs on Ca6, two were in the same gene (Ca_10411), encoding dynamin-related protein 1C. The MTA identified in the FarmCPU model, Ca1_16163105 on Ca1, was present in the Ca_06938 gene that encodes organic cation/carnitine transporter 4-like (OCT), which is involved in homeostasis in animals and has been well-studied in Arabidopsis thaliana—disruption of AtOCT1 affects root development (Lelandais-Brière et al., 2007). Four major pathways in Lupinus albus that contribute to PUE are carbon fixation, cluster-root formation, soil P mobilization, and cellular P reuse (Xu et al., 2020). A recent effort to understand the genetic basis of photosynthesis and PUE affecting yield reported that three major QTL (q14-2, q15-2, and q19-2) explained 6.6–58.9% of the phenotypic variation (Li et al., 2016). Furthermore, the gene that encodes purple acid phosphatase within the q19-2 region (Glyma.19G193900) is a potential candidate for regulating both soybean PUE and photosynthetic capacity. A total of 159 QTNs within 31 genomic regions and genes associated with photosynthesis-related traits under P stress conditions were genotyped in 2,019 soybean accessions using 292,035 high-quality SNPs and phenotyped under adequate- and low-P conditions for 2 years (Lü et al., 2018). The MTA (Ca1_12310101) for physiological PUE is present in a gene that encodes fasciclin-like arabinogalactan protein 12 (Ca_02521); in Arabidopsis, this gene is involved in cell wall biogenesis (Figure 2C). Further, this SNP locus was associated with total shoot P content and shoot DW in different models (Supplementary Table 1 and Figure 2). Such shared associations were also reported using different models in mungbean PUE (Reddy et al., 2020). A recent study reported that low expression of selected fasciclin-like arabinogalactan protein genes led to kernel abortion in maize (Zea mays) and Arabidopsis thaliana seeds (Cagnola et al., 2018). In cereals such as rice, SNP loci on chromosomes 1, 4, 11, and 12 are associated with PUE (Wissuwa et al., 2015).

Pang et al. (2018a) reported that root foraging and root physiology, such as the exudation of carboxylates into the rhizosphere, are important strategies for plant P acquisition efficiency. A positive correlation was also identified between mature leaf Mn concentration and rhizosheath carboxylate amount relative to root DW, and hence the carboxylate-releasing P-mobilizing strategy was proxied by foliar Mn concentration in a large set of chickpea germplasm under low P supply (Pang et al., 2018a). The MTA for Mn concentration in mature leaves identified on Ca4 was in gene Ca_14893 that encodes zinc finger BED domain-containing protein RICESLEEPER 2-like protein. Carboxylate exudation is an important physiological root trait that enables plants to mine soil P (Lambers et al., 2008, 2015; Richardson and Simpson, 2011). Carboxylate concentrations in the rhizosheath are positively correlated with shoot P content (Pang et al., 2018a). A SNP locus associated with carboxylate amount was identified in the rhizosheath (μmol g^–1 root DW) associated with a Ca_11019 gene that encodes for ABC1 family protein, and is an integral part of the membrane. The FarmCPU model identified two significant MTAs (Ca4_15651907 and Ca8_230753) in RhizoPH.

Conclusion

In summary, the SNP loci associated with more than one trait were identified. For instance, Ca1_12310101 on Ca1 is associated with three traits (i.e., physiological PUE, shoot DW, and shoot P content), Ca2_31290805 on Ca2 is associated with P utilization and total root length, Ca4_37796452 on Ca4 is associated with RMR and shoot P content, Ca3_6798755 on Ca3 is associated with RGR and TRLz1, and Ca7_5414752 on Ca7 is associated with leaf intracellular CO₂ concentration and WUE. The MTA for Mn in mature leaves identified on Ca4 was in gene Ca_14893 that encodes zinc finger BED domain-containing protein RICESLEEPER 2-like protein. The MTAs reported in this study can be used in chickpea breeding programs to enhance PUE.

Data Availability Statement

Publicly available datasets were analyzed in this study. This data can be found here: CNSA (https://db.cngb.org/cnsa/) of CNGBdb with accession code CNP0000370.

Author Contributions

KHMS and RKV conceived the study. MT and DK performed GWAS. YC and JP generated the phenotyping data. PB, MRo, AC, MRy, and HL contributed to the resources and the writing. MT, RKV, and KHMS prepared the manuscripts. All authors read and approved the final manuscript.

Funding

RKV greatly acknowledges financial support from the Department of Agriculture and Cooperation and Farmers Welfare, Government of India, UWA Research Support from the University of Western Australia, Australia and the Bill and Melinda Gates Foundation, United States.

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.

Supplementary Material

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

Supplementary Figure 1 | Example of spurious associations detected for taproot length (TRL, mm) in the GLM model, with no associations detected in the MLM and FarmCPU models.

Supplementary Figure 2 | Population structure. (A,B) Three sub-populations among 233 genotypes phenotyped for root and phosphorus-related traits, and (C,D) three sub-populations among 91 genotypes phenotyped for shoot/root morphological and physiological traits. Each colored vertical line represents proportions of ancestral populations (K) for each individual. Optimum K value determined using ADMIXTURE’s cross-validation procedure.

Supplementary Table 1 | Details of significant MTAs with different traits and functional annotation of SNP loci.

Footnotes

1. ^ https://www.uniprot.org/uniprot/Q8RXC8

References

Assefa, T., Zhang, J., Chowda-Reddy, R. V., Lauter, A. N. M., Singh, A., O’Rouke, J. A., et al. (2020). Deconstructing the genetic architecture of iron deficiency chlorosis in soybean using genome-wide approaches. BMC Plant Biol. 20:42. doi: 10.1186/s12870-020-2237-5

PubMed Abstract | CrossRef Full Text | Google Scholar

Bharadwaj, C., Tripathi, S., Soren, K. R., Thudi, M., Singh, R. K., Sheoran, S., et al. (2020). Introgression of the “QTL-hotspot” region enhances drought tolerance and grain yield in multiple genetic backgrounds of elite chickpea cultivars in India. Plant Genome 14:e20076. doi: 10.1002/tpg2.20076

PubMed Abstract | CrossRef Full Text | Google Scholar

Browning, S. R., and Browning, B. L. (2016). Genotype imputation with millions of reference samples. Am. J. Hum. Genet. 98, 116–126. doi: 10.1016/j.ajhg.2015.11.020

PubMed Abstract | CrossRef Full Text | Google Scholar

Cagnola, J. I., Dumont de Chassart, G. J., Ibarra, S. E., Chimenti, C., Ricardi, M. M., Delzer, B., et al. (2018). Reduced expression of selected fasciclin-like arabinogalactan protein genes associates with the abortion of kernels in field crops of Zea mays (maize) and of Arabidopsis seeds. Plant Cell Environ. 41, 661–674. doi: 10.1111/pce.13136

PubMed Abstract | CrossRef Full Text | Google Scholar

Chaudhary, R. S., Somasundaram, J., Mandal, K. G., and Hati, K. M. (2018). Enhancing water and phosphorus use efficiency through moisture conservation practices and optimum phosphorus application in rainfed maize–chickpea system in vertisols of Central India. Agric. Res. 7, 176–186. doi: 10.1007/s40003-018-0316-8

CrossRef Full Text | Google Scholar

Chen, X., Wang, H., Li, X., Ma, K., Zhan, Y., and Zeng, F. (2019). Molecular cloning and functional analysis of 4-Coumarate:CoA ligase 4 (4CL-like 1) from Fraxinus mandshurica and its role in abiotic stress tolerance and cell wall synthesis. BMC Plant Biol. 19:231. doi: 10.1186/s12870-019-1812-0

PubMed Abstract | CrossRef Full Text | Google Scholar

Chen, Y., Ghanem, M. E., and Siddique, K. H. M. (2017). Characterising root trait variability in chickpea (Cicer arietinum L.) germplasm. J. Exp. Bot. 68, 1987–1999.

Google Scholar

Cong, W.-F., Suriyagoda, L. D. B., and Lambers, H. (2020). Tightening the phosphorus cycle through phosphorus-efficient crop genotypes. Trends Plant Sci. 25, 967–975. doi: 10.1016/j.tplants.2020.04.013

PubMed Abstract | CrossRef Full Text | Google Scholar

Cosgrove, D. J. (2016). Plant cell wall extensibility: connecting plant cell growth with cell wall structure, mechanics, and the action of wall-modifying enzymes. J. Exp. Bot. 67, 463–476. doi: 10.1093/jxb/erv511

PubMed Abstract | CrossRef Full Text | Google Scholar

Diapari, M., Sindhu, A., Bett, K., Deokar, A., Warkentin, T. D., and Tar’an, B. (2014). Genetic diversity and association mapping of iron and zinc concentrations in chickpea (Cicer arietinum L.). Genome 57, 459–468. doi: 10.1139/gen-2014-0108

PubMed Abstract | CrossRef Full Text | Google Scholar

Farahani, S., Maleki, M., Mehrabi, R., Kanouni, H., Scheben, A., Batley, J., et al. (2019). Whole genome diversity, population structure, and linkage disequilibrium analysis of chickpea (Cicer arietinum L.) genotypes using genome-wide DArTseq-based SNP markers. Genes 10:676. doi: 10.3390/genes10090676

PubMed Abstract | CrossRef Full Text | Google Scholar

Gali, K. K., Sackville, A., Tafesse, E. G., Lachagari, V. B. R., McPhee, K., Hybl, M., et al. (2019). Genome-wide association mapping for agronomic and seed quality traits of field pea (Pisum sativum L.). Front. Plant Sci. 10:1538. doi: 10.3389/fpls.2019.01538

PubMed Abstract | CrossRef Full Text | Google Scholar

Giovannetti, M., Göschl, C., Dietzen, C., Andersen, S. U., Kopriva, S., and Busch, W. (2019). Identification of novel genes involved in phosphate accumulation in Lotus japonicus through genome wide association mapping of root system architecture and anion content. PLoS Genet 15:e1008126. doi: 10.1371/journal.pgen.1008126

PubMed Abstract | CrossRef Full Text | Google Scholar

Hammond, J. P., Broadley, M. R., White, P. J., King, G. J., Bowen, H. C., Hayden, R., et al. (2009). Shoot yield drives phosphorus use efficiency in Brassica oleracea and correlates with root architecture traits. J. Exp. Bot. 60, 1953–1968. doi: 10.1093/jxb/erp083

PubMed Abstract | CrossRef Full Text | Google Scholar

He, J., Jin, Y., Du, Y.-L., Wang, T., Turner, N. C., Yang, R.-P., et al. (2017). Genotypic variation in yield, yield components, root morphology and architecture, in soybean in relation to water and phosphorus supply. Front. Plant Sci. 8:1499. doi: 10.3389/fpls.2017.01499

PubMed Abstract | CrossRef Full Text | Google Scholar

Hiremath, P. J., Kumar, A., Penmetsa, R. V., Farmer, A., Schlueter, J. A., Chamarthi, S. K., et al. (2012). Large-scale development of cost-effective SNP marker assays for diversity assessment and genetic mapping in chickpea and comparative mapping in legumes. Plant Biotechnol. J. 10, 716–732. doi: 10.1111/j.1467-7652.2012.00710.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Jadhav, A. A., Rayate, S. J., Mhase, L. B., Thudi, M., Chitikineni, A., Harer, P. N., et al. (2015). Marker-trait association study for protein content in chickpea (Cicer arietinum L.). J. Genet. 94, 279–286. doi: 10.1007/s12041-015-0529-6

PubMed Abstract | CrossRef Full Text | Google Scholar

Jaganathan, D., Thudi, M., Kale, S., Azam, S., Roorkiwal, R., Gaur, P. M., et al. (2015). Genotyping-by-sequencing based intra-specific genetic map refines a “QTL-hotspot” region for drought tolerance in chickpea. Mol. Genet. Genom. 290, 559–571. doi: 10.1007/s00438-014-0932-3

PubMed Abstract | CrossRef Full Text | Google Scholar

Kaashyap, M., Ford, R., Kudapa, H., Jain, M., Edwards, D., Varshney, R. K., et al. (2018). Differential regulation of genes involved in root morphogenesis and cell wall modification is associated with salinity tolerance in chickpea. Sci. Rep. 8:4855.

Google Scholar

Kale, S. M., Jaganathan, D., Ruperao, P., Chen, C., Punna, R., Kudapa, H., et al. (2015). Prioritization of candidate genes in “QTL-hotspot” region for drought tolerance in chickpea (Cicer arietinum L.). Sci. Rep. 5:15296.

Google Scholar

Kaler, A. S., Gillman, J. D., Beissinger, T., and Purcell, L. C. (2020). Comparing different statistical models and multiple testing corrections for association mapping in soybean and maize. Front. Plant Sci. 10:1794. doi: 10.3389/fpls.2019.01794

PubMed Abstract | CrossRef Full Text | Google Scholar

Katuuramu, D. N., Hart, J. P., Porch, T. G., Grusak, M. A., Glahn, R. P., and Cichy, K. A. (2018). Genome−wide association analysis of nutritional composition−related traits and iron bioavailability in cooked dry beans (Phaseolus vulgaris L.). Mol. Breed. 38:44.

Google Scholar

Khazaei, H., Podder, R., Caron, C. T., Kundu, S. S., Diapari, M., Vandenberg, A., et al. (2017). Marker–trait association analysis of iron and zinc concentration in lentil (Lens culinaris Medik.) seeds. Plant Genome 10:lantgenome2017.02.0007. doi: 10.3835/plantgenome2017.02.0007

PubMed Abstract | CrossRef Full Text | Google Scholar

Kohli, P. S., Kumar Verma, P., Verma, R., Parida, S. K., Thakur, J. K., and Giri, J. (2020). Genome-wide association study for phosphate deficiency responsive root hair elongation in chickpea. Funct. Integr. Genomics 20, 775–786. doi: 10.1007/s10142-020-00749-6

PubMed Abstract | CrossRef Full Text | Google Scholar

Lambers, H., Hayes, P. E., Laliberté, E., Oliveira, R. S., and Turner, B. L. (2015). Leaf manganese accumulation and phosphorus-acquisition efficiency. Trends Plant Sci. 20, 83–90. doi: 10.1016/j.tplants.2014.10.007

PubMed Abstract | CrossRef Full Text | Google Scholar

Lambers, H., Raven, J. A., Shaver, G. R., and Smith, S. E. (2008). Plant nutrient-acquisition strategies change with soil age. Trends Ecol. Evol. 23, 95–103. doi: 10.1016/j.tree.2007.10.008

PubMed Abstract | CrossRef Full Text | Google Scholar

Lelandais-Brière, C., Jovanovic, M., Torres, G. A., Perrin, Y., Lemoine, R., Corre-Menguy, F., et al. (2007). Disruption of AtOCT1, an organic cation transporter gene, affects root development and carnitine-related responses in Arabidopsis. Plant J. 51, 154–164. doi: 10.1111/j.1365-313x.2007.03131.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, B.-B., Wang, X., Tai, L., Ma, T.-T., Shalmani, A., Liu, W.-T., et al. (2018). NAD kinases: metabolic targets controlling redox co-enzymes and reducing power partitioning in plant stress and development. Front. Plant Sci. 9:379. doi: 10.3389/fpls.2018.00379

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, H., Yang, Y., Zhang, H., Chu, S., Zhang, X., Yin, D., et al. (2016). A genetic relationship between phosphorus efficiency and photosynthetic traits in soybean as revealed by QTL analysis using a high-density genetic map. Front. Plant Sci. 7:924. doi: 10.3389/fpls.2016.00924

PubMed Abstract | 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. Bioinformatics 28, 2397–2399. doi: 10.1093/bioinformatics/bts444

PubMed Abstract | CrossRef Full Text | Google Scholar

Lü, H., Yang, Y., Li, H., Liu, Q., Zhang, J., Yin, J., et al. (2018). Genome-wide association studies of photosynthetic traits related to phosphorus efficiency in soybean. Front. Plant Sci. 9:1226. doi: 10.3389/fpls.2018.01226

PubMed Abstract | CrossRef Full Text | Google Scholar

MacDonald, G. K., Bennett, E. M., Potter, P. A., and Ramankutty, N. (2011). Agronomic phosphorus imbalances across the world’s croplands. Proc. Natl. Acad. Sci. 108, 3086–3091. doi: 10.1073/pnas.1010808108

PubMed Abstract | CrossRef Full Text | Google Scholar

Mamidi, S., Lee, R. K., Goos, J. R., and McClean, P. E. (2014). Genome−wide association studies identify seven major regions responsible for iron deficiency chlorosis in soybean (Glycine max). PLoS One 9:e107469. doi: 10.1371/journal.pone.0107469

PubMed Abstract | CrossRef Full Text | Google Scholar

Mannur, D. M., Babbar, A., Thudi, M., Sabbavarapu, M. M., Roorkiwal, M., Sharanabasappa, Y., et al. (2019). Super Annigeri 1 and improved JG 74: two Fusarium wilt-resistant introgression lines developed using marker-assisted backcrossing approach in chickpea (Cicer arietinum L.). Mol. Breed. 39:2.

Google Scholar

Mogollón, J. M., Beusen, A. H. W., van Grinsven, H. J. M., Westhoek, H., and Bouwman, A. F. (2018). Future agricultural phosphorus demand according to the shared socioeconomic pathways. Global Env. Change 50, 149–163. doi: 10.1016/j.gloenvcha.2018.03.007

CrossRef Full Text | Google Scholar

Motsara, M. R. (2002). Fertility status of Indian soils. Fertilizer News 47, 15–22.

Google Scholar

Nayak, S. N., Zhu, H., Varghese, N., Datta, S., Choi, H.-K., Horres, R., et al. (2010). Integration of novel SSR and gene-based SNP marker loci in the chickpea genetic map and establishment of new anchor points with Medicago truncatula genome. Theor. Appl. Genet. 120, 1415–1441. doi: 10.1007/s00122-010-1265-1

PubMed Abstract | CrossRef Full Text | Google Scholar

Ozkuru, E., Ates, D., Nemli, S., Erdogmus, S., Karaca, N., Yilmaz, H., et al. (2018). Association mapping of loci linked to copper, phosphorus, and potassium concentrations in the seeds of C. arietinum and C. reticulatum. Genomics 111, 1873–1881. doi: 10.1016/j.ygeno.2018.12.010

PubMed Abstract | CrossRef Full Text | Google Scholar

Pang, J., Zhao, H., Bansal, R., Bohuon, E., Lambers, H., Ryan, M. H., et al. (2018c). Leaf transpiration plays a role in phosphorus acquisition among a large set of chickpea genotypes. Plant Cell Environ. 41, 2069–2079.

Google Scholar

Pang, J., Bansal, R., Zhao, H., Bohuon, E., Lambers, H., Ryan, M. H., et al. (2018a). The carboxylate-releasing phosphorus-mobilizing strategy can be proxied by foliar manganese concentration in a large set of chickpea germplasm under low phosphorus supply. New Phytologist 219, 518–529. doi: 10.1111/nph.15200

PubMed Abstract | CrossRef Full Text | Google Scholar

Pang, J., Ryan, M. H., Lambers, H., and Sidddique, K. H. M. (2018b). Phosphorus acquisition and utilisation in crop legumes under global change. Curr. Opin. Plant Biol. 45, 248–254. doi: 10.1016/j.pbi.2018.05.012

PubMed Abstract | CrossRef Full Text | Google Scholar

Purushothaman, R., Krishnamurthy, L., Upadhyaya, H. D., Vadez, V., and Varshney, R. K. (2017). Root traits confer grain yield advantages under terminal drought in chickpea (Cicer arietinum L.). Field Crops Res. 201, 146–161. doi: 10.1016/j.fcr.2016.11.004

PubMed Abstract | CrossRef Full Text | Google Scholar

Reddy, V. R. P., Das, S., Dikshit, H. K., Mishra, G. P., Aski, M., Meena, S. K., et al. (2020). Genome-wide association analysis for phosphorus use efficiency traits in mungbean (Vigna radiata L. Wilczek) using genotyping by sequencing approach. Front. Plant Sci. 11:537766. doi: 10.3389/fpls.2020.537766

PubMed Abstract | CrossRef Full Text | Google Scholar

Richardson, A. E., and Simpson, R. J. (2011). Soil microorganisms mediating phosphorus availability. Plant Physiol. 156, 989–996. doi: 10.1104/pp.111.175448

PubMed Abstract | CrossRef Full Text | Google Scholar

Roorkiwal, M., Bharadwaj, C., Barmukh, R., Dixit, G. P., Thudi, M., Gaur, P. M., et al. (2020). Integrating genomics for chickpea improvement: achievements and opportunities. Theor. Appl. Genet. 133, 1703–1720. doi: 10.1007/s00122-020-03584-2

PubMed Abstract | CrossRef Full Text | Google Scholar

Sab, S., Lokesha, R., Mannur, D. M., Somasekhar, Jadhav, K., Mallikarjuna, B. P., et al. (2020). Genome-wide SNP discovery and mapping QTLs for seed iron and zinc concentrations in chickpea (Cicer arietinum L.). Front. Nutr. 7:559120. doi: 10.3389/fnut.2020.559120

PubMed Abstract | CrossRef Full Text | Google Scholar

Sabbavarapu, M. M., Sharma, M., Chamarthi, S. K., Swapna, N., Rathore, A., Thudi, M., et al. (2013). Molecular mapping of QTLs for resistance to Fusarium wilt (race 1) and Ascochyta blight in chickpea (Cicer arietinum L.). Euphytica 193, 121–133. doi: 10.1007/s10681-013-0959-2

CrossRef Full Text | Google Scholar

Sani, S. G. A. S., Chang, P. L., Zubair, A., Carrasquilla-Garcia, N., Cordeiro, M., Penmetsa, R. V., et al. (2017). Genetic diversity, population structure, and genetic correlation with climatic variation in chickpea (Cicer arietinum) landraces from Pakistan. Plant Genome 11:170067. doi: 10.3835/plantgenome2017.08.0067

PubMed Abstract | CrossRef Full Text | Google Scholar

Syers, J. K., Johnston, A. A., and Curtin, D. (2008). Efficiency of Soil and Fertilizer Phosphorus Use. Reconciling Changing Concepts of Soil Phosphorus Behaviour with Agronomic Information. Rome: Food and Agriculture Organization of the United Nations (FAO), 110.

Google Scholar

Thudi, M., Bohra, A., Nayak, S. N., Varghese, N., Shah, T. M., Penmetsa, R. V., et al. (2011). Novel SSR markers from BAC-end sequences, DArT arrays and a comprehensive genetic map with 1,291 marker loci for chickpea (Cicer arietinum L.). PLoS One 6:e27275. doi: 10.1371/journal.pone.0027275

PubMed Abstract | CrossRef Full Text | Google Scholar

Thudi, M., Khan, A. W., Kumar, V., Gaur, P. M., Katta, A. V. S. K., Garg, V., et al. (2016b). Whole genome re-sequencing reveals genome wide variations among parental lines of mapping populations in chickpea (Cicer arietinum). BMC Plant Biol. 16:10. doi: 10.1186/s12870-015-0690-3

PubMed Abstract | CrossRef Full Text | Google Scholar

Thudi, M., Chitikineni, A., Liu, X., He, W., Roorkiwal, M., Yang, W., et al. (2016a). Recent breeding programs enhanced genetic diversity in both desi and kabuli varieties of chickpea (Cicer arietinum L.). Sci. Rep. 6:38636.

Google Scholar

Thudi, M., Upadhyaya, H. D., Rathore, A., Gaur, P. M., Krishnamurthy, L., Roorkiwal, M., et al. (2014). Genetic dissection of drought and heat tolerance in chickpea through genome-wide and candidate gene-based association mapping approaches. PLoS One 9:e96758. doi: 10.1371/journal.pone.0096758

PubMed Abstract | CrossRef Full Text | Google Scholar

Upadhyaya, H. D., Bajaj, D., Narnoliya, L., Das, S., Kumar, V., Gowda, C. L. L., et al. (2016). Genome-wide scans for delineation of candidate genes regulating seed-protein content in chickpea. Front. Plant Sci. 7:302. doi: 10.3389/fpls.2016.00302

PubMed Abstract | CrossRef Full Text | Google Scholar

Upadhyaya, H. D., Dwivedi, S. L., Baum, M., Varshney, R. K., Udupa, S. M., Gowda, C. L. L., et al. (2008). Genetic structure, diversity, and allelic richness in composite collection and reference set in chickpea (Cicer arietinumL.). BMC Plant Biol. 8:106. doi: 10.1186/1471-2229-8-106

PubMed Abstract | CrossRef Full Text | Google Scholar

Varshney, R. K., Song, C., Saxena, R. K., Azam, S., Yu, S., Sharpe, A. G., et al. (2013b). Draft genome sequence of chickpea (Cicer arietinum) provides a resource for trait improvement. Nat. Biotechnol. 31, 240–246.

Google Scholar

Varshney, R. K., Gaur, P. M., Chamarthi, S. K., Krishnamurthy, L., Tripathi, S., Kashiwagi, J., et al. (2013a). Fast-track introgression of “QTL-hotspot” for root traits and other drought tolerance traits in JG 11, an elite and leading variety of chickpea. Plant Genome 6:lantgenome2013.07.0022. doi: 10.3835/plantgenome2013.07.0022

CrossRef Full Text | Google Scholar

Varshney, R. K., Thudi, M., Nayak, S. N., Gaur, P. M., Kashiwagi, J., Krishnamurthy, L., et al. (2014). Genetic dissection of drought tolerance in chickpea (Cicer arietinum L.). Theor. Appl. Genet. 127, 445–462.

Google Scholar

Varshney, R. K., Thudi, M., Roorkiwal, M., He, W., Upadhyaya, H. D., Yang, W., et al. (2019). Resequencing of 429 chickpea accessions from 45 countries provides insights into genome diversity, domestication and agronomic traits. Nat. Genet. 51, 857–864. doi: 10.1038/s41588-019-0401-3

PubMed Abstract | CrossRef Full Text | Google Scholar

Wissuwa, M., Kondo, K., Fukuda, T., Mori, A., Rose, M. T., Pariasca-Tanaka, J., et al. (2015). Unmasking novel loci for internal phosphorus utilization efficiency in rice germplasm through genome-wide association analysis. PLoS One 10:e0124215. doi: 10.1371/journal.pone.0124215

PubMed Abstract | CrossRef Full Text | Google Scholar

Xiao, S., Gao, W., Chen, Q. F., Chan, S. W., Zheng, S. X., Ma, J., et al. (2010). Overexpression of Arabidopsis acyl-CoA binding protein ACBP3 promotes starvation-induced and age-dependent leaf senescence. Plant Cell 22, 1463–1482. doi: 10.1105/tpc.110.075333

CrossRef Full Text | Google Scholar

Xu, W., Zhang, Q., Yuan, W., Xu, F., Aslam, M. M., Miao, R., et al. (2020). The genome evolution and low-phosphorus adaptation in white lupin. Nat. Commun. 11:1069.

Google Scholar

Zhou, H., Alexander, D., and Lange, K. (2011). A quasi-Newton acceleration for high-dimensional optimization algorithms. Stat. Comput. 21, 261–273. doi: 10.1007/s11222-009-9166-3

PubMed Abstract | CrossRef Full Text | Google Scholar