Association analysis of biotic and abiotic stresses resistance in chickpea (Cicer spp.) using AFLP markers

ABSTRACT

Chickpea suffers from different biotic and abiotic stresses, which has led to considerable decreases in seed yield. Since the economic value of a cultivar is determined by its phenotypic characteristics, good knowledge of the genetic structure and identification of genomic loci associated with desired traits will facilitate crop breeding. Amplified fragment length polymorphism (AFLP) markers associated with different agro-morphological characteristics, cold and drought tolerance as well as Ascochyta blight resistance were identified in chickpea by using 44 genotypes comprising cultigens, landraces, internationally developed improved lines and wild relatives, evaluated in a randomized complete block design with three replications at the Urmia rainfed research station. Using six AFLP primer combinations, 64 clear and reproducible markers were developed. The polymorphism information content values calculated for each primer pair ranged from 0.155 (EcoR1-ACC/Mse1-CAG) to 0.270 (EcoR1-ACC/ Mse1-CTG) with an average of 0.237. Analysis of molecular variance (AMOVA) revealed that 90% of the total variance was due to differences within populations and 10% due to differences among populations. Using a general linear model and applying multiple testing corrections, 24 AFLP markers associated with genes controlling the studied traits were identified. Some identified markers were associated with more than one trait. The identified markers could be of interest in marker-assisted selection in chickpea breeding programmes.

KEYWORDS:

• Chickpea lines
• genetic diversity
• linkage disequilibrium mapping
• molecular markers
• tolerance to stresses

Introduction

Chickpea (Cicer arietinum L.) as a third major cool-season grain legume in the world is grown in 13.9 million hectares with 13.7 million ton production as of 2014 (http://faostat.fao.org). Stresses caused by high or low temperatures lead to crop yield losses due to damage to reproductive organs [1], floral bud development [2], pollen viability and fertilization [1], pod production, seeds per plant, seed yield [3], seed filling [4] and seed composition [5]. Chickpea may be exposed to low temperatures and excess water early in the growing season (especially if it is planted in the fall) or high temperatures and inadequate water at the end of the growing season [6]. Tolerance to low temperatures in chickpea primarily depends on the weather conditions and the sowing date. High levels of cold tolerance are known in local cultivars and related species such as Cicer bijugum, Cicer reticulatum, Cicer echinospermum, Cicer pinnatifidum and Cicer judaicum [7]. Wery [8], by evaluating 29 chickpea genotypes for frost resistance in Montpellier, France, identified three frost-resistant genotypes, FLIP 81-293C, FLIP 82-127C and FLIP 82-128C. A key recurring problem in the screening for abiotic stress tolerance is environmental heterogeneity [6].

Drought, as another major abiotic stress in many parts of the world [9], causes heavy production losses, 30%–100% [10] in cool-season food legumes [6,11]. Stoddard et al. [6] reviewed the screening techniques in cool-season food legumes and listed 25 drought-resistant chickpea genotypes. Nevertheless, these genotypes are far from possessing the desirable characteristics to grow in farmer fields due to being small seeded and susceptibile to Ascochyta blight [12]. The quantitative inheritance of drought tolerance and its interaction with the environment are the major challenges that exist in identifying drought tolerance genotypes in chickpea [13].

Ascochyta blight (causal agent Ascochyta rabiei) is considered to be a major constraint for chickpea production in cool and wet regions of the world [14]. About 50% yield reduction has been reported to have occurred due to Ascochyta blight epidemics in 1979–1982 [15]. Most of the cultivated chickpea varieties are only partially resistant to this pathogen [16]. Pyramiding quantitative trait loci (QTLs) using marker-assisted selection (MAS) is generally considered an appropriate strategy for improving the tolerance levels of chickpea in breeding programmes. Comparing to other crop species, such as cereals, the advances in the development of genomic resources in the leading legume crops of the semi-arid tropics, especially chickpea (Cicer arietinum), have been very slow [17]. However, these programmes have significantly gained by integrating MAS within empirical selection [18].

Association analysis has been used for identification of markers associated with important traits in several plants species including wheat [19], maize [20], sunflower [21,22] and rice [23]. Darvishzadeh [22] by using 106 dispersed oily sunflower lines genotyped with 30 simple sequence repeat (SSR) and 28 inter-retrotransposon amplified polymorphism (IRAP) and retrotransposon microsatellite amplified polymorphism (REMAP) markers identified the genomic regions associated with agro-morphological traits in sunflower. Sukumaran et al. [24] identified single nucleotid polymorphism (SNP) in chromosome 5A and 6A of spring wheat that were significantly associated with yield and yield components. Chickpea is well suited to association analysis for its self-pollination mating system [25]. Association analysis has been performed in chickpea [26–31] and QTLs controlling different traits, such as flowering time [32], branch number [33], pod and branch number/plant and plant hairiness [34], 100-seed weight and seed coat colour [35,36], seed protein content [37], resistance to Fusarium wilt and Ascochyta blight [38], as well as drought tolerance [26], have been identified.

The power of association analysis to identify and characterize loci/genes associated with different complex traits is highly affected by admixtures of populations [39]. It is generally considered that good understanding of the population structure in an association panel is necessary to avoid identifying false positive correlations between markers and traits [40]. The objectives of the present study were to characterize the population structure within chickpea genotypes and to identify AFLP markers associated with different agro-morphological characteristics, cold and drought tolerance as well as Ascochyta blight resistance in chickpea using association analysis.

Materials and methods

Plant materials

A population of 44 chickpea genotypes comprising cultigens, landraces, internationally developed improved lines and wild relatives was used for association analysis (Table 1). The genotypes of the association panel were prepared from six geographically diverse regions including International Center for Agricultural Research in the Dry Areas (ICARDA) (28 genotypes), Turkey (five genotypes), Iran (eight genotypes), USSR (one genotype), Egypt (one genotype) and Spain (one genotype) with Latitude of 5°30 W, 40°02 N to 45°09 E, 37°12 N (Figure 1).

Table 1. Code and characteristics of 44 chickpea lines used in the present study.

code	Genotype	Pedigree	Type	Subgroup	DT	TAB	CT
g1	ILWC 81	C. reticulatum IG	wild	mixed	I	HT	HT
g2	ILWC 106	C. echinuspermum IG	wild	mixed	HS	HT	HT
g3	ILWC 139	C. judaicum IG	wild	blue	RS	HT	HT
g4	ILWC 181	C. reticulatum IG	wild	blue	RS	HT	HT
g5	ILWC 235	C. echinuspermum IG	wild	mixed	RT	HT	HT
g6	C. oxyodon L.-4	C. oxyodon L	wild	red	RT	HT	HT
g7	C. oxyodon L.-9	C. oxyodon L	wild	red	RT	HT	HT
g8	Shahindezh	Landrace	Desi	yellow	RS	T	HS
g9	Germezi Nukhud	Landrace	Desi	yellow	RT	T	HS
g10	Mahabad	Landrace	Desi	mixed	I	T	RT
g11	Kaka	Not traced	Desi	green	I	HT	RT
g12	Piruz	Not traced	Desi	green	I	HT	T
g13	Qazvin	local variety	Kabuli	mixed	I	HT	RS
g14	ILC 3279	Landrace	Kabuli	green	I	HT	HS
g15	ILC533	Not traced	Kabuli	mixed	I	S	SS
g16	ILC 8262	Not traced	Kabuli	yellow	I	HT	T
g17	ILC 8617	ILC482 (mutation)	Kabuli	yellow	I	HT	T
g18	X 96TH3 K3	FLIP 91-150C X FLIP 93-93C	Kabuli	yellow	I	HT	I
g19	X 95TH5 K10	FLIP 91-149C X FLIP 93-194C	Kabuli	yellow	RT	HT	RS
g20	Sel 93TH2 4416	ILC482 X NEWC 36	Desi	mixed	RT	RS	I
g21	Sel 93TH2 4477	ILC3470 x ILC 8617	Desi	yellow	I	RS	RT
g22	Sel 95TH1 716	ILC482 X ILWC182	Desi	mixed	RT	I	T
g23	Sel 95TH1 744	ILC482 X NEWC 36	Desi	green	RT	HS	I
g24	Sel 95TH1 745	ILC482 X NEWC 36	Desi	green	RT	HS	T
g25	Sel 96TH1 1403	ILC482 X NEWC 36	Desi	yellow	RT	HS	HT
g26	Sel 96TH1 1404	ILC482 x NEWC 36	Desi	green	RT	HS	T
g27	Sel 96TH1 1406	ILC482 x NEWC 36	Desi	mixed	RT	HS	T
g28	Sel 96TH1 1439	ILC482 x ILWC 182	Desi	yellow	RT	RT	HT
g29	Sel 96TH1 1484	ILC482 x NEWC 36	Desi	green	I	HS	T
g30	Sel 96TH1 1485	ILC482 x NEWC 36	Desi	green	RT	S	HT
g31	Sel 96TH1 1488	ILC482 x NEWC 36	Desi	mixed	RT	S	I
g32	Sel 98TH1 744	ILC482 x NEWC 36	Desi	red	I	S	T
g33	Sel 98TH1 1518	Not traced	Desi	mixed	I	S	T
g34	Azad	FLIP 93-93C	Kabuli	yellow	I	HT	I
g35	Flip 93-261C	X90TH566/FLIP87-69CXILWC219	Kabuli	red	T	RS	I
g36	Flip 98-258C	Not traced	Kabuli	mixed	I	HT	RS
g37	Flip 99-34C	X96TH8/FLIP91-150CXFLIP91-105C	Kabuli	green	T	HT	HS
g38	Flip 99-66C	X96TH61/(F91-159XILC1278)XF91-149C	Kabuli	mixed	T	HT	RT
g39	Flip 01-06C	X98TH23-1-BH-20/ILC1799XFLIP92-148	Kabuli	red	I	HT	HS
g40	Flip 03-46C	X99TH 13/FLIP93-62CX FLIP93-50C	Kabuli	green	I	HT	HS
g41	Flip 03-112C	X00TH51/FLIP 98-52C x FLIP 98-47C	Kabuli	green	I	HT	I
g42	x2001TH 45	X2000TH 66XS98008	Kabuli	mixed	I	RS	I
g43	x2003TH 21	(S 00 791 x FLIP 98-23C) x ICCV2	Kabuli	mixed	T	HT	HS
g44	x2003TH 164	FLIP 98-50C x FLIP 97-90C	Kabuli	green	T	HT	S

Note: R: representative of resistant; HT: high tolerant; T: tolerant; RT: relatively tolerant; I: intermediate; RS: relatively susceptible; S: susceptible; HS: high susceptible; SS: strongly susceptible.

Figure 1. Geographical distribution and origin map of 44 chickpea genotypes.

Phenotypic data evaluations

The selected set was planted in a randomized complete block design with three replications at the Urmia rainfed research station. The traits, including number of stems (NS), plant height (PH), days from the first effective raining after sowing to 50% flowering (DF), days from the first raining after sowing to 90% maturity (DM), 100-seed weight (100SW) and seed yield (YLD) as well as susceptibility to Ascochyta blight (TAB), and cold (CT) and drought tolerance (DT), were recorded (see the Online Supplemental Data). For the visual screening of cold tolerance, we used a scale of 1–9, where 1 denotes resistant, with no visible symptoms of damage; 2, highly tolerant, with slight foliar damage (1%–10% of leaflets show withering) and up to 10% of branches showing withering and drying; 3, tolerant, with slight foliar damage (11%–20% of leaflets show withering) and up to 20% of branches show withering and drying, no plant killing; 4, relatively tolerant, with foliar damage (21%–30% of leaflets show withering) and up to 30% of branches showing withering and drying, no plant killing; 5, intermediate, 41%–60% of leaflets and 21%–40% of branches show withering and drying, up to 5% plant killing; 6, relatively susceptible; 61%–79% of leaflets and 41–60% of branches show withering and drying, 6%–25% plant killing; 7, susceptible, with 81%–99% of leaflets and 61%–90% of branches showing withering and drying, 26%–50% plant killing; 8, highly susceptible, with 81%–99% of leaflets and 91%–99% of branches showing withering and drying, 51%–75% plant killing; and 9, strongly susceptible, with 100% plant killing [41]. For drought tolerance evaluation, also on the basis of ICARDA recommendation, we used the drought tolerance score, where 1 indicated free, early flowering, good early plant vigour, 100% pod setting; 2, highly tolerant, early flowering, good early plant vigour, 96%–99% pod setting; 3, tolerant, early flowering, good early plant vigour, 86%–95% pod setting; 4, moderately tolerant, early flowering, moderate early plant vigour, 76%–85% pod setting; 5, intermediate, medium flowering, poor early plant vigour, 51%–75% pod setting; 6, moderately susceptible, medium flowering, lack of early plant vigour, 26%–50% pod setting; 7, susceptible, late flowering, lack of early plant vigour, 11%–25% pod setting; 8, highly susceptible, late flowering, lack of early plant vigour, 1%–10% pod setting; 9, 100% plants killed, lack of early plant vigour, no flowering, no pod setting [42]. The blight severity was recorded on a 1–9 scale as suggested by Singh and Reddy [43], taking into consideration the extent of damage to foliage, including leaves, petioles and stems.

Genomic DNA isolation, PCR conditions and electrophoresis

Young leaf tissues of two-week-old seedlings were cut and were subjected to vacuum freeze-drying for dehydration. Genomic DNA was extracted from leaf tissues using the modified method described by Dellaporta et al. [44]. An association panel was fingerprinted with six AFLP primer combinations (Table 2). To identify the most informative AFLP primer combination, the means of the polymorphism information content value, the effective multiplex ratio and the resolving power of markers (Rp) were considered [45]. The DNA concentration in each sample was adjusted to 50 ng μL⁻¹. To 6.0 μL of diluted DNA in suitable 1.5-mL tubes, 0.25 μL EcoRI, 0.25 μL MseI, 4.0 μL 10X buffer and 9.5 μL ddH₂O were added. The tubes were incubated at 37 °C overnight to digest. The DNA digestion step was confirmed by running the products in a 1% agarose gel. In the next step, to 20 μL of digested DNA, 0.5 μL EcoRI adapter, 0.50 μL MseI adapter, 0.50 μL adenosine-5'-triphosphate (ATP) (10 mmol L⁻¹), 1.0 μL T4 DNA ligase, 2.0 μL T4 buffer and 0.5 μL distilled water were added. All tubes were incubated at 37 °C for 3 h. Selective amplification was performed according to the manufacturer's instructions with unlabelled EcoRI and MseI (Generay Biotech, Shanghai, China) primers. For 3.0 μL of pre-amplified products, 3.0 μL EcoRI, 3.0 μL MseI, 1.5 μL 10X buffer, 2.5 μL deoxyribonucleoside triphosphate mixture (dNTPs) (2 mmol L⁻¹), 0.2 μL Tag DNA polymerase (5U μL⁻¹) and 4.5μl H₂O were added. Touchdown polymerase chain reaction (PCR) was used for amplification (Eppendorf 5331 Gradient MasterCycler, Hamburg, Germany) as follows: 94 °C for 4 min, 1 cycle of 94 °C for 30 s, 65 °C for 30 s, 72 °C for 1 min and was followed by 10 cycles with a decrease of annealing temperature at 1 °C per cycle. This was followed by 23 cycles of 94 °C for 30 s, 56 °C for 30 s and 72 °C for 1 min.

Table 2. Description of the primers used for amplified fragment length polymorphism analysis.

Primer	Function	Sequence
EcoRI-1	Adapters	5′-CTCGTAGACTGCGTACC
		CATCTGACGCATGGTTAA-5′
MseI -1	Adapters	5′-TACTCAGGACTCATC
		GAGTCCTGAGTAGCAG-5′
EcoRI(E01)	Pre-selective amplification	5'-GACTGCGTACCAATTCA-3'
MseI(M02)	Pre-selective amplification	5'-GATGAGTCCTGAGTAAC-3'
EcoRI +3	Selective amplification	MseI +3
EcoRI (Primers +3)	Selective amplification	5'-GACTGCGTACCAATTCACC-3'
EcoRI (Primers +3)	Selective amplification	5'-GACTGCGTACCAATTCAAC-3'
MseI (Primers +3)	Selective amplification	5'-GATGAGTCCTGAGTAACAG-3'
MseI (Primers +3)	Selective amplification	5'-GATGAGTCCTGAGTAACTG-3'
MseI (Primers +3)	Selective amplification	5'-GATGAGTCCTGAGTAACAC-3'

E seq: GACTGCGTACCAATTC; MseI: GATGAGTCCTGAGTAA.

Four microliters of reaction products were mixed with 1 µL of formamide dye (98% formamide, 10 mmol L⁻¹ ethylenediaminetetraacetic acid (EDTA), 0.25% v/v bromophenol blue and 0.25% v/v xylene cyanol FF) and denatured at 95 °C for 4 min and quenched on ice. The products were resolved in a 1.2% (w/v) metaphor agarose gel (Sigma–Aldrich, St. Louis, MO) in 0.5X Tris-Borate-EDTA (TBE) buffer, stained with ethidium bromide (1.0 μg mL⁻¹) and photographed under ultraviolet (UV) light using a Gel-Doc image analysis system (Gel Logic 212 PRO, Woodbridg, NJ). The standard 50–900 bp AFLP DNA ladder (Fermentas, Carlsbad, CA) was used in the flanking lanes. Compared to polyacrylamide-gel electrophoresis, or automated analysis, agarose-gel electrophoresis is the most-appropriate and safe technology for routine analyses of these types of markers.

Population structure and association analysis

The AFLP data were scored as dominant markers for each locus. Population structure was analysed using a model-based Bayesian approach in the software package Structure 2.3.4 [46]. Five independent runs were performed, setting the number of sub-populations (K) from 1 to 10, the burn in time and MCMC (Markov Chain Monte Carlo) replication number both to 100,000, and a model for admixture and correlated allele frequencies. The K value was determined by the log likelihood for each K; Ln P (D) = L (K), and an ad hoc statistic ΔK based on the second-order rate of change in the likelihood (ΔK) [47]. The appropriate value of K was selected when the estimate of Ln P (D) reached a minimum stable value. Inferred ancestry estimates of individuals (Q-matrix) were derived for the selected sub-population [46]. Association analysis was performed to analyse marker–trait association by the structured association approach using ancestry coefficient (Q-values) estimates as covariate in a general linear model (GLM) function using TASSEL 2.1. Multiple testing corrections were performed by adjusting maker probability values for multiple test runs by a permutation test in the TASSEL 2.1 software.

Results and discussion

Morphological diversity

Analysis of variance revealed significant differences among the studied genotypes for most of the studied traits (Table 3). The experimental coefficient of variation (CV) ranged from 2.22 to 89.13. However, for the majority of traits, the values were less than 20%. The number of stems per plant ranged from 1 (C. oxyodon) to 7 (ILWC139); the days from the first effective raining after sowing to 50% flowering, from 168 (Sel95TH1) to 188 days (Flip99-34C and Flip03-46C); the plant height, from 16 (ILWC235) to 78 cm (C. oxyodon); the 100-seed weight, from 12 (Kaka) to 35.6 g (Flip03-46); the seed yield; from 31.8 (ILWC181) to 4190.4 kg ha⁻¹ (ILC8262); the resistance to Ascochyta blight disease, from 2 (wild types, FLIPs and Sel95THs genotypes) to 8 (Sel95THs genotypes); the cold tolerance, from 2 (wild types, FLIPs and Sel95THs genotypes) to 9 (ILC533). Of the studied genotypes, 42.8% were tolerant to drought stress. Most of the studied genotypes showed resistance to Ascochyta blight disease. The majority of chickpea landraces (83.3% and 50%) were susceptible to drought and cold stress, respectively. On the other hand, 78.6%, 50%, 12.1% and 50% of the genotypes from Sel95THs, ILCs and FLIPs groups, respectively, were tolerant to cold stress, whereas 71.4% and 37.5% of the genotypes from Sel95THs and ILCs groups, respectively, were tolerant to drought stress. Most of the FLIPs and ILC genotypes were resistant to Ascochyta blight disease but the majority of the Sel95THs genotypes (85.7%) were susceptible. The studied wild genotypes were tolerant to cold and Ascochyta blight disease. Comparing the studied chickpea types, 80.8% of desi type were cold tolerant, whereas just 16.7% of kabuli types were tolerant. Regarding the tolerance to Ascochyta blight disease, 50% and 83.3% of desi and kabuli types were tolerant, respectively. In terms of drought stress tolerance, 53.8% and 33.3% of the two types were tolerant, respectively (Online Supplemental Data).

Table 3. Mean of square and some statistical parameters for 44 genotypes of chickpea comprising cultigen, landraces, internationally developed improved lines and wild relatives.

Source of variation	NS	DF (days)	DM (days)	CT	TAB	PH (cm)	100SW (gr)	YLD (kg ha^−1)
Block	1.1031**	626.58**	267.96**	2.06**	5.58 ^ns	708.93**	95.82**	12,057,874.70**
Genotype	1.23**	42.38**	31.89^ns	1.67**	6.39**	364.27**	77.32**	661,609.50**
Mean	3.38	179.80	224.60	1.87	3.88	40.57	26.79	637.97
CV	19.58	2.48	2.22	19.75	39.09	15.60	11.92	89.13
Max	5.00	189.30	229.00	6.00	6.70	57.00	35.70	1663.30
Min	2.00	173.60	207.00	1.30	1.30	22.30	15.00	71.00
SE	±0.07	±0.54	±0.49	±0.07	±0.17	±1.11	±0.51	±70.26

Note: NS: number of stems; DF: days from the first effective raining after sowing to 50% flowering; DM: days from the first effective raining after sowing to 90% maturity; TAB: tolerance to Ascochyta blight; CT: cold tolerance; DT: drought tolerance; PH: plant height; 100SW: 100-seed weight; YLD: seed yield; CV: coefficient of variation; SE: standard error.

The potential seed yield of chickpea, about 6 t ha⁻¹, is far below the world average yield (0.98 t ha⁻¹) [48]. This is due to several complex field environmental factors such as the global climate change [49], rainfed sowing, and poor farm management, lack of well adapted cultivars with low levels of tolerance to different biotic and abiotic stresses. High genetic variability was observed among the chickpea genotypes for the studied agro-morphological traits as well as resistance to biotic and abiotic stresses. Our findings showed that considerable part of genotypes involved in the current minicore germplasm were drought tolerant. Toker et al. [50] reported that, under drought stress conditions, all accessions of perennial wild Cicer species were significantly superior to the annual wild and cultivated chickpeas, including the best tolerant chickpea, ICC4958. The results showed that the chickpea sowing time in dry lands is a major parameter in drought tolerance by accelerating the plant life cycle before dramatic water deprivation. Drought escape through short and/or better timed life cycles, matched with large total biomass, is the primary form of crop adaptation to low rainfall Mediterranean-type conditions [51]. Therefore, achieving successful production under terminal drought stress depends on the development of short-season (earliness) varieties that enable the crop to escape severe soil-water deficits [52]. Drought avoidance as an alternative mechanism corresponds to the maintenance of high tissue water potential and consists in mechanisms that both reduce water loss from plants, via stomatal control of transpiration, and also maintain water uptake through an extensive and prolific root system [53]. In the present study, we took advantage of both mechanisms, so the suitable drought tolerance varieties showed acceptable yield at stressed conditions.

The results showed that Cicer genotypes are well adapted to germination and seedling establishment in moderate conditions, followed by vegetative growth in cooler conditions which are typical following autumn or winter sowing in temperate Mediterranean and cool West Asia regions. The cold tolerance scores of Cicer species indicated that the level of cold tolerance in all accessions of C. reticulatum and C. echinospermum were superior to those of well-known cold tolerant varieties, ILC 8262 and ILC 8617. These results were in agreement with the data reported by Robertson et al. [54] and Singh et al. [55]. They found that most of the accessions of C. bijugum, C. echinospermum and C. reticulatum are tolerant to cold and have significantly higher levels of cold tolerance than cultivars. The results indicated that some desi types, such as Sel 95TH 1716 and Sel 95TH 1744, as well as kabuli types, such as FLIP 93-261C and x03TH21, showed higher cold tolerance during the early seedling stage and withstood temperatures of −15.6 °C without snow covering, under dry cold conditions. The results also indicated that, for colder areas of highlands, genotypes adapted for winter or autumn sowing could be obtained by the use of genes for cold tolerance from the wild Cicer species. In a preliminary screening for cold tolerance, Toker [7] found that 11 accessions of C. reticulatum and 4 accessions of C. echinospermum were highly cold tolerant and superior than that the best cold tolerant cultigens, ILC 8262 and FLIP 93-53C. This was supported by our findings, and the findings of Singh et al. [56].

The results showed that all wild Cicer genotypes and landraces as well as some kabuli types were resistant to Ascochyta blight. Most of desi types showed susceptibility to disease. Reddy and Sing [15] evaluated 9574 desi and 3836 kabuli germplasm accessions from the world collection of chickpea for resistance to Ascochyta blight both in vegetative and podding stages and reported that 11 kabuli and 6 desi accessions were resistant. A major part of reports focus on strategies to minimize the damage caused by Ascochyta blight [16]. Factors such as using Ascochyta blight-free seeds, one to two years of non-host crops for warm and wet areas, 3–4 years crop rotation for cold and dry areas, optimum sowing date, deep sowing, optimizing plant density, balanced nutrition and alternative sowing patterns should be considered as useful means of reducing this disease [16].

Molecular diversity and population structure

Using 6 AFLP primer combinations, 64 clear bands were generated. The number of polymorphic markers varied from 6 (EcoRI-ACC/ MseI-CAG) to 46 (EcoRI-AAC/ MseI-CAG) for primer combinations. Detailed information on the levels of genetic diversity of genotypes is available elsewhere [45].

Analysis of population structure was carried out by using a model-based clustering approach that assigns individuals into one or more sub-populations probabilistically based on the allele frequencies detected at different loci [46]. The structure clustering results from K = 2–10 are shown in Figure 2. As indicated in Figure 2, the K values (number of sub-populations) steadily kept on increasing until K = 4, indicating that the 44 genotypes in this study are clustered into four sub-populations. Genetic variation among the four identified sub-populations was tested using F-statistics, estimated from pairwise comparisons as a measure of genetic distance between sub-populations. The F-statistics values between sub-populations were significant (P = 0.01) and ranged from 0.02 to 0.96, supporting the existence of genetic structure. Of all genotypes, 65.9% were assigned into the corresponding subgroups, and the remaining ones were categorized into the ‘mixed’ genotypes based on their Q-values (Table 1 and Figure 3). In our association panel based on AFLP markers, the r² values ranged from 0.082 to 0.333 with an average of 0.027 and D' ranged from 0.001 to 0.999 with an average of 0.356 (Figure 4). Of 576 AFLP marker pairs, 7.76% showed a significant level of LD (P < 0.05).

Figure 2. Bilateral charts to determine the optimal number of K identified by Structure program.

Figure 3. Genetic relatedness of 44 genotypes of chickpea comprising cultigen, landraces, internationally developed improved lines and wild relatives with 6 AFLP primer combinations as analysed by the Structure program.

Figure 4. LD plot generated by AFLP marker pairs. The upper diagonal shows the linkage disequilibrium rates (D´) among each pair of markers and the below diagonal corresponds to their P-values.

Association analysis

By applying the GLM procedure and multiple testing corrections in the TASSEL software, one AFLP marker was identified to be associated with genes controlling the stem number; two markers, with the days from the first effective raining after sowing to 50% flowering and two ones, to 90% maturity; two markers, with partial tolerance to Ascochyta blight; one marker, with the 100-seed weight and two ones, with the yield (Table 4). None of identified markers for plant height, cold and drought tolerance were statistically significant. Among the identified AFLP markers, six markers were common, whereas 12 ones were specific.

Table 4. QTLs detected for the studied agro-morphological traits in chickpea using the GLM procedure.

Trait	Locus	F_Marker	P_Marker	P-adj_Marker	R^2
NS	M4	6.4293	0.0156	0.2867^ns	0.1319
NS	M16	4.1074	0.0499	0.8531^ns	0.0890
NS	M19	10.1468	0.0029	0.0270*	0.1918
NS	M40	4.2585	0.0461	0.7742^ns	0.0920
NS	M59	5.1173	0.0297	0.7493^ns	0.1083
DF	M4	5.8454	0.0207	0.3387^ns	0.1335
DF	M6	8.5599	0.0058	0.0549*	0.1838
DF	M7	10.2297	0.0028	0.0170**	0.2119
DF	M9	8.3008	0.0066	0.0689^ns	0.1793
DF	M16	5.4390	0.0252	0.5035^ns	0.1254
DF	M40	7.1448	0.0111	0.1718^ns	0.1584
DF	M51	6.7238	0.0135	0.2667^ns	0.1505
DM	M7	6.0390	0.0188	0.3477^ns	0.1351
DM	M9	5.6133	0.0232	0.4466^ns	0.1268
DM	M16	9.4508	0.0040	0.0280*	0.1959
DM	M17	10.5953	0.0024	0.0150**	0.2143
DM	M20	4.6801	0.0370	0.7682^ns	0.1081
DM	M27	3.8713	0.0566	0.9311^ns	0.0912
DM	M51	6.0959	0.0183	0.2408^ns	0.1362
TAB	M6	4.4174	0.0421	0.8192^ns	0.0913
TAB	M7	23.0527	0.0000	9.99E−04**	0.3334
TAB	M9	11.7478	0.0015	0.0070**	0.2077
TAB	M10	4.4693	0.041	0.7433^ns	0.0923
TAB	M50	6.8767	0.0124	0.2218^ns	0.1345
CT	M6	7.03	0.0115	0.1818^ns	0.1445
CT	M40	7.5266	0.0091	0.1568^ns	0.1530
CT	M62	5.649	0.0225	0.4715^ns	0.1197
PH	M5	4.2697	0.0459	0.9211^ns	0.0905
PH	M6	4.4885	0.0409	0.7053^ns	0.0946
PH	M37	4.4862	0.0410	0.7722^ns	0.0946
PH	M40	3.8020	0.0588	0.9441^ns	0.0815
PH	M43	5.3061	0.0270	0.5455^ns	0.1097
PH	M44	4.2210	0.0470	0.8811^ns	0.0896
PH	M51	7.6077	0.009	0.1179^ns	0.1492
100SW	M3	4.9304	0.0326	0.7333^ns	0.1006
100SW	M8	4.7211	0.0363	0.7542^ns	0.0968
100SW	M20	9.8208	0.0034	0.0420 *	0.1794
100SW	M44	4.7499	0.0357	0.7752^ns	0.0973
100SW	M63	6.9094	0.0124	0.2677^ns	0.1346
YLD	M5	13.95	6.31E-04	9.99E-04**	0.2256
YLD	M8	4.6475	0.0377	0.7333^ns	0.0919
YLD	M12	15.879	3.05E-04	9.99E-04**	0.2474
YLD	M37	4.6179	0.0383	0.7113^ns	0.0914
YLD	M54	5.006	0.0314	0.6603^ns	0.0982

Note: R² is the per cent of individual phenotypic variance explained by each QTL.

The significant level of LD found in the present study supports the idea of using natural populations for genomic studies. Being a predominantly self-pollinating species, chickpea is expected to show higher LD extent, 150–200 Kb [36] and 500–600 Kb [57], than outcrossing species like maize with (≈2 Kb) LD extent [58]. Several factors which cause LD content increasing include low recombination rate, small population size, population admixture, inbreeding, genetic isolation between lineages, population subdivision, genetic drift and epistasis [22,59]. The admixture of populations is the most important factor influencing marker–traits associations [39]. Therefore, good understanding of the population structure is a prerequisite in association mapping and necessary to avoid identifying false-positive correlations between markers and traits [40,60]. Using polymorphic bands produced by six AFLP primer combinations, the association panel subdivided into seven sub-populations. The GLM was used to identify genes controlling some agro-morphological traits in the studied chickpea germplasm mini-core collection. Heretofore a considerable number of QTLs have been identified for desirable agronomic traits in chickpea on different linkage groups (LG 2, 3, 4, 6 and 8) [61]. For instance, Iruela et al. [62] identified three QTLs for resistance to Ascochyta blight on linkage group 4 (LG4) (QTL AR1 and QTL AR2) and 2 (QTL AR3). Cho et al. [63] reported two major QTLs on LG 2 close to the GA16 and TA37 loci controlling the resistance to Ascochyta rabiei pathotype I, as well as two QTLs for resistance to pathotype II on LG 4 such that one is linked to CaETR or GAA47 and the other is linked to TA72/ScY17. Bian et al. [64], by comparing three chickpea LGs comprising QTLs for tolerance to A. rabiei, found that LG3 correspond to the sub-centromere region of chromosome C in chickpea. LG8 is located on the long arm of chromosome H (W-Ca-LG8) and LG4 is located in the sub-centromere region of chromosome B (W-Ca-LG4). Based on bibliography, the majority of AB-resistance QTLs were reported mainly on two LGs, CaLG02 and CaLG04 [61,65]. Two important QTLs (Q3-1 and Q1-1) were reported for drought tolerance from ILC 588 × ILC 3279 population which was located on LG3 and LG1 [66,67]. In order to understand the complex nature of drought tolerance in chickpea, Thudi et al. [68] carried out an association mapping study and identified 15 markers significantly associated with five root traits (root dry weight, root length density, root surface area, root volume and rooting depth). Among them, seven markers showed significant relation with a single trait (rooting depth) and two markers (NCPGR7 and DR_237) showed associations with more than one trait. They concluded that these two markers were associated with co-localized/pleiotropic QTLs. In the case of cold tolerance, we were not able to identify any reported QTLs. The association studies in chickpea have revealed that the linkage-based QTL analyses can be complemented by LD-based association studies [68]. The markers identified in this study could be transferred to sequence characterized amplified region (SCAR) markers to increase the efficiency of MAS in chickpea breeding programmes for biotic and abiotic stresses.

Conclusions

High genetic diversity was observed for the studied agro-morphological traits in the association panel. Using a GLM, several AFLP markers were detected for the studied traits. The identified markers can be transferred to sequence characterized amplified region markers for facilitating MAS. Detection of molecular markers associated with genes controlling different traits could increase the efficiency of MAS in chickpea breeding programmes for biotic and abiotic stresses.

Acknowledgments

The authors are grateful to Dr Razban A. and the Institute of Biotechnology, Urmia University, Iran for the scientific support of this work.

Disclosure statement

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.