Investigation of the domestication of common bean (Phaseolus vulgaris) using multilocus sequence data

Additional keywords: ABC approach, association mapping, bottleneck, demography, genepools, linkage disequilibrium, races.

Domestication is a complex process in which human usage of plant and animal species has led to morphological and physiological changes that made them genetically different from the wild types and better adapted to different agro-ecosystems (Glémin and Bataillon 2009). Beginning in the Epipalaeolithic and extending into the Neolithic period (13 000–11 000 years ago), cultivation started with just a few plant species as food sources (Fuller 2007). Many morphological and physiological changes were associated with the process of domestication and are termed the ‘domestication syndrome’ (Glémin and Bataillon 2009). The study of domestication as an evolutionary model can identify events associated with the origins of crop species and describe the selective pressures experienced by domesticated taxa.

Common bean (Phaseolus vulgaris L.) is the most important dietary legume in the world because of its high concentrations of protein, fibre and complex carbohydrates. It is an especially important food for many developing countries in Latin America, Asia and Africa (Graham and Vance 2003). It is estimated that the global harvest is ~18.7 million tons and is grown on 27.7 million ha in ~148 countries (Gepts et al. 2008).

Based on the discovery of wild common bean in Argentina (Burkart and Brücher 1953) and Guatemala (McBryde 1947) and archaeological remains found in the Americas (Kaplan and Kaplan 1988; Kaplan and Lynch 1999), common bean is commonly thought to have originated in the Americas. Two large gene pools of wild types were identified based on phaseolin seed protein variation (Gepts et al. 1986; Gepts 1990), DNA marker diversity (Becerra Velasquez and Gepts 1994; Sonnante et al. 1994; Freyre et al. 1996; Tohme et al. 1996), morphology (Evans 1976; Gepts and Debouck 1991), isozymes (Koenig and Gepts 1989), mitochondrial DNA variation (Khairallah et al. 1992) and amplified fragment length polymorphism (AFLP) (Rossi et al. 2009) and short sequence repeats (SSR) (Kwak and Gepts 2009) marker data. The Middle American (MA) gene pool extends from Mexico through Central America and into Venezuela, whereas the Andean gene pool is found in Peru, Chile, Bolivia and Argentina. Recently, multilocus sequence data considered demographic events in wild common bean (S. Mamidi, M. Rossi, D. Annam, S. M. Moghaddam, R. K. Lee, R. Papa, P. E. McClean, unpubl. data) and determined that the two gene pools diverged ~110 000 years before present (BP) followed by a geographical adaptation bottleneck in each wild gene pool.

The discovery of landraces in archeological sites dating from 10 000 years BP in Argentina and 7000 years BP from Mexico (Kaplan et al. 1973) suggests that common bean was domesticated early in Middle and South America. Accelerator mass spectrometry (AMS) analyses provide evidence of cultivation of common bean before ~2500 BP in Tehuacan valley, 1300 BP in Tamaulipas, 2100 BP in Oaxaca and 4400 BP in the Peruvian Andes (Kaplan and Lynch 1999). Linguistic evidence suggests the presence of bean 3400 BP in Middle America (Brown 2006). Mexico is suggested as one centre of domestication of common bean (Gepts et al. 1986; Gepts 1988; Smith 1995; Piperno and Flannery 2001; Doebley 2004) and another centre of domestication is suggested in the Andes (Gepts et al. 1986; Gepts 1998). Apparently, the divergence of ancestral wild common bean at 110 000 years BP provided the genetic basis for the domestication within the two gene pools.

The presence of distinct groups of landraces has been described at both the morphological (Singh et al. 1991a) and molecular level (Gepts et al. 1986; Singh et al. 1991b; Becerra Velasquez and Gepts 1994; Freyre et al. 1996). Singh et al. (1991a) classified the MA cultivars into three races – lowland race Mesoamerican (M) and highland races Durango and Jalisco. Race Durango occupies the semi-arid northern highlands of Mexico. The pinto, great northern, medium red and pink market classes are assigned to this race. Race Jalisco overlaps the southern distribution of the race Durango. Race Mesoamerica, the third MA race, occupies the lowlands of Latin America from Mexico to northern Colombia and Venezuela. Black, navy and small red market classes represent this race. The Andean gene pool is divided into three races based on morphological and ecological criteria – Nueva Granada, Peru and Chile (Singh et al. 1991a). Race Nueva Granada is the most widely cultivated Andean race and includes the majority of the commercial large seeded cultivars. It is grown at mid-altitudes of the Andes and Africa, in warm lowland environments of Brazil, Mexico and the Caribbean and in the temperate climates of North America and Europe. The dark and light red kidney, cranberry and most horticultural snap bean market classes are found within this race. Race Peru is found in the Andean highlands, whereas race Chile occupies the southern Andes.

Although it has been established that domestication was an independent event in each gene pool, the number of domestication events in each gene pool is debated. A single MA domestication event (Gepts et al. 1986; Papa and Gepts 2003; Kwak et al. 2009, Rossi et al. 2009) would imply that divergence into races followed the domestication process instead of resulting from separate domestications. Beebe et al. (2000) used random amplified polymorphic DNA (RAPD) data to suggest that the MA gene pool was the result of two distinct domestication events. This was further supported by chloroplast data (Chacón et al. 2005) and sequence diversity data for DFR and CHI introns (McClean et al. 2004; McClean and Lee 2007).

Archeological and multilocus sequence data provide complementary information to understand the demographic and evolutionary events associated with the domestication of a species (Doebley et al. 2006; Burke et al. 2007). Recent population genomics studies suggest that domestication affects the entire genome and that selection acts on a large number of loci (Wright et al. 2005; Caicedo et al. 2007), so a multilocus approach is appropriate to study the effects of domestication and selection. Similar results were obtained by Papa et al. (2007) on P. vulgaris using 2509 AFLPs. Computationally, approximate Bayesian computation (ABC) has emerged as a preferred approach to simulate many models consisting of various combinations of demographic parameters (some derived from archeological data). Subsequent statistical analyses can then select the model(s) that best fits the observed summary population genetic data obtained from sequencing multiple loci. The ABC simulation method considers the population summary data and make inferences with less computational time than when all available data are analysed in detail. This method was used successfully to untangle many evolutionary processes in humans (Fagundes et al. 2007; Patin et al. 2009; Scheinfeldt et al. 2009; Ray et al. 2010; Batini et al. 2011) and plants (Ross-Ibarra et al. 2007; Ingvarsson 2008; François et al. 2008). The advantages and disadvantages of the methodology are reviewed extensively by Bertorelle et al. (2010), Lopes and Beaumont (2010) and Csilléry et al. (2010). The objectives of this research were to collect multilocus sequence data from domesticated P. vulgaris landraces, evaluate the nucleotide variation within this collection, consider different models of domestication using ABC and estimate the domestication model and population parameters that best describe domestication within each of the two common bean gene pools. The results were compared to previous data at the same loci for wild relatives to assess the genome wide effects of domestication. To our knowledge, this is the first crop species where multiple demographic models were tested to find the best model of domestication.

A collection of 24 landraces were analysed (Table 1). Based on DFR (McClean et al. 2004) and CHI intron-3 (McClean and Lee 2007) haplotypes, unique landraces were selected to represent races within each of the two common gene pools. Phaseolus coccineus L. genotypes PI 325 589 and PI 325 599 were used as out group members. Thirteen nuclear genes were selected for sequencing with at least one locus on each Phaseolus vulgaris (L.) linkage group (Table 2; McConnell et al. 2010). Two additional Pv08 loci, g776 and D1468, were included. D1468 is associated with several domestication traits including number of pods and days to flowering and maturity (Koinange et al. 1996). Locus g776 maps < 2cM from D1468 (McConnell et al. 2010).

Table 1.  Common bean landraces used for the study of domestication
Accession numbers are from the National Plant Germplasm System

Table 2.  Summary of common bean genes analysed for domestication

DNA was extracted from young leaves (Brady et al. 1998) and fragments from the 13 loci were amplified using standard PCR conditions. The amplified fragments were sequenced from both directions using a Beckman CEQ 2000XL DNA analysis system (Beckman Coulter Inc., Brea, CA). The DNA sequence chromatograms were analysed using the Staden package (Staden 1996). Gene annotation and structure were determined by blastx against the Viridiplantae database at NCBI (http://blast.ncbi.nlm.nih.gov/Blast.cgi, accessed 2 January 2009).

Population structure was determined using the STRUCTURE 2.2 software (Pritchard et al. 2000). Combined data for all loci was used for this analysis. The sequence files were converted into STRUCTURE input format using xmfa2struc (Didelot and Falush 2007). We set k (the number of subpopulations) from 1 to 6 and performed 10 runs for each k value. For each run, a burn in of 100 000 iterations was followed by an additional 500 000 iterations. The Δk method proposed by Evanno et al. (2005) was used to choose the best k value. The assignment of an individual to a subpopulation was based on subpopulation probability values estimated in STRUCTURE. Individuals were assigned to a subpopulation based on a coefficient q_i ≥ 0.7. To further differentiate the subpopulations, a neighbour joining (NJ) tree using the combined loci was built in ClustalX (Larkin et al. 2007) and bootstrapped over 1000 replicates.

Unless otherwise stated, population genetic statistics were calculated using DnaSP 4.90 (Rozas and Rozas 1999). Population differentiation was described using F_st (Hudson et al. 1992) and Hudson’s S_nn (Hudson 2000) statistic with significance determined by 10 000 permutations. The number of shared (S_s) and fixed (S_f) silent sites between gene pools and the number of unique polymorphic sites (S_And and S_MA) in each population were determined. Nucleotide diversity for synonymous and non-coding silent sites were estimated using (1) the Watterson’s estimator (θ_w = 4N_eµ), (2) the average number of pairwise differences per site between sequences in a sample (π), (3) the number of segregating sites (S), (4) the number of haplotypes (H) and (5) haplotype diversity (Hd). Haplotype diversity is a measure of uniqueness of a haplotype in a population. To test the departure from the neutral equilibrium model of evolution, Tajima’s D (D_T; Tajima 1989) was estimated. ZnS, the average R² (linkage disequilibrium (LD) coefficient) over all pairwise comparisons was also calculated. The expected decay of intergenic LD with physical distance was estimated as described by Remington et al. (2001) and Pyhajarvi et al. (2007) by fitting the data a nonlinear regression equation using the NLIN procedure in SAS 9.2 (SAS Institute, Cary, NC, USA).

Six domestication models (Fig. 1) were simulated using Hudson’s ms (Hudson 2002) to find the best domestication model in each gene pool. Models 1 and 2 describe a single domestication event in each gene pool. Models 3 and 4 describe the presence of two domestication events, two races together as an event and the other race as another domestication event. Therefore, models 3 and 4 each have three submodels accounting for all three possible combinations in each of the gene pool. Models 5 and 6 indicate the presence of three domestication events, one for each race. All the models consist of a bottleneck during the start of domestication. Models 2, 4 and 6 are characterised by an exponential growth of the landrace population after the bottleneck, whereas models 1, 3 and 5 are characterised by instantaneous expansion of population after the bottleneck. Since population size has little effect on the simulation results (Tenaillon et al. 2004; Wright et al. 2005), ancestral population size (N_w), present population size (N_L) and effective population size (N_e) were assumed to be 100 000 individuals. For the MA gene pool models, domestication began at 10 000 years BP and ended before 2500 years BP. For the Andean models, domestication began at 10 000 years BP and ended at 4000 years BP. The ending dates were based on the results by Kaplan and Lynch (1999). A variable mutation rate (µ, based on a uniform distribution of 1 × 10^–10 to 1 × 10^–6 substitutions per synonymous site per year) and a symmetric migration rate (m, between 1 × 10^–2 to 10 individuals per generation) were included in the simulations. The bottleneck population size is assumed to be 0.0001 to 100% of effective population size, equivalent to 1–100 000 individuals. For each locus within the model, we simulated 1 000 000 priors. The simulation results were piped into the msstats 0.2.9 software (available at http://molpopgen.org/) to obtain the summary statistics for each of the simulation. For each simulation, the Euclidean distance between the simulated and observed summary statistics (segregating sites (S), number of haplotypes (H) and nucleotide polymorphism (π)) was calculated. We accepted the top 10 000 simulations for each model with a Euclidean distance <0.1. The best model was selected by combining the accepted simulations across all models and estimating the posterior probability of the model in the top 5000 simulations as described by Pritchard et al. (1999), Estoup et al. (2004) and Ray et al. (2010).

Fig. 1.  Models of domestication tested in each genepool. Models 1 and 2 describe a single domestication event, Models 3 and 4 describe the presence of two domestication events, two races together as an event and the other race as another domestication event and so has three submodels accounting for all three possible combinations, Models 5 and 6 indicate the presence of three domestication events, one for each race. Models 2, 4, and 6 are characterized by an exponential growth of the landrace population after the bottleneck, whereas model 1, 3, 5 are characterized by instantaneous expansion of population after the bottleneck. N_A is ancestral population size, N_w is the present size of wildtype populations, N_b, N_b1, N_b2, N_b3, represent the bottleneck population sizes, and N_L, N_L1, N_L2, N_L3 represent size of present day population size of landraces. Time t₂ is the start of domestication, and t₁, t_1a, t_1b, t_1c represent the ending of domestication.

For the best model, the summary statistics for the accepted simulations were reduced in dimensionality using principal component analysis (PCA) using PRINCOMP in SAS 9.2. This also helped produce a set of uncorrelated transformed statistics (Ray et al. 2010). Then the parameters were estimated from the accepted 10 000 simulations for the best model using a general linear model (GLM), described by Leuenberger and Wegmann (2010) using the GLM procedure in SAS 9.2. For estimating the goodness of fit of our model and the parameters, we compared the mean of observed statistics with a posterior distribution of summary statistics of the accepted simulations as described by Pascual et al. (2007) and Ingvarsson (2008).

STRUCTURE analysis defined two subpopulations: one composed of Durango, Jalisco and Mesoamerica landraces; and a second represented by Chile, Peru and Nueva Granada landraces. A similar result was observed with a NJ tree where a 100% bootstrap value supported two gene pools (Fig. 2). There was significant genetic differentiation between the two gene pools with an average F_st = 0.38 (Table 3). This was further confirmed with a significant S_nn value for 10 loci. The number of shared sites is 38 and there were six fixed sites between the two gene pools. The unique polymorphic sites in MA gene pool (32) were greater than in Andean (18).

Fig. 2.  Neighbor joining tree for the combined loci of the landraces under study built in ClustalX bootstrapped over 1000 replicates. Only bootstrap supports >70% are represented here. The tree shows the presence of two differentiated populations.

Table 3.  Tests of population differentiation between the two common bean gene pools at silent sites
Statistically significant differences are indicated: ns, not significant; *, 0.01 < P < 0.05; **, 0.001 < P < 0.01; ***, P < 0.001

The MA gene pool had 69 segregating sites and the Andean gene pool had 56 segregating sites (Table 4). The average number of MA haplotypes (4.15) and haplotype diversity (0.562) was slightly higher than that of Andean (3.08 and 0.438 respectively). The level of nucleotide diversity was heterogeneous among loci. On average, the polymorphism within landraces was θ_sil = 0.0088 for MA and θ_sil = 0.0080 for Andean. For nearly all loci, the ratio of π_nonsyn/π_syn was <1 in both gene pools. Based on D_T, most loci showed a significant departure from neutral equilibrium expectations in the Andean gene pool (Table 5). D_T values were negative for most Andean loci.

Table 4.  Diversity and neutrality estimates for the silent sites at each common bean locus studied for domestication
Note: And, Andean; MA, Middle American; n, sample size; S, number of segregating sites; H, number of haplotypes; Hd, haplotype diversity; θ, Watterson estimator; π, average number of pairwise differences per site between sequences in a sample; π_nonsyn/π_syn, ratio of pairwise differences at non synonymous sites to pairwise differences at synonymous sites for the entire sequence. Statistically significant differences are indicated: *, P < 0.05; **, P < 0.01

Table 5.  Posterior probabilities of models tested
For Andean a = Chile + Nueva Granada; Peru, b = Chile + Peru; Nueva Granada, c = Nueva Granada + Peru; Chile. For Middle America (MA) a = Durango + Jalisco; Mesoamerica, b = Durango + Mesoamerica; Jalisco. c = Jalisco + Mesoamerica; Durango

The value of Zns is 0.430 and 0.627 for MA and Andean landraces. In the MA gene pool, LD decayed at 500 bp (r² < 0.1) whereas in the Andean gene pool, the decay was within 100 bp (r² = 0.1 (Fig. 3). For the landraces, 42 and 47% of pairwise comparisons were significant in MA and Andean gene pools, respectively, and only 2% were significant in MA and Andean when the Bonferroni correction was applied.

Fig. 3.  Linkage disequilibrium decay of the landraces under study. The distance between two polymorphic sites in bp is presented on X-axis and linkage disequilibrium coefficient (r²) is presented on the Y-axis. The scattered dots represent the r² values for individual pairwise combination. The smoothened line represents the non-linear regression line.

In both gene pools, model 2 performed best with a posterior probability of at least 0.26 in MA and 0.21 in Andean (Table 5). For model 2, seven and three principal components explained 98% of the variability of the summary statistics. For model 2, MA domestication started at 8160 years BP and ended at 6260 years BP with a bottleneck size (in % effective wild types population) of ~48% for the MA gene pool (Table 6; Fig. 4). For model 2 within the Andean gene pool, domestication began 8500 years BP and ended 7012 years BP with a bottleneck size of ~47% (Table 6; Fig. 4). Based on the effective population size estimates of 292 362 and 137 248 for the MA and Andean wild gene pools (S. Mamidi, M. Rossi, D. Annam, S. M. Moghaddam, R. K. Lee, R. Papa, P. E. McClean, unpubl. data), it is estimated that the domestication bottleneck population sizes were 139 310 and 72 827 respectively. The estimated migration rate between the wild and domesticated population in each gene pool was ~0.5 migrants per generation.

Table 6.  Parameter estimates for the top two models of common bean domestication in each genepool
Values in parenthesis are 95% CI

Fig. 4.  Domestication parameters of the best model. The population sizes of wild types are adapted from S. Mamidi, M. Rossi, D. Annam, S. M. Moghaddam, R. K. Lee, R. Papa, P. E. McClean (unpubl. data). The population sizes of wildtypes (N_w) and bottleneck (N_b) are presented as number of effective individuals. Time is in years before present. Times t₂ and t₁ represent the start and ending times of bottleneck.

Several simulation summary statistics were compared with the observed mean values to assess the validity of model 2. A reasonable fit was observed for theta (θ_w) (95% CI = 0–0.324 for MA and 0–0.353 for Andean) and D_T (95% CI = 1.451–1.442 for MA and –1.128–0.014 for Andean). The observed means in all cases fell within the 95% distribution of the accepted simulations.

Common bean has a high level of phenotypic diversity illustrated by its wide geographical distribution from northern Mexico to northern Argentina and its adaptation to tropical and temperature environments. Wild types as a whole have the largest level of diversity available (McClean and Lee 2007) and domesticated landraces selected for many important agronomic traits have arisen from the two main wild gene pools. These landraces contain much of the diversity that has been captured in production of cultivars that were mostly developed through hybrid breeding. In this study we evaluated multilocus sequence data from a diverse group of landrace genotypes to assess nucleotide diversity, population differentiation and demography of the species. These results have implications for understanding the history of domestication of common bean, which itself contributes to our understanding of the origin and development of modern cultivation and agronomy (Guo et al. 2010). This study also enhances our understanding of the factors that contribute to LD in present day cultivars which, in turn, has implications for association mapping studies.

The landraces of common bean were divided into two gene-pool specific subpopulations. This observation is consistent with earlier research (Becerra Velasquez and Gepts 1994; Gepts et al. 1986; Singh et al. 1991b; Freyre et al. 1996; McClean et al. 2004; McClean and Lee 2007; Kwak and Gepts 2009). For the landraces we analysed, the nucleotide diversity of MA was higher than that the Andean gene pool, similar to earlier studies (Cattan-Toupance et al. 1998; Beebe et al. 2001; McClean et al. 2004; McClean and Lee 2007; Kwak and Gepts 2009; Rossi et al. 2009). High F_st and significant S_nn values and the presence of fixed sites between the two geographically separated gene pools further support their reproductive isolation. The greater haplotype diversity at silent MA sites suggests a larger effective population size for that gene pool.

We found high variance among our D_T estimates, with both positive and negative values for the MA subpopulation, which may be due to the influence of evolutionary processes on nucleotide variation (Wright and Gaut 2005; Moeller et al. 2007). Also, this could be due to the initial period of positive D_T possible after a bottleneck due to accumulation of intermediate frequency variants (Maruyama and Fuerst 1985; Depaulis et al. 2003). Since the MA gene pool is more diverse than the Andean (Cattan-Toupance et al. 1998; Beebe et al. 2001; Galván et al. 2001; Rossi et al. 2009; Kwak and Gepts 2009), a question for subsequent research is whether selection acted differentially on the MA gene pool that directly influenced neutrality estimates. D1468, a locus previously mapped to a QTL for pod number, days to flowering and maturity (Koinange et al. 1996), had zero diversity in both the gene pools. This may suggest that this locus was selected during the domestication process. However, we do not have the appropriate data to test for the effect of selection at this locus.

The structure analysis of 21 wild types belonging to two gene pools (S. Mamidi, M. Rossi, D. Annam, S. M. Moghaddam, R. K. Lee, R. Papa, P. E. McClean, unpubl. data) along with the landraces under study here suggested the presence of two subpopulations with a probability of assignment of each individual to a group of q_i > 0.7. The landraces and wild types of each gene pool were grouped together using STRUCTURE and NJ procedures (data not shown). The close relationship of these is further confirmed by the low F_st values between the wild and landraces within each gene pool (0.15 in MA and 0.04 in Andean). The two major subpopulations identified here are consistent with earlier work of Rossi et al. (2009), Kwak and Gepts (2009) and McClean et al. (2011). This data clearly suggests that landraces within a gene pool arose by a domestication event specific to that gene pool. The hypothesis of independent domestication events in each gene pool was suggested previously (Gepts et al. 1986; Papa and Gepts 2003; Kwak and Gepts 2009; Rossi et al. 2009).

Since the domestication event is independent in each gene pool, model selection and parameter estimates were performed separately in each gene pool. Also, since each gene pool is further defined by a specific race structure (Singh et al. 1991a; Beebe et al. 2000), models with 1–3 domestication events were tested using the ABC approach. Even though computationally intensive, this approach has successfully described evolutionary events in many species (Fagundes et al. 2007; Ingvarsson 2008; Patin et al. 2009; Ray et al. 2010; Batini et al. 2011). All models tested included a population bottleneck imposed on a founding population that subsequently led to the landrace population. This is an additional example of another crop species that experienced a bottleneck during domestication (Eyre-Walker et al. 1998, Wright et al. 2005; Zhu et al. 2007). We also evaluated models that included both exponential and instantaneous expansion of the landrace populations. The best models selected in each gene pool consisted of a single domestication event. This single domestication model showed similar posterior probabilities for the two gene pools. These results support hypotheses proposed in common bean by previous researchers (Gepts et al. 1986; Papa and Gepts 2003; Kwak and Gepts 2009; Rossi et al. 2009).

The parameters estimated here suggest a short duration with a large founding population. This result is contrary to other species where the founding population was small and the bottleneck duration was long (Tenaillon et al. 2004; Hamblin et al. 2006; Zhu et al. 2007). Unlike most crops, wild-type common beans belong to the same species and are members of the same gene pool as related landraces. Based on their highly similar phenotypes found for landraces and their wild relatives, along with similar levels of nucleotide diversity, it appears that the domestication bottleneck may have only involved a limited number of genes. The strong population structure of common bean and phenotypic similarity between wild types and landraces could explain the short bottleneck duration.

When wild type and landraces were evaluated together, more shared sites than fixed sites were observed (29 and 0 in MA, 36 and 0 in Andean). This a strong indicator of gene flow and is consistent with the earlier results by Papa and Gepts (2003) and Papa et al. (2005) for common bean, as well the bidirectional gene flow observed in 12 of the 13 most important food crops (Ellstrand et al. 1999; Stewart et al. 2003). Moreover, for 7 of these 13 crops, introgression of domesticated traits increased the competitiveness of the related weed species (Ellstrand et al. 1999; Ellstrand and Schierenbeck 2000). In these cases, gene flow may simply reflect the close proximity of wild types and cultivated landraces in production fields where wild genotypes grow next to cultivated fields in native domains (Papa et al. 2005).

As observed for other plants, a bottleneck increases population structure and reduces within subpopulation diversity. Within each landrace gene pool, population differentiation was higher than that observed for wild types (F_st = 0.29; S. Mamidi, M. Rossi, D. Annam, S. M. Moghaddam, R. K. Lee, R. Papa, P. E. McClean, unpubl. data). This is further supported by presence of more fixed sites in the landraces compared to the wild types (S_F = 0; S. Mamidi, M. Rossi, D. Annam, S. M. Moghaddam, R. K. Lee, R. Papa, P. E. McClean, unpubl. data). The higher F_st values for landraces than wild types appear to reflect the cumulative effect of gene flow in wild types over thousands of years, compared to landraces, which are of a much more recent origin. In addition, after the initial domestication event, genotypes became adapted to specific environmental conditions (Kwak and Gepts 2009) and further differentiated into races.

Landrace nucleotide variation was slightly lower than the estimates for the wild types (Table 4). Typically, during domestication, important agronomic characters are selected which result in a genome-wide reduction of genetic diversity in the domesticates (Tanksley and McCouch 1997; Eyre-Walker et al. 1998; Buckler et al. 2001; Diamond 2002; Clark et al. 2004; Papa et al. 2007, Pozzi et al. 2004; Otero-Arnaiz et al. 2005; Vasemagi et al. 2005; Wright et al. 2005; Doebley et al. 2006; Kilian et al. 2007; Zhu et al. 2007). Earlier studies in bean using DNA marker, protein and morphological variation, determined that domesticated landraces indeed contain a subset of the variability found in wild beans (Gepts and Bliss 1986; Gepts et al. 1986; Debouck et al. 1993; Kami et al. 1995; Tohme et al. 1996; Beebe et al. 2001; Chacón et al. 2005). Here we observed that the amount of variability/polymorphisms retained is similar to the diversity estimates found in wild types. These estimates are higher than those proposed by Buckler et al. (2001), who suggested that 60–80% of variability is retained in the domesticated crops. The loss of variability in other crops is within that range: grasses (66%, Buckler et al. 2001), soybean (50%, Hyten et al. 2007; Guo et al. 2010) and rice (70%, Li et al. 2009). In the present study, a high level of loss of diversity is not evident from nucleotide diversity or SNP density data. The SNP density in wild types are one per ~41 bp and ~50 bp for MA and Andean wild types respectively. For landraces, these values are one per ~46 bp and one per ~57 bp. SNP density is higher in other species. For example, wild rice has one SNP per ~19 bp and cultivated rice has one per ~40 bp (Zhu et al. 2007). Similarly there is one SNP per ~19 bp in wild sunflower and one SNP per ~39 bp in cultivated sunflower (Liu and Burke 2006) and one per ~22 bp in teosinte and one per ~30 bp in maize (Tenaillon et al. 2004). This implies that in common bean domestication may have been less severe than in other species. However, we need to consider these results along with other diversity studies using molecular markers (Papa and Gepts 2003; Papa et al. 2007; Kwak et al. 2009; Rossi et al. 2009) where a much larger effect of domestication on the genetic diversity of the common bean was found, particularly in MA genotypes. More analysis at the nucleotide level on larger samples is needed in order to better estimate the selection intensity and the drift associated with domestication in common bean.

In general, it would be expected that LD for self-pollinated species like P. vulgaris would extend to the kilobase level as observed for Arabidopsis thaliana (L.) Heynh (Nordborg et al. 2002), rice (Garris et al. 2003) and soybean (Zhu et al. 2003). The lack of LD decay in Andean landraces is reflective of the low diversity levels and low population differentiation within the gene pool. The mean LD coefficient (Zns) estimate is higher for Andean landraces than MA landraces, which is consistent with earlier results by Rossi et al. (2009).

A higher level of LD in landraces than wild populations is likely due to lower diversity and the short time frame to accumulate recombination events among the domesticated genotypes (Morrell et al. 2005; Caldwell et al. 2006; Rostoks et al. 2006; Hyten et al. 2007). The higher level of LD in Andean compared to MA gene pool in both wild types (0.27 in MA and 0.46 in Andean) and landraces suggests that the higher LD in the Andean gene pool originated before domestication (Rossi et al. 2009) and is suggested to be the result of migration, genetic drift and selection (Rossi et al. 2009). LD decay in landraces is more than that of wild types, consistent with previous estimates that suggested increase in LD decay distance after a bottleneck (Flint-Garcia et al. 2003; Gupta et al. 2005; Li et al. 2009). Similar estimates have been observed in other crops (Morrell et al. 2005; Caldwell et al. 2006; Liu and Burke 2006; Zhu et al. 2007).

Wild common bean is divided into two gene pools (Singh et al. 1991a, 1991b; Blair et al. 2006; Díaz and Blair 2006; McClean et al. 2004; McClean and Lee 2007; Kwak and Gepts 2009; Rossi et al. 2009) that appear to have arisen from a common ancestor 110 000 years BP (S. Mamidi, M. Rossi, D. Annam, S. M. Moghaddam, R. K. Lee, R. Papa, P. E. McClean, unpubl. data). Both wild gene pools arose via a bottleneck probably associated with regional adaptation. Only one domestication event is suggested in each landrace gene pool that is consistent with suggestions from previous research (Kwak and Gepts 2009; Rossi et al. 2009). The domestication events in each gene pool were characterised by a bottleneck of ~50% effective population size and a bottleneck length of ~2000 and 1500 years in the MA and Andean gene pools respectively. After the bottleneck, diversification of landraces into races occurred which differ by morphological and physiological characteristics.

Association mapping uses the linkage disequilibrium in a population of choice to discover QTL for various traits of importance. The major advantage of association mapping vs bi-parental mapping is that it samples more recombination events are available than in a single pair-wise cross. This presumably will lead to higher mapping resolution (Myles et al. 2009). As we discuss below, the results here have important implications for association mapping in common bean.

Different demographic factors influence LD in different ways. Bottlenecks reduce genetic variation and change the gene frequency spectrum by removing low-frequency alleles (Hamblin et al. 2011). The extent of LD increases due to the elimination of a subset of recombination events. Selection also increases LD distance that may extend beyond the average for the whole genome (Myles et al. 2009). These factors can lead to an extensive haplotype structure which is more pronounced in self-pollinating crops (Hamblin et al. 2011). In P. vulgaris, the two gene pools diverged, with a reduction in diversity due to bottlenecks, in wild types at 110 000 years BP. The domestication bottleneck within each gene pool, the subsequent differentiation into races and selection by breeders to develop cultivars led to an increase in LD. Finally, as a result of the its self-pollinating nature, it is difficult to break up the LD generated by these factors in common bean (Myles et al. 2009; Hamblin et al. 2011) and mapping resolution would generally be low. A final concern is the high degree of population structure observed in common bean. As has been documented, population structure can result in spurious associations (false positives) between phenotypes and unlinked markers (Knowler et al. 1988; Cardon and Palmer 2003). Even though mixed linear models, which account for population structure and relatedness can minimise the discovery of false positives, this is at the expense of reducing the power to detect true positives (Zhao et al. 2007; Brachi et al. 2010).

To counteract these effects, we propose that common bean association mapping should be performed with populations consisting of individuals from within a single gene pool. First, this will greatly reduce the population structure problem often observed for common bean populations (Kwak and Gepts 2009; Rossi et al. 2009; McClean et al. 2011). This will be most beneficial for mapping in the Andean gene pool because the very low genetic differentiation among the races and the number of markers required would be high due to its low LD decay. Depending on the specific population selected, it might be possible, as with humans, to perform the mapping without correcting for population structure. As for mapping in the MA gene pool, it may even be of further benefit to consider populations derived from the Durango and Jalisco races as a pool and Mesoamerican races as a second pool. Population structure analyses consistently define Durango and Jalisco landraces as a single subpopulation and Mesoamerican genotypes as a second subpopulation (McClean et al. 2011). Due to a significant LD decay distance, fewer markers would be necessary. Given the low diversity within each of these subpopulations, it will be necessary to use a much larger core set of SNP markers to discover those polymorphic within these distinct subpopulations. Soon this should not be a concern, given the low cost of discovering SNPs using massively parallel next generation techniques.