Genetic Dissection and Identification of Candidate Genes for Salinity Tolerance Using Axiom®CicerSNP Array in Chickpea

Globally, chickpea production is severely affected by salinity stress. Understanding the genetic basis for salinity tolerance is important to develop salinity tolerant chickpeas. A recombinant inbred line (RIL) population developed using parental lines ICCV 10 (salt-tolerant) and DCP 92-3 (salt-sensitive) was screened under field conditions to collect information on agronomy, yield components, and stress tolerance indices. Genotyping data generated using Axiom®CicerSNP array was used to construct a linkage map comprising 1856 SNP markers spanning a distance of 1106.3 cM across eight chickpea chromosomes. Extensive analysis of the phenotyping and genotyping data identified 28 quantitative trait loci (QTLs) explaining up to 28.40% of the phenotypic variance in the population. We identified QTL clusters on CaLG03 and CaLG06, each harboring major QTLs for yield and yield component traits under salinity stress. The main-effect QTLs identified in these two clusters were associated with key genes such as calcium-dependent protein kinases, histidine kinases, cation proton antiporter, and WRKY and MYB transcription factors, which are known to impart salinity stress tolerance in crop plants. Molecular markers/genes associated with these major QTLs, after validation, will be useful to undertake marker-assisted breeding for developing better varieties with salinity tolerance.

as a better strategy to identify QTLs and develop crop varieties for enhanced stress tolerance [41]. Such relative performance has practical relevance since genotypes with low yield potential under control conditions often show higher tolerance to stress than high-yielding genotypes. Therefore, it is important to identify QTLs governing salinity stress tolerance for yield and yield components that can be used effectively in marker-assisted breeding in chickpea.
With an objective to dissect the genetic basis of salt tolerance in chickpea, the present study has been undertaken with following objectives: (a) develop a dense genetic map for ICCV 10 (salt stress-tolerant) × DCP 92-3 (salt stress-sensitive) population, (b) identify QTLs for salinity stress tolerance in chickpea, and (c) identify genes underlying major QTLs associated with salinity tolerance component traits in chickpea.

Phenotypic Variation for Salinity Tolerance Component Traits
RILs, along with parents, were analyzed for their phenotypic performance in control and salt-stressed environments. Salt stress significantly reduced yield in the salinity-sensitive parent (DCP 92-3), being seven-fold more than the salinity-tolerant parent (ICCV 10) across seasons (Table 1). During the 2015-16 and 2016-17 crop seasons, a normal frequency distribution was observed for all traits analyzed, except for stress tolerance index (STI) and stress susceptibility index (SSI) for yield per plant (Supplementary Figures S1 and S2).
The salinity tolerance component traits were subjected to Pearson's correlation analysis to identify relationships between them. For simplicity, we classified significant correlation values into three different degrees, viz. high, moderate, and low. For instance, a high degree of correlation represented correlation coefficient values between ±0.50 and ±1.0 and was considered to be a strong correlation. A moderate degree of correlation represented coefficient values between ±0.30 and ±0.49, while coefficient values below ±0.29 were classified as low degree of correlation. In this context, we evaluated the correlations between traits that were measured during the 2015-16 season ( Table 2, Supplementary Figure S3a). A strong positive correlation was observed between SSI for yield per plant and yield per plant evaluated under control conditions, and between STI for yield per plant and yield per plant under stress conditions. In contrast, SSI for yield per plant and yield per plant measured under salinity stress showed a strong negative correlation. Furthermore, STI for yield per plant displayed a medium positive correlation with pods per plant, 100-seed weight and plant height, all evaluated under salinity stress conditions. STI for 100-seed weight showed a small positive correlation with STI for yield per plant and yield per plant assessed under stress conditions.
We also assessed the correlation between traits evaluated during the 2016-17 crop season (Table 2, Supplementary Figure S3b). Here, STI for yield per plant displayed a strong positive correlation with yield per plant measured under salinity stress. Although SSI for yield per plant was strongly and positively related to yield per plant measured under control conditions, it showed a strong negative correlation with yield per plant measured under stress conditions. Furthermore, STI for yield per plant showed a positive and moderate correlation with pods per plant and yield per plant, both measured under control scenarios. The SSI for yield per plant was positively and moderately correlated with pods per plant and plant height under control conditions. Furthermore, a small and positive correlation was observed between STI for yield per plant and plant height measured under control, while SSI for yield per plant showed a small and negative correlation with 100-seed weight measured under stress conditions.

High-Density Linkage Map
Of the 50,591 SNP probes on the Axiom ® CicerSNP array, 5123 SNPs were polymorphic between both parents and displayed segregation within the population. In the constructed linkage map, 1856 markers were assigned to eight linkage groups. The integrated map had a total length of 1106.3 cM with an average distance of 0.59 cM between adjacent markers. The number of loci onto eight linkage groups ranged from 56 (CaLG04) to 487 (CaLG07) and the length ranged from 45.6 cM (CaLG08) to 270.7 cM (CaLG06) with a mean value of 138.2 cM (Table 3). The marker density of the linkage groups varied significantly, with the highest marker density on CaLG01 with 327 markers distributed over a distance of 147.7 cM, while CaLG04 had the lowest density with 56 markers distributed over a distance of 87 cM (Table 3). Table 3. Distribution of markers on the eight linkage groups (LGs) of the chickpea genetic map for the ICCV 10 × DCP 92-3 RIL population.

S. No.
LGs

QTLs for Salinity Tolerance Component Traits
The phenotypic and genotypic data were analyzed for identification of QTLs to understand the genetic basis of salinity tolerance in the RIL mapping population derived from ICCV 10 × DCP 92-3. The QTL analysis was performed for seven yield and yield-related traits, namely, pod number per plant (PPP), 100 seed weight (100SW), yield per plant (YPP), SSI for yield per plant (SSI_YP), SSI for 100 seed weight (SSI_100SW), STI for yield per plant (STI_YP), and STI for 100 seed weight (STI_100SW), and two agronomic traits, i.e., plant height (PH) and branch number per plant (NB). The QTLs that contributed >10% of the phenotypic variation explained (PVE) were considered as major QTLs. Furthermore, if QTL for a given trait in a particular treatment appeared in both years, it was considered as a consistent QTL [32].

Pod number per plant (PPP):
In the case of PPP, although no major QTL was identified, one minor QTL (qPPP8.1; PVE 8.1%) was identified on CaLG08 in the salinity treatment in 2015-16, and one minor QTL with 9.1% PVE was identified on CaLG07 in the control (Table 4).  Table 4).

Number of branches per plant (NB):
One major QTL (qNBC8.1; PVE 12.7) was identified for NB on CaLG08 flanked by AX-123638389 and AX-123638445 markers in 2015-16 under control treatment. This QTL was considered as a consistent QTL as it also appeared in 2016-17 season (Table 4). Two minor QTLs (qNBC8.2; PVE 8.9% and qNBC8.3; PVE 6.1%) were also identified for NB under control condition on CaLG08 in 2015-16 and 2016-17, respectively. No QTL was identified for NB under the salinity treatment across the seasons (Table 4).

Candidate Genes for Salinity Tolerances
The genes located in the genomic regions for the identified QTLs were extracted. A total of 1121 genes were present in all QTL regions, based on the reference chickpea genome [1]. Of the 1121 genes, 136 putative candidate genes (88 on CaLG03, 25 on CaLG06, 14 on CaLG02, six on CaLG07, and one each on CaLG04 and CaLG08) belong to several gene families, including kinases, genes encoding transcription factors, and ion channels encoding proteins in response to salinity stress, among others. These genomic regions on CaLG03 and CaLG06 are of great interest as they hold QTLs for traits related to yield under salinity (Supplementary Table S1, Figure 1, and Supplementary Figure S4). The major QTL intervals AX-123622602 and AX-123622699 on CaLG03 and AX-123640392 and AX-123655575 on CaLG06 harbored 113 key genes that are reportedly involved in salinity and other abiotic stress tolerances in chickpea and other crops (Supplementary Table S2). protein, which reportedly improve and enhances salt tolerance in various crops, was identified 235 [45,46]. An abscisic acid insensitive (ABI) 5, a typical subfamily of proteins belonging to the basic 236 domain/leucine zipper (bZIP) transcription factors, was recently reported to contribute to salt stress 237 tolerance in Arabidopsis thaliana via abscisic acid signaling [47]. Notably, the ABI5 gene was also found 238 in the predicted QTL region. The stress-responsive NAC transcription factor family of proteins has 239 been studied extensively for its role in salinity stress tolerance in transgenic plants [48,49]. We 240 identified the presence of the NAC-domain-containing protein in the QTL region. In addition to the 241 above genes, some trait-specific candidate genes, such as E3-ubiquitin-protein ligase, WRKY 242 transcription factor, and zinc finger proteins, were also identified. The gene encoding E3 ubiquitin 243 ligase activity has been shown to regulate grain width and weight in rice [50]. Furthermore, while the 244 differential expression of the WRKY transcription factor gene alters seed size in soybean [51], a zinc 245 finger CCHC-type protein is involved in regulating seed size in Medicago truncatula [52]. However, 246 further fine mapping and validation of these genes would be needed to pinpoint the candidate genes 247 regulating traits of interest.

248
Analysis of the QTLs from the CaLG06 genomic region revealed that QTLs for three traits, been demonstrated to play important roles in conferring salinity tolerance in Chenopodium quinoa [55].
In this region, the presence of transmembrane domain proteins is anticipated to be involved in A detailed analysis of QTLs from CaLG03 showed that QTLs for four traits (STI_100SW, SSI_100SW, 100SW, SSI_YP) were flanked by AX-123622602 and AX-123622699 markers. Here, the~3.3 Mb region between these two markers was selected to identify candidate genes. The~3.3 Mb region contained 763 predicted genes. Functional annotation of the candidate genes revealed their roles in various biotic and abiotic stress responses. For instance, a MYB transcription factor, which is mainly involved in the regulation of plant growth and development under abiotic stress, such as salinity tolerance, was identified in this region [42][43][44]. Similarly, a C3HC4-type RING zinc finger protein, which reportedly improve and enhances salt tolerance in various crops, was identified [45,46]. An abscisic acid insensitive (ABI) 5, a typical subfamily of proteins belonging to the basic domain/leucine zipper (bZIP) transcription factors, was recently reported to contribute to salt stress tolerance in Arabidopsis thaliana via abscisic acid signaling [47]. Notably, the ABI5 gene was also found in the predicted QTL region. The stress-responsive NAC transcription factor family of proteins has been studied extensively for its role in salinity stress tolerance in transgenic plants [48,49]. We identified the presence of the NAC-domain-containing protein in the QTL region. In addition to the above genes, some trait-specific candidate genes, such as E3-ubiquitin-protein ligase, WRKY transcription factor, and zinc finger proteins, were also identified. The gene encoding E3 ubiquitin ligase activity has been shown to regulate grain width and weight in rice [50]. Furthermore, while the differential expression of the WRKY transcription factor gene alters seed size in soybean [51], a zinc finger CCHC-type protein is involved in regulating seed size in Medicago truncatula [52]. However, further fine mapping and validation of these genes would be needed to pinpoint the candidate genes regulating traits of interest.
Analysis of the QTLs from the CaLG06 genomic region revealed that QTLs for three traits, namely, STI_YP, SSI_YP, and YPP, were flanked by AX-123640392 and AX-123655575 markers spanning a~0.1 Mb region. The predicted~0.1 Mb region encompassed 155 genes. Analysis of the genes underlying the QTL interval identified several potential candidates with a predicted role in salinity stress tolerance. For example, overexpression of cation/proton antiporter genes encoding cellular Na + /H + exchanger proteins can upregulate plant performance under salinity stress [53]. A cation/H + antiporter 4-like gene was identified in this region. Further, a plant calmodulin-binding family protein, predicted to be involved in salinity stress responses in plants during early germination, was also detected [54]. Membrane proteins containing a transmembrane domain have been demonstrated to play important roles in conferring salinity tolerance in Chenopodium quinoa [55]. In this region, the presence of transmembrane domain proteins is anticipated to be involved in conferring salinity tolerance in the genotypes. Notably, histidine kinases act as sensory molecules and are involved in the transduction of environmental signals in plants, fungi, and prokaryotes. A previous study identified the involvement of histidine kinases in perception or transduction of salt stress signals in Synechocystis sp. [56]. The presence of a histidine kinase gene underlying the QTL region is thought to be involved in the perception of salt stress signals and the regulation of salt-inducible genes. We anticipate that further fine mapping and cloning of the candidate genes underlying the major QTL regions will unravel the salinity tolerance mechanism in chickpea.

Discussion
Recent technological advances have led to the evolution of breeding methodologies that have the potential to accelerate the breeding process [57]. Plant response to abiotic stresses, such as salinity, is a complex trait that is governed by many genes; therefore, breeding approaches for such complex stress tolerance and crop stability have been challenging. The use of marker-assisted selection/marker-assisted breeding in this regard has helped to simplify things to a certain extent. Marker-assisted breeding has successfully produced crops with tolerance to biotic stresses and has increased yield in many crops [57][58][59].
Salinity stress appears to be a polygenic and quantitative trait that is regulated by several genes under diverse environments [60]. To understand the genetic basis of salinity stress, significant efforts have been made to identify QTLs associated with salinity tolerance, with several QTLs mapped on different chromosomes of chickpea by multiple research groups in the last decade. In a study by Vadez et al. [36], QTLs for seed weight, pod number, and harvest index in the salinity treatment were obtained on CaLG07 using a RIL population derived from a cross between salt-tolerant JG62 and salt-sensitive ICCV 2. Similarly, a minor QTL for yield that explained 8% PVE on CaLG07 was identified in chickpea [34]. In another study, two key genomic regions on CaL05 and CaL07, harboring QTLs for yield and other salinity-associated traits, were reported in a RIL population derived from the cross ICCV 2 × JG 11, using SSR and SNP markers [35].
Recent advances in genome-based approaches have endorsed the development of high-throughput approaches for genotyping, enabling the identification of and access to desirable alleles, with different QTLs having the potential to affect desired responses. In the current study, using a reasonably large RIL population, high-density genetic map, and phenotyping under controlled conditions, we identified 28 QTLs for nine traits across seasons and treatments in chickpea. Importantly, two genomic regions on CaLG03 and CaLG06, harboring 10 QTLs for salinity stress indices and yield-related traits, were detected.
Identification of QTLs for relative performance of genotypes under stress and controlled conditions has significant over identification of QTLs based on phenotypic performance in the stress environment alone [41]. QTLs regulating salinity stress tolerance by utilizing stress tolerance indices, which differentiate the performance of crops both under control and stress conditions, have been reported in several crops, but not in chickpea. In rice, genomic regions governing heat stress tolerance [40] and salinity stress tolerance have been successfully mapped using stress indices [37][38][39]61]. In the present study, major QTLs for salinity stress indices SSI_YP, SSI_100SW, STI_100SW, and yield component (100SW) under salinity have been identified. These QTLs, after validation, may play a key role in marker-assisted breeding programs to produce saline-tolerant lines in chickpea. Salt-tolerant genotypes will have lower SSI values, indicating the smaller difference in yield between control and stress treatment, with the reverse being true for susceptible genotypes [37,62].
Identifying the genes involved in the salinity stress response is the first step towards gaining the necessary knowledge for genetics, physiology, and breeding. The application of genomics-based knowledge with breeding platforms is expected to provide useful insights into the molecular responses of plants to salinity. Plants alleviate the salinity stress effects of increased Na + in cells by excluding and sequestering Na + . Increased cytosolic Na + concentration is detected and, in turn, the response pathways to stress ion channels, calcium-dependent kinases, histidine kinases, and membrane receptors, etc., activate to regulate gene expression, resulting in the synthesis of compatible solutes to abate the destructive effects of increased Na + ions in the cytoplasm [63]. Plant hormones also play a vital role in regulating plant responses to salinity. Plants induce ABA synthesis to mediate the plant responses through the ABA-dependent signal transduction pathway. To maintain cellular ionic and osmotic homeostasis, ABA induces osmolyte synthesis and stomatal closure in the guard cells, and other related mechanisms leading to the maintenance of cellular ionic and osmotic homeostasis [64][65][66].
The putative candidate genes found in the QTL regions in this study have been experimentally validated for their role in the salinity stress response in several earlier studies in different plant species. For example, genes identified on CaLG03 and CaLG06 that encode several kinases, calcium-dependent protein kinases (CDPKs), mitogen-activated protein kinases (MAPKs), histidine kinases (HKs), sucrose non-fermenting related kinases (SnRK1), transcription factors such as WRKY, basic leucine zipper (bZIP), MYB/MYC, and cation calcium exchanger, were reported to be playing a vital role in the salinity stress response [67][68][69][70]. Similarly, candidate genes coding for proteins related to cation calcium exchanger, cation antiporter 4 (regulates plasma membrane antiporter activity), and potassium channel AKT1 (involved in regulating K + /Na + ratio) led to salinity stress tolerance in Arabidopsis [71,72]. The transmembrane protein-encoding genes found in the QTL region on CaLG03 are candidate genes for seed weight in chickpea [73].
Among the 113 putative candidate genes found in the QTL regions on CaLG03 and CaLG06, most were involved in osmoregulation, which help plants to cope with salinity and other biotic stresses (Table S3). The identification of probable candidate genes for salinity tolerance on small genomic regions on CaLG03 and CaLG06 make these regions promising for future use in genomics-assisted breeding (GAB) to improve salt and other abiotic stress tolerance in chickpea.

Development of RIL Population
Two diverse chickpea parental lines DCP 92-3 and ICCV 10 that differ in salinity tolerance were used to develop a chickpea mapping population comprising of 201 RILs (F 8 ) at the ICAR-Indian Institute of Pulses Research, Kanpur. ICCV 10 is a desi-type genotype and is highly tolerant to salinity stress, while DCP 92-3 is susceptible to salinity. The population was advanced using the single seed descent method.

Phenotyping for Salinity Stress
The RIL population was screened for two consecutive years from November to April (2015- 16 and 2016-17) at Karnal, Haryana, India located between 29 • 43 N longitude and 76 • 58 E latitude. Karnal has an average elevation of 240 m from mean average sea level (m.a.s.l) and received total rainfall of 85.8 mm during January, 7.8 mm in March, and 3.4 mm in April, during the crop season 2015-16. The average temperature during cropping season ranged between 11.92 • C and 26.85 • C with minimum temperature of 6.8 • C in December 2015. The initial humidity during seedling stage was 91.87% (max) and 46.06% (min). During 2016-17 season, the average temperature was between 11.82 • C and 27.36 • C, with a minimum temperature of 6.78 • C in the month of January 2017. The total rainfall received during January was 85.8 mm, 7.8 mm in March, and 3.4 mm in April, 2017. The average humidity of 86% (max) and 11% (min) was observed. Overall, the two cropping seasons had cool environmental conditions suitable for crop growth and expression. The experiment was done under control conditions; hence, it did not affect the level of salinity stress in micro-plot. The mapping population was screened using the micro-plot screening method described below.
A total of 201 RILs, along with the parents and Karnal Chana-1 as a positive check, were sown in micro-plots (2 × 2 m 2 ) top covered with 2 mm polycarbonate sheet to protect it from unseasonal rains. The experiment was organized in an augmented design to evaluate seed yield and other yield contributing traits in control and saline (ECiw 6 dS/m) conditions at ICAR-Central Soil Salinity Research Institute (CSSRI), Karnal. Initial soil salinity ranged from ECe 0.92-0.97 dS/m. Saline irrigation of ECiw 6 dS/m was applied at 30, 60, and 90 day intervals after sowing. Soil samples were taken from time to time to measure the salinity build-up of ECiw 6 dS/m in the soil.
The genotype response to salt stress was expressed as the stress susceptibility index (SSI; [61]) and stress tolerance index (STI; [37]) using the formulae: SSI = (1 − Ys/Yp)/SI;SI = 1 −Ȳs/Ȳp [61] and stress tolerance index STI = Yp × Ys/Ȳp 2 where Ys and Yp are the yield of genotypes evaluated under saline (stress) and non-saline (nonstress) condition, andȲs andȲp are the mean yield of all genotypes evaluated under stress and non-stress conditions.
The SSI is a measure of reduction in yield caused by a stress, relative to a favorable environment. Salt-tolerant genotypes would have a lower SSI value, indicating a smaller yield difference with the control than the susceptible genotypes with a larger yield difference. The STI identifies genotypes that produce higher yields in the control and under stress conditions (higher STI value), which is more desirable than only under stress conditions.

High Density Genotyping
Genomic DNA was collected from pooled young leaves of RILs and both parents using the CTAB method [74]. DNA quality was tested in 1% agarose gel and quantified using NanoDrop 8000 spectrophotometer (Thermo Fisher Scientific Inc. Waltham, MA, USA). DNA concentration was adjusted to a minimum of 50 ng/µL. RILs were genotyped using an Axiom ® CicerSNP array with 50,590 probes distributed on all eight linkage groups as described in Roorkiwal et al. [30].

Genetic Map Construction and QTL Analysis
A total of 5123 polymorphic SNP markers were used to construct a genetic linkage map using Join Map v. 4.1 (https://www.kyazma.nl/index.php/JoinMap) [75] as described in [31]. To find QTLs responsible for salinity tolerance, we used multi-seasonal field phenotyping data for yield and agronomical traits associated with salinity tolerance of parental accessions and 201 RILs from the mapping population derived from DCP 92-3 × ICCV 10. SNP genotyping data, obtained from Axiom ® CicerSNP array-based SNP genotyping, was correlated with the aforementioned phenotyping data using Windows QTL Cartographer ver. 2.5 [76] (composite interval mapping at significant LOD (logarithm of odds) threshold >3.0 at a p < 0.05) as per earlier studies [77,78].

Mining of Candidate Genes
To identify candidate genes, present within the QTL region, flanking markers were subjected to BLAST against chickpea reference genome assembly, and candidate genes in the QTL regions were retrieved.

Statistical Analysis
Pearson's correlation analysis was conducted using the "Hmisc" package in R. Frequency distributions of the component traits within the RIL population were analyzed and plotted with "UsingR" package within the R environment.

Conclusions
The present study, using a DCP 92-3 × ICCV 10 RIL mapping population and 50K SNP genotyping array, identified two genomic regions that harbor QTLs for salinity tolerance in a~3.3 Mb region on CaLG03 and~0.1 Mb region on CaLG06 with major QTLs for yield and salinity tolerance. The genes present in the identified QTL regions in this study reportedly play a key role in salinity tolerance. Further sequencing and functional validation are required to identify the candidate genes in these QTL regions responsible for salt tolerance. Nonetheless, the high-density linkage map, major QTLs, and genomic regions identified in this study can be used in GAB, after suitable markers are developed and validated, to breed high-yielding salinity-tolerant chickpea varieties.
Supplementary Materials: Supplementary materials can be found at http://www.mdpi.com/1422-0067/21/14/ 5058/s1. Table S1: Summary of major and minor QTLs identified on CaLG03 and CaLG06. Table S2: List of genes present in identified QTL regions. Table S3: List of putative candidate genes associated with salinity stress response on CaLG03 and CaLG06. Figure