A three-population wave-of-advance model for the European early Neolithic

Ancient DNA studies have shown that early farming spread through most of Europe by the range expansion of farmers of Anatolian origin rather than by the conversion to farming of the local hunter-gatherers, and have confirmed that these hunter-gatherers continued to coexist with the incoming farmers. In this short report, I extend a previous three-population wave-of-advance model to accommodate these new findings, and derive the conditions supportive of such a scenario in terms of the relative magnitudes of the parameters. The revised model predicts that the conversion rate must, not surprisingly, be low, but also that the hunter-gatherers must compete more strongly with the converted farmers than with the alien farmers. Moreover, competition with the hunter-gatherers diminishes the speed of the wave-of advance of the farmers. In addition, I briefly consider how the wave-of-advance approach may contribute to interpreting the results of archaeological studies using the summed probability distribution of radiocarbon dates.


Introduction
A longstanding question with regard to the spread of early farming in Europe is whether it occurred by the range expansion of farmers of Near Eastern origin or by the iterative conversion to farming of the local hunter-gatherers. Genetic studies of ancient DNA have produced strong evidence in favor of the former proposal. Thus, the broad pattern appears to be one in which farmers migrated in from Anatolia and their descendants colonized Europe with minimal en route interaction or admixture with the indigenes [1]. More specifically, a principal component analysis shows that the early Neolithic individuals from central and western Europe including Spain, and possibly also southern Scandinavia, cluster together with the Anatolian farmers [2][3][4]. Similarly, Valdiosera et al. [5] estimate the genetic contributions of western hunter-gatherers, early Anatolian farmers, and Yamnaya (a later incursion) by a three-way admixture analysis and show that early Neolithic individuals from the north and south of Spain have greater than 95% Anatolian ancestry, with a small remaining contribution from the western hunter-gatherers. Equally important is the observation that the resident hunter-gatherers continued to coexist with the incoming farmers [6][7][8][9][10].
The seminal theoretical studies on the spread of early farming in Europe by Ammerman and Cavalli-Sforza [11,12] (see also [13,14]) were judiciously equivocal as to the two a1111111111 a1111111111 a1111111111 a1111111111 a1111111111 possibilities. On the one hand, these authors invoked the wave-of-advance model [15] to explain the approximately constant rate of spread of early farming. On the other hand, they argued that the frequency gradients of classical genetic markers currently seen in the southeast to northwest direction were formed during this process, which is more consistent with the second theory, as it may imply a gradual dilution of the Near Eastern genetic contribution. Note that the one-population Fisher model is not equipped to predict the formation of ancestry gradients. The principal components study by Menozzi et al. [16] that appeared to support this latter claim has subsequently been questioned on statistical grounds [17], and Arenas et al. [18] have suggested that these gradients were formed during the Palaeolithic. It has also been shown that the spread of early farming "was not . . . regular across Europe . . . but proceeded in leaps" [19]. However, approximate rate constancy does hold [20], and the truth may lie somewhere in between [21].
The above perspective on the initial phase of the Neolithic transition in Europe-farming was spread by populations of Anatolian descent at an approximately constant rate-needs to be qualified by the recognition of inevitable regional variation (e.g. [9,10,21,22]). Doubts have even been raised on the applicability of the wave-of-advance (reaction-diffusion) approach [23]. Hence, what follows and in particular the relevance of the predictions made in this paper are highly dependent on how these issues are perceived.
Aoki et al. [24] proposed a three-population wave-of-advance model describing the joint dynamics of the incoming Near Eastern farmers ("initial" farmers), resident hunter-gatherers, and hunter-gatherers that have adopted farming ("converted" farmers). This model incorporated conversion of the hunter-gatherers to farming by horizontal cultural transmission and competition between the initial and converted farmers ("intra-subsistence" competition). Within this relatively simple setting, we identified four classes of traveling wave solutions. The solution we judged to be in qualitative agreement with the then accepted observations is reproduced in Fig 1. It shows a wave front of the converted farmers moving at constant speed and that, in its wake, reciprocal gradients form in the densities of the initial and converted farmers. None of the four classes of solutions predicted the steady range expansion of the initial farmers of Near Eastern descent that is coupled with the persistence of the hunter-gatherers, which is what is suggested by the ancient DNA data. This model has subsequently been extended by Mimura and colleagues to allow for differential mobility in hunter-gatherers and farmers and then rigorously analyzed from a mathematical standpoint [25,26]. Ackland et al. [27] have generalized the model by the addition of various interaction terms and a geographical barrier to dispersal; their numerical simulations suggest that the converted farmers may often mediate the spread of farming beyond the barrier.
The main purpose of the present paper is to revise the model of Aoki et al. [24] to be qualitatively consistent with our current understanding of the Neolithic transition in Europe and to address the limited question of the wave profile or composition of the traveling wave solutions, specifically whether the wave front of the spread of farming will be dominated by the initial farmers, and whether the hunter-gatherers will persist in the wake of the wave. More sophisticated models and analyses have been applied to the speed of the wave-of-advance [28,29]. Finally, I suggest that the relevant traveling wave solution may provide a suitable null hypothesis for interpreting the results of archaeological studies using the so-called summed probability distribution of radiocarbon dates, e.g. [30,31].
respectively, at position x and time t in a linear habitat. The revised model is described by the equations (The subscripts t and xx indicate the time derivative and the second spatial derivative, respectively.) The original model [24] contained seven parameters: the diffusion constant D, assumed to be the same for all three populations; the intrinsic growth rates r f , r c , and r h of the initial farmers, converted farmers, and hunter-gatherers; the carrying capacities K and L of the farmers and hunter-gatherers; and the conversion rate e of hunter-gatherers to farming. Here I add three terms, multiplied by parameters w, p, and q, which represent competition between the two types of farmers and the hunter-gatherers ("inter-subsistence" competition). Note the initial and converted farmers are assumed to have the same competitive effect, q, on the hunter-gatherers, and the intra-subsistence competition coefficients have been set to 1. Hence, the model is not completely general. All parameters in the revised model are positive. It is convenient to work with the following non-dimensional equations where F, C, and H have been normalized by their carrying capacities; and a = r f /r c , b = r h /r c , s = eL/r c , g = eK/r h . The meanings of the symbols are summarized in Table 1.
Numerical solutions of Eq 2 are computed on the finite interval 0�x�l with reflecting boundary conditions and from the initial conditions Fðx; Eq 3B entails that before the spread of early farming European hunter-gatherers were at carrying capacity (unity in the non-dimensional variables). An anthropological study based on the ages at death of skeletons in Mesolithic cemeteries [32] suggests that this was approximately the case; farmers far behind the wave front were likely also at carrying capacity (also unity in the non-dimensional variables). Typically, I take l = 800 and d = 80. (The code written in an advanced form of BASIC is available upon request.) Table 1. List of symbols.

Variables (except in Eq 1)
x rescaled position in linear space

Ad hoc analysis
The non-dimensional equations that ignore spatial structure, Eq A1, and their equilibria are given in Appendix 1. The traveling wave solutions of interest are those that connect a stable equilibrium E5 whereĈ < <F to equilibrium E2, as such solutions imply a wave front dominated by the initial farmers and persistence of the hunter-gatherers in its wake. Note and A monostable wave connects a stable equilibrium to an unstable one. It is shown in Appendix 2 that E5 exists and that the necessary conditions for linear stability are met if These conditions can be rewritten in terms of the original parameters by noting s =eL/r c and g = eK/r h . Fig 2 illustrates a monostable wave generated under these conditions. I have been unable to obtain the sufficient conditions for linear stability of E5. An example of a numerically-simulated monostable wave connecting E5 to E2 that is predicted for parameter range p>w+s, w<1, q+g<1. Specifically, a = 1, b = 1, s = 0.01, g = 0.15, w = 0.3, p = 1.1, q = 0.1. In particular, the values of s and g entail that the carrying capacity of farmers is 15 times that of hunter-gatherers [33]. Solid grey, broken black, and solid black curves give the normalized densities of initial farmers F, converted farmers C, and hunter-gatherers H, respectively, on the interval 0�x�800. It is worth noting that F and H are both below carrying capacity behind the wave front. The black arrow marks the wave front. The numerically estimated speed is 1.62, which is in good agreement with the heuristic minimum speed v ¼ 2 ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi að1 À wÞ p = 1.67 (Appendix 3). https://doi.org/10.1371/journal.pone.0233184.g002

PLOS ONE
A bistable wave connects two stable equilibria. Appendix 2 shows that E5 and E2 cannot coexist as stable equilibria. Hence, there are no bistable waves that connect E5 and E2.
In terms of the original parameters, the revised model predicts a minimum speed of either v ¼ 2 ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi Dr f ð1 À wÞ or v ¼ 2 ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi D½r c ð1 À pÞ þ eL� whichever is larger, for a monostable wave connecting E5 to E2, and incidentally also for a monostable wave connecting E4 to E2 (see Fig 3 and Appendix 3). These minimum speeds are slower than those proposed by Aoki et al. [24], because of the presence of the competition terms w and p. In particular, in the special case of r f �r c (a�1), the speed of a monostable wave connecting E5 to E2 is necessarily given by Eq 7A, which entails that the speed, as well as the density (see Eq 5 and Appendix 2), of the invading farmers is monotone decreasing in the competition coefficient w. On the other hand, when Eq 7B applies, the speed is monotone increasing in the conversion rate, e, in agreement with Fort [29].

Discussion
Recent genetic studies of ancient DNA have shown that early farming spread through Europe principally by the demic expansion of farmers of Anatolian origin rather than by cultural diffusion. They have also confirmed that the indigenous hunter-gatherers were not replaced, but continued to coexist for a while with the incoming farmers. In this paper, I extend the threepopulation wave-of-advance model of Aoki et al. [24] to accommodate these new findings. The three populations are the farmers of Anatolian descent (initial farmers, F), the farmers of European descent that have converted from hunting and gathering (converted farmers, C), and the hunter-gatherers that have maintained their lifestyles (H). I formulate the conditions for the scenario suggested by the ancient DNA studies to be realized, in terms of the

PLOS ONE
competition coefficients between the three populations and the non-dimensional conversion rate of hunter-gatherers to farming, which are the relevant parameters of the model. Depending on the parameter values, various traveling wave solutions are possible. Among these, the solution that is qualitatively most compatible with the ancient DNA observations is the one that connects the stable equilibrium E5 (Eq 5) to the unstable equilibrium E2 (Eq 4). Provided certain conditions are met, this solution entails that the wave front will be dominated by the initial farmers and that the hunter-gatherers will persist behind the wave front (Fig 2). The inequalities listed in Eq 6 indicate that such a solution is more likely to exist when the conversion rate (the normalized parameters s = eL/r c and g = eK/r h are proportional to the conversion rate e) is relatively low, the competitive effect of hunter-gatherers on the initial farmers (w) and of the farmers on the hunter-gatherers (q) are relatively low, and the competitive effect of hunter-gatherers on the converted farmers (p) is relatively high. In particular, p > w + s must be satisfied (Eq 6). When this inequality is reversed, the traveling wave solution may connect equilibrium E4 to equilibrium E2, which entails that the spread of farming is mediated by the converted farmers (see Appendix 1 and Fig 3).
Compatibility with the ancient DNA data also requires that the density of converted farmers behind the wave front be relatively low, which is predicted for small values of s, g, w, and q, and large values of p (Appendix 2). From a mathematical standpoint, these conditions are all intuitively reasonable. It is an interesting open question whether the hunter-gatherers (H) would have behaved less competitively toward the phenotypically divergent initial farmers (F) [34,35] than toward the phylogenetically closer converted farmers (C); i.e. w < p. Conversion is modeled here as an instantaneous process, but may actually have occurred gradually [36,37]. Populations in transition may have encroached on the resources used by the pristine hunter-gatherers, eliciting opposition from the latter.
In order that the model described in this paper apply to the Neolithic transition in Europe, the conversion rate of hunter-gatherers to farming must be low. An ethnological study shows that hunter-gatherers may have profitably coexisted with farmers, for example by trading animal protein and labor for carbohydrates and "luxury" items [38]. An archaeological study of central and western Europe shows that farmers of the Linearbandkeramik culture and huntergatherers may have coexisted by spatial exclusion: "People of the LBK settled in exactly those areas only marginally exploited by hunter-gatherers and not . . . with . . . more intense huntergatherers exploitation" [39]. A definitive study combines information from ancient DNA and dietary stable isotope data to show "persuasive evidence for the prolonged coexistence of genetically distinct hunter-gatherers and farming groups over the course of the Neolithic in Central Europe" [7]. Thus, various lines of evidence argue against the ready conversion of hunter-gatherers to farming. Moreover, Bowles [40] compares the productivity of foraging with early farming and concludes that, if hunter-gatherers converted to farming, they did not do so "because cultivation of crops was simply a better way to make a living." The motivation for conversion may have been social competition, perhaps for marriage partners [41]. Prolonged coexistence provides opportunities for conversion, but also entails day-to-day competition unless there is active avoidance.
Fort [29] provides valuable estimates of the conversion rate expressed as "the average number of hunter-gatherers converted by each farmer per generation." These estimates (e.g. 2.34 .7) are calculated from ethnographic records of Christian missionaries in contact with a larger number of hunter-gatherers. Since only "successful" interactions appear to be included, the true conversion rate may be lower. In addition, it is not clear whether relations between missionaries and hunter-gatherers can be regarded as typical of what transpired in prehistoric encounters between farmers and hunter-gatherers.
There was apparently a considerable time delay between first contact of the indigenous hunter-gatherers with the incoming farmers and conversion of the former to farming; conversion may have occurred in the wake of the traveling wave rather than at the wave front. If true, there are at least two implications worth considering. First, it may be interesting to introduce into the model a time delay between contact and conversion, similar to the fixed time interval between migration events assumed by Fort and Méndez [28] and Fort [29]. Second, although conversion and intermarriage are not synonymous, if intermarriage also occurred well behind the wave front, then the "surfing" effect resulting in the introgression of European huntergatherer genes into the farming populations of Anatolian descent might have been considerably weakened [42] (see also Edmonds et al. [43]).
The model proposed in this paper predicts that the farmers of Anatolian descent will invade at below their intrinsic carrying capacity (i.e. the carrying capacity in the absence of competition with the converted farmers and hunter-gatherers, which equals unity in the non-dimensional variables; see Eq 5 and Fig 2) and, depending on the parameter values, well below it (Figs 3 and 4). With the exception of this difference and possibly a slower predicted speed (see Eq 7A), an equally good qualitative description of the wave-of-advance of the incoming farmers is provided by the original Fisher model [15]. Since a low conversion rate of hunter-gatherers to farming is required by the three-population model for it to apply, the motivation for invoking it in lieu of the simpler one-population Fisher model may be questioned. However, the speed of the spread of farming may in fact be slower [20,21] than originally thought (~1 km/yr [11]), which in terms of the present model may perhaps be interpreted as due to the competitive effect, w, of hunter-gatherers on the alien farmers. In addition, if future empirical work can show that the first European farmers were well below their intrinsic carrying capacity because of competition with the hunter-gatherers, the present approach-or at the very least a two-population Lotka-Volterra competition model for the alien farmers and resident huntergatherers-may be justified.
Changes in the relative sizes of populations over time have been investigated using the method of summed probability distributions. This method assumes that the number of reported radiocarbon dates from all known sites in a given area and falling within a given time window can serve as an approximate proxy for the relative population size of that space-time sector [44,45]. The summed probability distribution is a plot of these numbers over the period of interest. Applying this method to the Neolithic transition in Europe, Shennan et al. [30] identified a boom-bust temporal pattern at regional levels, which they tentatively attributed to endogenous causes such as population growth to unsustainable levels. On the other hand, an exponential growth curve apparently gives a good fit to the summed probability distribution for the whole of western Europe. Silva and Vander Linden [31] obtain similar results but interpret them differently.
The traveling wave solution pictured in Fig 2 suggests that the arrival of the wave-ofadvance in any region entails exponential growth followed by demographic stasis of the incoming farmers. In other words, the summed probability distribution is expected to be flat after the initial exponential phase. It is not clear why a bust should follow the boom, but it may be as Shennan et al. [30] suggest that resources are depleted. Note this effect is not incorporated in the current wave-of-advance model, which does not address events subsequent to the initial spread of farming. At the pan-European level, the summed probability distribution should be monotone increasing as indeed it is, but I suggest that the null hypothesis for dependence on time should be closer to quadratic than exponential. This is because, given a constant radial rate of expansion in a planar habitat-as opposed to the linear habitat assumed in the present model-the total surface area occupied by the farmers should increase at a quadratic rate, whereas the density will be uniform throughout the occupied area.
Eq A6A. If a(1−w)<1+s−p on the other hand, a similar argument using Eq A6B shows that v ¼ 2 ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffiffi 1 þ s À p p may be the minimum speed.