Demographic history of european populations of Arabidopsis thaliana

Abstract

The model plant species Arabidopsis thaliana is successful at colonizing land that has recently undergone human-mediated disturbance. To investigate the prehistoric spread of A. thaliana, we applied approximate Bayesian computation and explicit spatial modeling to 76 European accessions sequenced at 876 nuclear loci. We find evidence that a major migration wave occurred from east to west, affecting most of the sampled individuals. The longitudinal gradient appears to result from the plant having spread in Europe from the east approximately 10,000 years ago, with a rate of westward spread of approximately 0.9 km/year. This wave-of-advance model is consistent with a natural colonization from an eastern glacial refugium that overwhelmed ancient western lineages. However, the speed and time frame of the model also suggest that the migration of A. thaliana into Europe may have accompanied the spread of agriculture during the Neolithic transition.

Figure 1. Bayesian clustering.

(A) Membership coefficients in K _max = 5 putative populations, computed using the average values over the 10 TESS runs with the smallest values of the deviance information criterion from a total of 100 runs. Similar results were obtained with other values of K _max from 4 to 10. (B) Interpolated membership coefficients in the three apparent subpopulations: western cluster, eastern cluster, and northern cluster.

Figure 2. Diversity regressed on geographic distance.

Correlation (R) map for the linear regression of expected heterozygosity on great circle distance. We used 300×180 points on a two-dimensional lattice covering Europe, and we computed distances from each lattice point considered as a potential source. The dots represent the centers of the 7 population samples used in the regression analysis.

Figure 3. Bayes factors.

The 4 demographic scenarios (Models A–D) and their associated Bayes factors. Model A is the model with constant population size, N ₀. Model B is a model with an exponentially growing population size (present size, N ₀, ancestral size, N ₁, time since the onset of expansion, t ₀). In Model C, the growth is exponential between two periods with constant size (present size, N ₀, ancestral size, N ₁, time since the onset of expansion, t ₀, time since the end of expansion, t ₁). Model D is similar to Model B, but it includes an ancient bottleneck before expansion. Variants of these 4 models, including variable mutation rates across loci, are considered here. The Bayes factors (top boxes) correspond to the ratio of the weight of evidence of each model to the weight of evidence of Model B. Two window sizes, δ _0.01 and δ _0.05, were used when computing the Bayes factors. These window sizes correspond to the 1% and 5% quantiles of the distance between the values of the summary statistics obtained under Model B and the observed values of the summary statistics. The Bayes factors were identical for the 2 window sizes and for values rounded for one decimal place, except for Model C, for which a minor difference was observed (1.8 for δ _0.05 instead of 1.9).

Figure 4. Onset and duration of the demographic expansion.

Plot of the joint posterior distribution for the time of onset of the expansion, t ₀, and the length of the expansion, t ₀−t ₁. Computations were performed under demographic Model C, in which the population was initially constant, then grew exponentially until t ₁, and then remained constant until the present. Percentages represent the cumulative probabilities under the density curve. The straight line indicates that the duration of expansion cannot be longer than the time elapsed since the onset of expansion.

Figure 5. Number of distinct and private haplotypes.

The mean number of distinct haplotypes and the mean number of private haplotypes of the central European population and the northern European population as functions of sample size. Vertical bars show standard error.

Figure 6. Estimation of the splitting time between the northern and central European populations of A. thaliana.

The mean number of distinct haplotypes and the mean number of private haplotypes of two simulated populations, as functions of sample size. The dark orange lines show the simulation results for a population of size 135,000, and the dark green lines show the simulation results for a population of size 135,000×1/4. The top panel shows the case when the split time is 0. Below follow the results for increasing split times. No migration is assumed. The split time T is given in units of population size. The fit of the simulated data to the observed data was evaluated by the mean across the 100 simulations of the sum of squared differences (SSD) between each simulated data set and the observed data.

Figure 7. Estimation of the migration rate between the northern and central European populations of A. thaliana.

The mean number of distinct haplotypes and the mean number of private haplotypes of two simulated populations as functions of sample size, shown for 100 replicates. The dark orange lines show the simulation results for a population of size N_CE = 135,000, and the dark green lines show the results for a population of size 135,000×1/4, when T = 13,500 years. The top panel shows the case when the migration rate, m, equals 0, and then follow the cases with m = 3 and m = 6 (normalized by N _CE). The results from the observed populations are also plotted for comparison (lighter orange and green lines).

Figure 8. Chi-square statistic maps for spatial range expansion.

(A) χ ² distances between the simulated and the empirical folded frequency spectra as a function of the time of onset of the expansion. The other parameters were fixed at m = 0.25, r = 0.6–1.2, and N ₁ = 10,000. The origin was placed north of the Black Sea (48°N, 35°E). The horizontal line corresponds to the 95% rejection interval of the χ ² test (df = 3, see Methods). (B) Interpolated map of χ ² distances between simulated and empirical folded spectra for 24 potential origins (black dots). The time of onset was fixed at 9,000 years BP, and the other parameters were fixed as in (A).

Figure 9. Frequency spectrum in actual and simulated data.

Minor allele frequency spectra of empirical data and data simulated under the best-fitting model of spatial range expansion. Population growth followed the logistic model within each deme (see text for the other parameter settings). The solid line (grey) corresponds to the neutral folded frequency spectrum. (A) The empirical folded spectrum was computed from the 648 inter-genic and non-coding sequences. (B) The simulated spectrum was computed using the same number of neutral nucleotides as in the data. In simulations, expansion started 9,000 years ago from a potential origin north of the Black Sea (48°N, 35°E). Other locations from a large region around this potential origin yielded very similar simulated spectra.