Selection and adaptive introgression guided the complex evolutionary history of the European common bean

Abstract

Domesticated crops have been disseminated by humans over vast geographic areas. Common bean (Phaseolus vulgaris L.) was introduced in Europe after 1492. Here, by combining whole-genome profiling, metabolic fingerprinting and phenotypic characterisation, we show that the first common bean cultigens successfully introduced into Europe were of Andean origin, after Francisco Pizarro's expedition to northern Peru in 1529. We reveal that hybridisation, selection and recombination have shaped the genomic diversity of the European common bean in parallel with political constraints. There is clear evidence of adaptive introgression into the Mesoamerican-derived European genotypes, with 44 Andean introgressed genomic segments shared by more than 90% of European accessions and distributed across all chromosomes except PvChr11. Genomic scans for signatures of selection highlight the role of genes relevant to flowering and environmental adaptation, suggesting that introgression has been crucial for the dissemination of this tropical crop to the temperate regions of Europe.

Fig. 1. Population structure of common bean in America and Europe.

a Admixture analysis (K = 2) showing inferred ancestry in the American (AM; left) and European (EU; right) accessions, with the identification of two gene pools (identified as clusters 1 and 2) that show correspondence to the two main common bean gene pools based on our passport data (cluster 1 – Andean, cluster 2 – Mesoamerican), and several intermediates and admixed genotypes in Europe. The chloroplast ancestry assignment is shown for the accessions (triangles below the chart) when not consistent with the nuclear assignment b Admixture plots for the AM Mesoamerican accessions (K = 2) grouped by geographic origin (i.e., latitude and state), which identifies two main subgroups (M1 and M2). c Admixture plots for the AM Andean accessions (K = 4) grouped by geographic origin (i.e., latitude and state), which identifies three Andean genetic subgroups (A1, A2 and A3). A fourth cluster in four accessions, based on the whole-set ADMIXTURE analysis (K = 2), was induced by the occurrence of Mesoamerican alleles with AM_M1/AM_M2 components (see also supplementary Note 15). d Neighbour-joining tree and seed pictures of the 66 pure American accessions. e Spatial interpolation of the geographic distributions of the EU Mesoamerican (M1 and M2) and EU Andean (A1, A2 and A3) ancestry components in Europe, as inferred by ChromoPainter analysis. Maps were designed using the map tools implemented in different R packages, such as spatial, maps, fields, maptools, raster, rgdal. f Primary and secondary domestications of Mesoamerican and Andean genetic groups/races in America. Loss of photoperiod sensitivity during secondary domestication was a relevant factor for the introduction of the Andean A1/Nueva Granada and A3/Chile and for the Mesoamerican M2/Mesoamerica ancestries in Europe (solid arrow). Genetic group M1 (Durango-Jalisco race) was successfully introduced into Europe after introgression from other genetic groups characterised by little or no photoperiod sensitivity (dashed arrow). Genetic group A2 (Peru race) was not introduced into Europe due to its high photoperiod sensitivity (discontinuous and truncated line). Source data are provided as a Source Data file.

Fig. 2. Genetic structure, molecular phenotyping and flowering data.

a, b Molecular phenotypes (PCA1 from 1493 putative secondary metabolites, H² > 0.65 over the entire dataset) of 94 American accessions and 96 European accessions confirm the subdivision into the two main groups based on the admixture coefficient (derived from nuclear genomic data, K = 2). Intermediate phenotypes and genotypes are seen in Europe. c Violin plots showing the distribution of PCA1 values related to secondary metabolites showing high heritability (H² > 0.65) by genetic subgroups in the American and European accessions. PCA1 was used as a representative molecular phenotype, and it explains 25.7% of the total variance for these traits. N. biologically independent samples, AM_A1 (19), AM_A2 (15), AM_A3 (8), AM_M1 (20), AM_M2 (31), EU_A1 (39), EU_A3 (17), EU_M1 (31), EU_M2 (4), EU_MIX (5). d Violin plots showing the distribution of PCA1 values related to the days to flowering (DTF) and photoperiod sensitivity (PS) by genetic subgroups in the American and European accessions. PCA1 was used as a representative phenotypic trait for DTF and PS, and it explains 68.8% of the total variance for these traits. N. biologically independent samples, AM_A1 (20), AM_A2 (18), AM_A3 (8), AM_M1 (22), AM_M2 (31), EU_A1 (40), EU_A3 (18), EU_M1 (33), EU_M2 (4), EU_MIX (5). e Proportions of the genetic memberships – P(AM_A1), P(AM_A2), P(AM_A3), P(AM_M1), P(AM_M2), P(SAND), and P(SMES) – inferred from the donor accessions and composing the American and European accessions (grouped as mainly AM_A1, AM_A2, AM_A3_AM_M1, AM_M2, EU_A1, EU_A3, EU_M1, EU_M2, and EU_MIX) are shown in the pie charts below the corresponding groups and flowering data (number and percentage of individuals with delayed or no flowering) in northern and southern Europe, related to the corresponding groups. c, d box plots represent minimum, first quartile, median, third quartile and maximum. Source data are provided as a Source Data file.

Fig. 3. Mapping introgression in the European common bean using ChromoPainter.

a Proportion of introgressed genome in the Mesoamerican (EU_MES; n = 43) and Andean (EU_AND; n = 71) groups. b, c Boxplots showing the median length of the introgressed blocks identified in each of the EU_AND and EU_MES accessions across all of the chromosomes and the median length of the introgressed blocks identified in each of the EU_AND and EU_MES individuals by chromosome. b, c box plots represent minimum, first quartile, median, third quartile and maximum. Sample size (N. accessions), (b, c; EU_AND = 71, EU_MES = 43). Source data are provided as a Source Data file.

Fig. 4. Boxplots of θπ averaged over 100 kb non-overlapping sliding windows, linkage disequilibrium decay and inter-chromosomal linkage disequilibrium and genetic load.

a Genetic diversity computed using whole chromosomes and the unmasked dataset. b Genetic diversity computed after the admixture masking process using whole chromosomes and linkage disequilibrium (LD) decay according to the physical distance. c Comparative LD decay in the American and European accessions. d Genome-wide measure of genetic load in the American and European accessions. The ratios are shown for missense (up) and loss-of-function (down) over synonymous mutations in the different groups. AM_M* and AM_A* are the admixed American accessions (not pure American individuals). e Private inter-chromosomal LD in American and European accessions (left), in the Mesoamerican and Andean European accessions (middle), and considering genomic regions under selection (S) in the Mesoamerican and Andean European accessions (right) a,b, box plots represent minimum, first quartile, median, third quartile and maximum. a, b, N = 5134 genomic windows. Source data are provided as a Source Data file.

Fig. 5. Candidate genes for adaptation.

Schematic representation of the regulatory networks underlying the four major flowering pathways in Arabidopsis thaliana. The genes involved in the photoperiod, vernalisation, autonomous and gibberellin pathways that lead to the transition from vegetative to flowering are shown below the corresponding pathway. Additional genes belonging to secondary pathways that interact with the main regulatory flowering networks are shown in italic. Orthologues genes in common bean showing signatures of adaptive introgression, and those located in GWA peaks in our study are highlighted as follows: yellow hexagons – common bean orthologues of LHY (Phvul.009G259400, Phvul.009G259650) and VRN1 and RTV1 (Phvul.011G050600) showing private inter-chromosomal linkage disequilibrium (LD) in the EU_M pool; pink hexagons – common bean orthologues of LD (Phvul.001G204600), NUC (Phvul.001G204700), CGA1 and GNC (Phvul.003G137100) showing private inter-chromosomal LD in the EU_M pool; red outlines – at least one orthologous gene in common bean showing signature of selection, introgression and with a significant differentiation (F_ST index) between American and European accessions (p < 0.05); orange outlines – at least one orthologous gene in common bean showing a signature of selection with no significant F_ST (p < 0.05); blue asterisks – at least one orthologous gene in common bean showing a signature of introgression; dashed blue outlines – at least one orthologous gene in common bean located within 50 kb centered on a significant GWA peak for days to flowering; dashed green outlines – at least one orthologous gene in common bean located within 50 kb centered on a significant GWA peak for growth habit; arrows – positive regulation of gene expression; truncated arrows – repression of gene expression; solid lines – direct interactions; dashed lines – indirect interactions in A. thaliana. Candidate genes for adaptation or post-domestication of the common bean in Europe, orthologous to those involved in flowering-related pathways, are shown in parentheses: UBP12/13 (Phvul.007G234000); LHY (Phvul.009G259400, Phvul.009G259650); LUX (Phvul.011G062100); PIL5 (Phvul.001G168700); CIB2 (Phvul.008G133600); LRB1 (Phvul.006G109600); DRIP1/2 (Phvul.001G157400, Phvul.007G177500); VRN1, RTV1 (Phvul.011G050600); UBC1/2 (Phvul.003G191900); LD (Phvul.001G204600); TFL1, ATC (Phvul.001G189200); GA2OX4 (Phvul.006G120700); CGA1, GNC (Phvul.003G137100); GAI, RGA1, RGL1, RGL2 (Phvul.001G230500); LMI1 (Phvul.001G184800, Phvul.001G184900); SIC (Phvul.008G182500); CRP (Phvul.008G142400); MYB30 (Phvul.008G041500); NUC (Phvul.001G154800, Phvul.001G204700, Phvul.011G074100); SUC9 (Phvul.004G085100, Phvul.004G085400, Phvul.004G085594); CYP715A1 (Phvul.007G071500).