Vertical Cultural Transmission Effects on Demic Front Propagation: Theory and Application to the Neolithic Transition in Europe

It is shown that Lotka-Volterra interaction terms are not appropriate to describe vertical cultural transmission. Appropriate interaction terms are derived and used to compute the effect of vertical cultural transmission on demic front propagation. They are also applied to a specific example, the Neolithic transition in Europe. In this example, it is found that the effect of vertical cultural transmission can be important (about 30%). On the other hand, simple models based on differential equations can lead to large errors (above 50%). Further physical, biophysical, and cross-disciplinary applications are outlined.


A. Demic diffusion
The motion of individuals is an efficient mechanism to propagate ideas or other cultural traits.Such a mechanism is called demic diffusion [13].In contrast to demic diffusion, in cultural diffusion the mechanism leading to the spread of ideas is not the motion of individuals but imitation or interbreeding (as explained in detail in the next section).
The earliest demic diffusion models [13,14] were based on the equation [6,15] p(x,y,t + T ) = +∞ −∞ +∞ −∞ p(x + x ,y + y ,t) × φ( x , y )d x d y + R T [p(x,y,t)] − p(x,y,t), (1) where p(x,y,t + T ) is the population density (number of individuals per unit area) at the location (x,y) and time t + T .The time interval T is that between two subsequent dispersal events or jumps, i.e., one generation [6] (defined as the mean age difference between an individual and her or his children).In Eq. (1), which is useful in several anthropological [6,13,14] and biophysical systems [4], the generation time T plays a role analogous to that of a delay time in many physical systems [16].
The dispersal kernel φ( x , y ) appearing in Eq. ( 1) is the probability per unit area that the children of an individual located at (x + x ,y + y ,t) become adults at (x,y,t + T ).The last two terms, R T [p(x,y,t)] − p(x,y,t) correspond to net reproduction (births minus deaths) during T .This is usually well described by the logistic growth [15,17] R T [p(x,y,t)] = p(x,y,t)p max e aT p max + p(x,y,t)(e aT − 1) where a is called the initial growth rate and p max is the carrying capacity.
For isotropic kernels, Eq. ( 1) can be Taylor expanded up to second order in space and first order in time to yield [6] where is the diffusion coefficient and F (p) is the time derivative of p(x,y,t) due to net reproduction, i.e., R T [p(x,y,t)] − p(x,y,t) = T F + T 2 2! ∂F ∂t + • • • [6,15].Equation ( 3) is called Fisher's equation [13,18,19] if a logistic growth function is assumed, which implies that, in the absence of diffusive effects [i.e., without the first term in the right-hand side of Eq. (3)], the population density p(x,y,t) evolves logistically [i.e., according to Eq. ( 2)].
A major drawback of Eqs. ( 1)-( 3) is the following.According to Eq. ( 1), newborn individuals can appear at (x,y) (last two terms), while their parents migrate away from (x,y) (first term on the right-hand side).In other words, in such models, some parents leave their newborn children alone.However, newborn humans can not survive away from their parents.This inconsistency can be avoided by using so-called cohabitation models, which are based on the following equation instead of Eqs.(1) or (3) [20]: where R T [p(x,y,t)] is given by Eq. (2).Note that, according to Eq. ( 6), children can not appear away from their parents' location, so that the cohabitation Eq. ( 6) is more reasonable than the noncohabitation equations ( 1)-( 3), as discussed in detail elsewhere [20].
In order to derive analytical front speeds, it is sufficient to consider the leading edge of the front.There, the population density p(x,y,t) is low, so Eq. ( 2) can be linearized into R T [p(x,y,t)] p(x,y,t)e aT , and Eq. ( 6) becomes where is called the net fecundity or reproductive rate (this is again a net rate, i.e., it includes the effect of mortality).In some papers [21,22], instead of the logistic growth (2), one assumes a linear cutoff, i.e., Eq. ( 7) unless it yields p(x,y,t + T ) > p max , and p(x,y,t + T ) = p max otherwise.Both approaches yield the same front speed [23], but the linear cutoff is simpler when including cultural transmission (see Sec. II) [24].

B. Cultural diffusion
An interesting application of the demic diffusion models summarized above is the Neolithic transition in the Near East and Europe.This is defined as the transition from huntergathering (Mesolithic) to agricultural (Neolithic) economics, which took place first in the Near East (about 10 000 years ago) and later spread westward and northward across Europe (during about 5000 years).Equations ( 3) and (7) have been applied to the Neolithic transition in Refs.[13] and [21], respectively (for some recent reviews on these and other models, see Refs.[9,11]).Such models are typically based on the assumption that the Neolithic transition was caused by an invasion of farmers, i.e., by demic diffusion.However, other possible mechanisms could have played a role in the propagation speed of the Neolithic transition in addition to demic diffusion, such as vertical cultural transmission (interbreeding between farmers and hunter gatherers), horizontal cultural transmission (imitation of the behavior of farmers by hunter gatherers, without mating with farmers), and local innovation (transition of hunter gatherers into farmers without any influence from farmers) [13].A crucial, long-standing problem is that (in contrast to demic diffusion) well-established mathematical models of the latter three mechanisms still do not exist, despite 40 years of research on reaction-diffusion models of the Neolithic transition [14].
There is an extensive literature on mathematical models of cultural transmission (see Refs. [25,26] for two seminal works).However, cultural transmission has been included only in a few mathematical models to compute its effect on the Neolithic front speed [22,[27][28][29][30][31], and always on the basis on Lotka-Volterra-type models.These are of two kinds.
(i) Sets of differential equations [27][28][29][30][31].These can be derived by extending the integrodifference Eq. ( 1) to two interacting populations y + y ,t) φ P ( x , y )d x d y + R T P [p P (x,y,t)] − p P (x,y,T ) + I P , (9) where the subscript N stands for Neolithic farmers and the subscript P for Mesolithic (i.e., late Palaeolithic) hunter gatherers.For simplicity, we have assumed the same generation time T for both populations, which is a good approximation according to anthropological data [32].The terms I N and I P correspond to the interaction between both populations (they depend on the population densities p N and p P , and will be carefully discussed below).As in the derivation above of Eq. ( 3), if the kernels φ N ( x , y ) and φ P ( x , y ) are isotropic, Eqs. ( 9) can be also Taylor expanded up to second order in space and first order in time to yield which generalize Fisher's equation (3) to two interacting populations.We will refer to Eqs. (10) as the differential-equations (DE) model.Classical diffusive Lotka-Volterra systems are a special case of this model.They are based on the assumption that the interaction between both populations [last terms in Eqs.(10)] is properly described by [17] I N = p N p P , I P = − p N p P .(11) Equations ( 10) and ( 11) and similar Lotka-Volterra differentialequations models have been applied to the Neolithic transition in several papers [27][28][29][30][31].However, in Sec.II, it is shown that, for vertical cultural transmission, more appropriate terms than Eq. ( 11) can be derived.
(ii) Sets of cohabitation equations for interacting species have been also applied to the Neolithic transition [22] p N (x,y,t y + y ,t) + I P (x + x ,y + y ,t)] φ P ( x , y )d x d y , (12) which generalize the cohabitation Eq. ( 7) to two interacting populations.If the Lotka-Volterra interaction (11) is assumed, Eqs.(12) become In Sec.II, it is shown that, in fact, the Lotka-Volterra interaction (11) is inappropriate for vertical cultural transmission.
The purpose of this paper is to develop a mathematical model that takes into account vertical cultural transmission (i.e., interbreeding between farmers and hunter gatherers) in addition to demic diffusion.The need for such a model has been also pointed out by several authors [33,34], but it does not exist so far.This paper aims to solve this problem.Before presenting this model, let us briefly describe the other two mechanisms mentioned above.A mathematical model of the third mechanism (horizontal cultural transmission [27,28,31,35]) will not be tackled in this paper because it would require substantial additional calculations since (i) it is well known from empirical data that human individuals often use the frequency of a trait as an indirect evaluation of its merit (biased transmission [26]), and (ii) horizontal transmission takes place at much longer distances than vertical transmission, even for pre-industrial populations [36], and this will lead to a distance-dependent probability of transmission.Finally, the fourth mechanism (local innovation) will not be analyzed here either, for the following reasons: (i) such local processes are probably not amenable to mathematical modeling, and (ii) a gradual front propagation in the Southeast-Northwest direction, with fairly constant speed, is observed for the Neolithic transition in Europe at the continental scale according to the archaeological data [37,38], thus, local innovations probably did not have an important effect on the global spread rate.
The plan of this paper is as follows.The derivation of the interaction terms (appropriate to vertical cultural transmission) is presented in Sec.II, considering homogeneous systems for simplicity.Section III contains the extension of the results to nonhomogeneous systems, and the analysis of their effect on the speed of front propagation.Section IV discusses the implications for the Neolithic transition in Europe.Section V is devoted to concluding remarks and further applications of these interaction terms.

II. VERTICAL CULTURAL TRANSMISSION
In this section, we consider homogeneous systems for simplicity, i.e., the population densities per unit area (p N and p P ) are independent of position.Nonhomogeneous extensions will be analyzed in the next section.Let P N (t) = p N (t)S and P P (t) = p P (t)S stand for the numbers of neolithic farmers (N ) and of mesolithic hunter gatherers (P ), respectively, in a system with area S at time t.
We assume that when a farmer mates a hunter gatherer, their children are all farmers.This is reasonable from the perspective that, if agriculture appeared at all, it was because some hunter gatherers became farmers, i.e., there was a tendency of hunter gatherers to become farmers.It is also in agreement with anthropological observations [13,39].However, it is very easy to generalize the model below to the case in which not all children from cross-matings are farmers [40].
We define I N 0 as the number of new farmers (N ) per generation due to vertical cultural transmission, and I P 0 as the corresponding decrease in the number of hunter gatherers (P ) per generation.Because in this section we neglect the dispersion of the populations for simplicity, the evolution equations are Note that the cohabitation equations ( 12) without dispersal [i.e., with φ N ( x , y ) = δ(0) and φ P ( x , y ) = δ(0) Dirac deltas centered at ≡ √ y =0] reduce to Eqs. ( 14) [41].In order to derive front speeds (next section), first we need to derive the dependence of I N and I P on the population numbers P N and P P for vertical cultural transmission.In the following two sections, we present two different derivations, both yielding the same results.However, first it will be useful to write down the conditions that we expect I N and I P to satisfy.
As explained in the previous section, by definition, vertical transmission is due to interbreeding between populations N and P [25].Intuitively, we expect that the interaction term I N should satisfy the following conditions (I P should also satisfy them, but at this point it is enough to consider I N ).
(i) I N should vanish if P N = 0 and/or P P = 0.
(ii) If P N is increased but P P remains the same, then the number of cross matings (and thus I N ) should become larger [unless P N P P ; see condition (iii)].(iii) In the case P N P P , it is reasonable to expect that, on average, P individuals will have reached their maximum possible number of social contacts (encounters per unit time, personal relationships, etc.) with N individuals.Thus, if P N is further increased but P P remains the same, we expect no additional cross matings, i.e., I N should remain approximately the same.
(iv) Similar to condition (ii), if P P is increased but P N remains the same, then I N should become larger [unless P P P N , see condition (v)].(v) Similar to condition (iii), in the case P P P N , if P P is further increased but P N remains the same, I N should remain approximately the same.

A. Phenomenological approach
First, we present a phenomenological or macroscopic approach that will very quickly lead us to the correct dependence of I N and I P on P N and P P .A more rigorous approach (which yields the same results, and also a uniqueness theorem) will be presented in the next sections.
Let us discuss the interaction term I N (the term I P is analogous).The Lotka-Volterra term (11), namely, is widely used for predator-prey and competition systems in ecology [17].However, this term does not satisfy conditions (iii) and (v).Therefore, it does not seem reasonable for vertical cultural transmission.In fact, all power-law interaction terms with the form have the same problem.On the other hand, a conditional term with the form does not satisfy conditions (ii) and (v).However, it is easy to see that the term I N = P N P P P N + P P (18) satisfies all of the five conditions above for the interaction term I N .Thus, this is the term that we will use for cultural vertical transmission.Similarly, we obtain, for I P , A rigorous derivation of Eqs. ( 18) and ( 19) is presented in the next section.Note that the interaction terms ( 18) and ( 19) are similar to the well-known interaction I N = P N P P c+P P , I P = − P N P P c+P P (with c a constant), which is called the Holling type-II interaction [42], Michaelis-Menten kinetics [43], or the quasi-steadystate approximation [44].However, both pairs of interaction terms are substantially different for the purposes of this paper because their front speed is different (see Sec. III).

B. Mating-frequencies approach
Aside from the phenomenological or macroscopic approach discussed in the previous section, a more detailed approach is to use microscopic models, which are based on mating frequencies [25].Here, we show that, in this way, the same result (18) as above can be obtained.This second approach is of utmost interest in its own right because (i) it makes it possible to show that Eq. ( 18) is in fact the only form for I N that satisfies all of the five conditions above (see Sec. II C), and (ii) it makes it possible to estimate the values of and from cross-mating probabilities [p(u) and p (u) below]. Let v = P P P N + P P (21) stand for the frequencies of individuals N and P , respectively.Obviously, Under random mating, marriage can occur with equal probability between any members of the population.In that special case, the probabilities of matings NN, P P , NP , and P N would obviously be u 2 , v 2 , uv, and uv, respectively.In other words, the mating matrix would be [25] Mating between members of two populations with substantially different cultural traits (farmers and hunter gatherers in our case) is obviously not random (e.g., NP and NN matings are not equally likely).Let p(u) stand for the probability that a N individual mates a P individual (i.e., the number of cross matings divided by P N ).Let p (u) stand for the probability that a P individual mates a N individual (i.e., the number of cross matings divided by P P ).Obviously, p(u)P N = p (u)P P . ( Let both probabilities p(u) and p (u) depend on frequency u (and, thus, on v = 1 − u).Indeed, I shall show that they must be frequency dependent in order for the five conditions below Eqs. ( 14) to be satisfied.
As mentioned at the beginning of Sec.II, in agreement with anthropological observations [13,39], we assume that farmers (N ) tend to become hunter gatherers (P ), but not the other way around, in such a way that only children from P P matings remain hunter gatherers (P ).This means that children from NN, NP , and P N matings are farmers (N ) [40].In other words, the probability matrix that a N child results from the mating ij [25] is, in our case, With this information, we can compute the probability that children from a mating are farmers (N ) as We neglect the individuals who do not mate and have children because they will not affect the front behavior (this is also done in mathematical models of vertical cultural transmission [25] because they do not affect the evolution of cultural traits).Then, the total number of matings in the population is simply (P N + P P )/2.By multiplying this by the probability (27) and by the average number of children per mating (namely, 2R 0N , where R 0N is the average number of children per farmer), we obtain that the number of farmers in the next generation is where time is measured in generations and P i stands for P i (t), with i = N,P .
Finally, subtracting to Eq. ( 28) the increase in the number of farmers due to net reproduction, i.e., R 0N P N [first term in the right of Eq. ( 14)], we obtain the increase in the population number of farmers due to vertical cultural transmission, By using Eqs.( 20) and ( 21), we can simplify this result into which simply states that the increase in the number of farmers due to vertical cultural transmission is the number of hunter gatherers who mate with farmers [P P p (u)] times the reproductive rate of farmers (R 0N ).
It is easy to repeat the previous steps to derive the decrease in the number of hunter gatherers.For this purpose, instead of the matrix (26), we have to consider the probability c ij that a P child results from the mating ij , in our case [40] and instead of Eq. ( 27), we now compute the probability that children from a mating are hunter gatherers (P ) as Following the same steps as above, we obtain that the analog to Eqs. ( 28)-( 30) for hunter gatherers instead of farmers is respectively.This last result simply states that the decrease in the number of hunter gatherers due to vertical cultural transmission is the number of hunter gatherers who mate with farmers [P P p (u)] times the reproductive rate of hunter gatherers (R 0P ).Again, this was to be expected intuitively.But, it is very important to stress that, in Eq. ( 35), the fecundity of hunter gatherers (R 0P ) appears, whereas in Eq. ( 30), the fecundity of farmers (R 0N ) appears instead.Therefore, we conclude that, in contrast to many ecological models [e.g., the Lotka-Volterra-type equations ( 10)-( 12)] for vertical cultural transmission, the interaction terms I N and I P in Eqs. ( 14) must be different from each other (not only in sign but also in magnitude).
In order to write the interaction terms I N and I P with explicit dependencies on the population numbers P N and P P , we need first to relate the cross-mating probability p (u) to the N-population frequency u.In order to do so, it will be very useful to consider several cases.
First, consider the special case in which p (u) is frequency independent, say, Then, we can rewrite Eq. ( 30) as which corresponds (for P N = 0) to the case given by the phenomenological equation ( 17) (with = R 0N a 0 ).However, as explained below Eq. ( 17), such an interaction term is not acceptable for vertical cultural transmission [because it is inconsistent with the conditions (ii) and (v)].Thus, the mating probability p (u) must be frequency dependent.Second, consider the case in which the probability p (u) that a P individual mates a N individual is proportional to the frequency of N individuals u, say, Then, making use of Eq. ( 20), the number of farmers in the next generation [Eq.(30)] is This result agrees with that of the phenomenological approach [Eq.( 18)] (with = R 0N η), which satisfies all of the five conditions stated below Eq. ( 14).Note from Eqs. ( 14) and ( 39) that η is dimensionless.On the other hand, use of Eqs.(38) and ( 20) into (24) yields where we have applied Eq. ( 22).Therefore, not only p (u) but also p(u) must be frequency dependent.Equations ( 30) and ( 35) imply that and, by using Eq. ( 39), the interaction term for hunter gatherers due to vertical cultural transmission is finally which also agrees with the phenomenological result (19) (with = R 0P η).We stress that it is very easy to generalize the model above to the case in which not all children from cross matings are farmers [40], but this is not necessary for the purposes of this paper.
In Secs.III and IV, these results will be applied to compute the effect of vertical cultural transmission on the Neolithic front speed.However, before doing so, the next section briefly tackles the question of whether it is possible to find some other forms for the interaction terms I N and I P such that they are also compatible with vertical cultural transmission.

C. Uniqueness of the interaction terms for vertical cultural transmission
In this section, we show that it is not possible to find any other forms for I N and I P aside from those given by Eqs. ( 39) and ( 42) such that they satisfy the five conditions stated below Eq. ( 14) for vertical cultural transmission.In order to do so, let us write the probability p (u) as a Taylor series and note the following points.
(i) The first term in the series (43) can not be accepted because [as explained below Eq. ( 37)], when used in Eq. ( 30), it gives rise to the term R 0N a 0 P P , which makes I N inconsistent with condition (v) below Eq. ( 14) [45].
(ii) The third term in Eq. ( 43) can not be accepted neither because, when used in Eq. ( 30), it will give rise to a term R 0N a 2 P 2 N P P /(P N + P P ) 2 in I N , which obviously makes it impossible for I N to fulfill condition (v) [45].
(iii) Similarly, all other terms a i u i (with i = 1) in Eq. ( 43) give rise to terms R 0N a i P i N P P /(P N + P P ) i in I N , which are also inconsistent with condition (v) [45].
Therefore, Eqs. ( 39) and ( 42) seem to be the only forms for I N and I P suitable to describe vertical cultural transmission.Using them in Eqs. ( 14) yields the nonspatial model of vertical cultural transmission, namely, P N (t)+P P (t) , P P (t + T ) = R 0P P P (t) − R 0P η P N (t)P P (t) P N (t)+P P (t) . ( In the next section, this model will be generalized to spatial systems (in order to take care of population dispersal).
It is very important to stress that terms with the Lotka-Volterra form (15) are widely applied in ecological models of interacting species for the so-called predator-prey (I P = −I N = − P N P p ) and the competition (I P = I N = − P N P p ) cases [17].However, Eqs. ( 39) and ( 42) definitely show that such Lotka-Volterra terms can not be applied to vertical cultural transmission (interbreeding).Indeed, not only the dependencies of I P and I N on the population densities (P N and P p ) are completely different than in the Lotka-Volterra case, but even their magnitudes differ from each other (|I P | = |I N |).We conclude that previous models of the effect of interbreeding on Neolithic front dynamics, which were based on Lotka-Volterra interaction terms [22,29,30], should be replaced by models based on the interaction terms ( 39) and (42) derived above.This is done in the following sections.

III. FRONT PROPAGATION
In order to generalize the nonspatial approach in the previous section, the following two spatial models will be analyzed.
(i) Sets of cohabitation equations.As explained above Eq.( 7), cohabitation equations should be preferred for biological reasons [20].They are given by Eqs. ( 12) with the interaction terms valid for vertical cultural transmission [i.e., using Eqs.( 39) and ( 42) instead of the Lotka-Volterra interaction (11)] p p (x+ x ,y+ y ,t)+p N (x+ x ,y+ y ,t) p P (x + x ,y + y ,t) φ P ( x , y )d x d y .(45) Let us recall that, for the case in which there is no dispersal [φ N ( x , y ) = φ P ( x , y ) = δ(0), i.e., the Dirac function centered at ≡ √ this set of cohabitation equations reduces to the nonspatial model in Sec.II [Eqs.(14)] with I N and I P given by Eqs.(39) and (42).Recall, also, that Eqs.(45) are strictly valid only if they yield population densities below the corresponding saturation values (see the last paragraph in Sec.I A), but here this point will not be discussed further because it would not change the result that we are interested in, namely, the front speed [24].
(ii) Sets of differential equations (DE model).It is included here only for comparison purposes, namely, to determine whether it is a reasonable approximation or not [it should be stressed again that this model (ii) is less reasonable than the cohabitation model (i) for biological reasons, see the text above and below Eq. ( 7) as well as Refs.[20]].This model (ii) is given by the classical diffusive Eqs.(10) with the vertical-transmission interaction terms [i.e., Eqs. ( 39) and ( 42) rather than Lotka-Volterra terms (11)] where the net reproduction rate of species i = N,P is given by the logistic growth function (5), In order to compare the predictions of Eqs. ( 46) and ( 47) to those of Eqs. ( 45), we note that the reproductive parameters appearing in Eqs. ( 47) and ( 45) are related according to Eq. ( 8), namely [46],

A. Front speed in the cohabitation-equations model
This model corresponds to Eqs. (45).The Neolithic transition is widely regarded as an invasion of Neolithic farmers N into a space populated by indigenous hunter gatherers [13].Thus, we can assume that the invasion front of farmers (with density p N ) spreads in a region where the density of hunter gatherers (with density p P ) is initially equal to its maximum possible value (p max P ) [29].Then, in the leading edge of the invasion front, p N is very small whereas p P is close to saturation: where O(2) stands for second-and higher-order terms ε(x,y,t) p max N and δ(x,y,t) p max P .By performing a two-variable Taylor expansion up to first order, it is easy to see that η p p (x,y,t) p N (x,y,t) By using this result, the first of Eqs. ( 45) can be linearized and the dependent variable p P (x,y,t) no longer appears in it: For isotropic kernels, invasion front speeds in two dimensions (2D) can be found by choosing the x axis parallel to the local velocity of the front [47].Let c ≡ |c x | stand for the front speed.Assuming constant-shape solutions with the form p N = p 0 exp[−λ(x − ct)] as x − ct → ∞ and an isotropic kernel φ( ), we obtain, from Eq. ( 51), where θ ≡ tan −1 y x .Finally, we need an expression for the kernel φ( ) in order to integrate Eq. ( 52).Kernels appropriate to Neolithic human dispersal have been recently discussed [15].However, here we want to focus our attention on the effect of vertical cultural transmission.Therefore, we analyze a simple model in which individuals either remain at rest (with probability p e , which is called the persistence [21]) or move a distance r (with probability 1 − p e ): φ N ( ) = p e δ (2) (0) + (1 − p e )δ (2) (r) = p e δ (1) (0) 2π where δ (2) (r) and δ (1) (r) are the 2D and 1D Dirac deltas centered at r, respectively.After integrating Eq. (52), we assume as usual that the minimum speed is that of the front [48] and obtain where is the modified Bessel function of the first kind and order zero.
In the limit in which the invading farmers N do not interbreed at all with the invaded hunter gatherers, there is no vertical cultural transmission (η → 0) but only demic diffusion, and Eq. ( 54) becomes the speed corresponding to the cohabitation Eq. (7) [see Eq. ( 14) in Ref. [21]], as it should [because then the first of Eqs.(45) reduces to Eq. ( 7)].

B. Front speed in the differential-equations model
As explained above Eq.( 7), cohabitation equations should be preferred for biological reasons [20].Therefore, the set of differential equations ( 46) is at most of approximate validity.But, can it really yield an approximately valid description of the effects of vertical cultural transmission on demic front propagation?In order to answer this question, we linearize Eqs. ( 46) and ( 47) and use Eq. ( 50) to obtain, instead of Eq. ( 51), and following the same approach as in the former section, the speed of front solutions to this equation is easily seen to be where we have used Eq. ( 48), and the subscript DE stands for the differential-equations model.Obviously, Eq. ( 57) also follows directly from Fisher's speed 2 √ a N D N if we note that Eq. ( 56) is Fisher's equation with a modified growth rate a N = a N + R 0N η/T .Finally, it is worthwhile to emphasize two points: (i) the cohabitation result ( 54) is more reasonable biologically for human populations than Eq.(57) [see the text above Eq.( 7) and the references in [20]]; and (ii) the simple formula ( 57) is considered here only in order to determine whether it is a useful approximation to the exact result (54).
Note that, for the Holling type-II interaction mentioned at the end of Sec.II B, the front speeds ( 54) and ( 57) would be different because the right-hand side in Eq. ( 50) would be replaced by ηp N (x,y,t)cp max P + O(2).

IV. APPLICATION TO THE NEOLITHIC TRANSITION IN EUROPE
Application of the models above to the Neolithic transition requires parameter values appropriate for pre-industrial populations.Population number data sets in previously empty islands are available [49] and the average value implied is R 0N = 2.2 [21].The generation time is T = 1 generation = 32 yr [50].The mean value of the population persistence is p e = 0.38 [21], which was estimated from the mobility data in Ref. [51].The value r = 50 km is estimated directly from those of the persistence and the mean-squared displacement 2 = 1544 km 2 [52].Finally, the diffusion coefficient is D N = 2 /(4T ) = 12 km 2 /yr.
In Fig. 1, we show the speed predicted by the cohabitation model (full curve) [Eq.(54)] [53].The front speed predicted by Eq. ( 54) increases with increasing values of the dimensionless interaction parameter η, as it was to be expected because, according to Eq. ( 39), higher values of η correspond to a stronger effect of vertical cultural transmission (i.e., to more cross matings between farmers and hunter gatherers per generation).It is remarkable that the observed speed range (hatched region in Fig. 1) severely constrains the possible range for the dimensionless vertical diffusion parameter as η < 0.1.Therefore, an important conclusion is that the observed Neolithic speed range can be used to constrain the possible range of a fundamental interbreeding parameter (η, introduced in this paper).Moreover, according to Fig. 1, vertical cultural transmission has an important effect on the front speed, up to 30% relative to the speed predicted by the purely demic model (i.e., η = 0).Figure 1 also includes the predictions of the approximate DE model (dashed curve) [Eq.( 57)].By comparing the dashed and solid curves in Fig. 1, it is found that the DE model makes an error relative to the cohabitation model that can be very important (up to 56%) if the effect of FIG. 1. Predictions for the Neolithic front speed as a function of the dimensionless vertical cultural transmission parameter η, according to the cohabitation model [Eq.( 54)] (solid curve).For comparison, the results for the approximate, noncohabitation differential-equation model [Eq.(57)] are also shown (dashed curve).The corresponding demic models (circle and square, respectively) are obtained if the effect of vertical cultural transmission (interbreeding) is neglected (η = 0).The hatched region corresponds to the observed speed range for the Neolithic transition in the Near East and Europe, as implied by the archaeological data using calibrated dates (0.6-1.1 km/yr) [37].vertical cultural transmission is included (as compared to an error of only 4% in the purely demic case η = 0).

V. CONCLUDING REMARKS
This paper contains two different derivations of the interaction terms under vertical cultural transmission (Sec.II).Probably the most important result is that, in this case, the interaction terms are not of Lotka-Volterra type [I N = p N p p and I P = − p N p p , see Eq. ( 11)] but have another form derived here [I N = R 0N η p N p P p N +p P and I P = −R 0P η p N p P p N +p P ; see Eqs. (39) and (42)].These interaction terms have also been used to compute the effect of vertical cultural transmission on demic front propagation (Sec.III) and applied to a specific example, the Neolithic transition in Europe (Sec.IV).The effect of vertical cultural transmission can be important (about 30%).Simple models based on differential equations [Eqs.(46)] can lead to large errors [above 50%, as compared to the integrodifference cohabitation model of vertical cultural transmission, Eqs.(45)].
This paper deals with the front speed problem, but the results should be applied in the future to another problem, namely, the shape of genetic clines (i.e., space-dependent gene frequencies).Simulations of genetic clines have been performed using interaction terms of the Lotka-Volterra type [i.e., similar to I N = p N p p and I P = − p N p p ; see Eq. ( 11)] [27,54].However, according to this paper, interbreeding interaction terms of the type derived here should be used instead [i.e., I N = R 0N η p N p P p N +p P and I P = −R 0P η p N p P p N +p P ; see Eqs. ( 39) and (42)].Aside from Lotka-Volterra terms, interaction terms of the type p N p P (p N +p P ) 2 were also applied to simulate vertical cultural transmission in Ref. [54], but such terms were taken from a horizontal transmission model [35] and are not valid for vertical transmission, according to the present paper.
Moreover, a combination of cultural transmission and demic diffusion (as modeled in this paper) is the leading explanation for some genetic clines observed in present European populations [55].This may lead in the future to another application of the work reported here because genetic studies of the Neolithic transition have sometimes led to contradictory results concerning the importance of interbreeding between farmers and hunter gatherers [56][57][58][59].This paper shows that the maximum value of a dimensionless parameter that quantifies the intensity of interbreeding (η, introduced here ) can be estimated directly from the observed value of the front speed (Fig. 1).
The approach reported in this paper could be applied in the future to similar human spread phenomena (e.g., the Austronesian Neolithic expansion [60], preceramic dispersals of maize and root crops into Panama [61], the diffusion of maize to the southwestern United States [62], etc.).
Additionally, in future work, it would be of interest to analyze whether these interaction terms (42) could be useful even for some systems without vertical cultural transmission, but such that it is reasonable to expect that the five conditions stated below Eq. ( 14) could be satisfied, e.g., predatorprey interactions [17], ecological competition (see Ref. [30], Sec.4.2), language competition [63], the spread of diseases [64] and epidemics [17], and virus infection fronts [4].
In addition to biophysical and interdisciplinary applications, the results in this paper may also have purely physical applications because the interaction term given by Eq. ( 18) [i.e., I N = p N p P p N +p P ] exhibits the following interesting saturation effect.If p P is constant, we have I N → p P for p N → ∞.Similarly, if p N is constant, then I N → p N for p P → ∞.Such a saturation effect of a reaction rate that depends on the concentrations of two species (p N and p P ) could be useful in the analysis of front propagation in systems in which the reactive process is locally limited due to impurities, porosity, etc.In other words, assuming, e.g., that p P has a low local value, then if the concentration of the other species (p N ) is increased, we should expect intuitively the reactive rate to stop increasing at some point (this will correspond to the limit I N → p P for p N → ∞).Thus, this saturation effect would correspond to the fact that the reaction is limited in some spatial sites due to the fact that a reactive component has a lower-than-average initial concentration, e.g., because such sites are (partially) filled up with additional, nonreactive substances (impurities).In future work, this could be applied to study the effect of impurities and/or porosity on chemical [65] and combustion [66] front propagation.Even in the absence of impurities and porosity, these interaction terms could be useful.Indeed, a saturation of the reaction rate should be expected intuitively whenever the concentration of a species is increased, if that of the the other one is fixed and substantially smaller.