## Info

Figure 1.3. Typical examples of the two types of waves of pursuit given by wavefront solutions of the predator (v)-prey (u) system (1.3) with negligible dispersal of the prey. The waves move to the left with speed c. (a) Oscillatory approach to the steady state (b, 1 — b), when a > a*. (b) Monotonie approach of (u, v) to (b, 1 — b) when a < a*.

Figure 1.3. Typical examples of the two types of waves of pursuit given by wavefront solutions of the predator (v)-prey (u) system (1.3) with negligible dispersal of the prey. The waves move to the left with speed c. (a) Oscillatory approach to the steady state (b, 1 — b), when a > a*. (b) Monotonie approach of (u, v) to (b, 1 — b) when a < a*.

the p(k; a = 0) curve. Since the local extrema are independent of a, we then have the situation illustrated in the figure. For 0 < a < a* there are 2 negative roots and one positive one. For a = a* the negative roots are equal while for a > a* the negative roots become complex with negative real parts. This latter result is certainly the case for a just greater than a* by continuity arguments. The determination of a* can be carried out analytically. The same conclusions can be derived using the Routh-Hurwitz conditions (see Appendix B, Volume I) but here if we use them it is intuitively less clear.

The existence of a critical a* means that, for a > a*, the wavefront solutions (U, V) of (1.6) with boundary conditions (1.7) approach the steady state (b, 1 — b) in an oscillatory manner while for a < a* they are monotonic. Figure 1.3 illustrates the two types of solution behaviour.

The full predator-prey system (1.3), in which both the predator and prey diffuse, also gives rise to travelling wavefront solutions which can display oscillatory behaviour (Dunbar 1983, 1984). The proof of existence of these waves involves a careful analysis of the phase plane system to show that there is a trajectory, lying in the positive quadrant, which joins the relevant singular points. These waves are sometimes described as 'waves of pursuit and evasion' even though there is little evidence of prey evasion in the solutions in Figure 1.3, since other than quietly reproducing, the prey simply wait to be consumed.

### Convective Predator-Prey Pursuit and Evasion Models

A totally different kind of 'pursuit and evasion' predator-prey system is one in which the prey try to evade the predators and the predators try to catch the prey only if they interact. This results in a basically different kind of spatial interaction. Here, by way of illustration, we briefly describe one possible model, in its one-dimensional form. Let us suppose that the prey (u) and predator (u) can move with speeds ci and c2, respectively, that diffusion plays a negligible role in the dispersal of the populations and that each population obeys its own dynamics with its own steady state or states. Refer now to

Figure 1.4. (a) The prey and predator populations are spatially separate and each satisfies its own dynamics: they do not interact and simply move at their own undisturbed speed q and C2. Each population grows until it is at the steady state (us, vs) determined by its individual dynamics. Note that there is no dispersion so the spatial width of the 'waves' wu and wv remain fixed. (b) When the two populations overlap, the prey put on an extra burst of speed hi vx, hi > 0 to try and get away from the predators while the predators put on an extra spurt of speed, namely, -h2ux, h2 > 0, to pursue them: the motivation for these terms is discussed in the text.

Figure 1.4. (a) The prey and predator populations are spatially separate and each satisfies its own dynamics: they do not interact and simply move at their own undisturbed speed q and C2. Each population grows until it is at the steady state (us, vs) determined by its individual dynamics. Note that there is no dispersion so the spatial width of the 'waves' wu and wv remain fixed. (b) When the two populations overlap, the prey put on an extra burst of speed hi vx, hi > 0 to try and get away from the predators while the predators put on an extra spurt of speed, namely, -h2ux, h2 > 0, to pursue them: the motivation for these terms is discussed in the text.

Figure 1.4 and consider first Figure 1.4(a). Here the populations do not interact and, since there is no diffusive spatial dispersal, the population at any given spatial position simply grows or decays until the whole region is at that population's steady state. The dynamic situation is then as in Figure 1.4(a) with both populations simply moving at their undisturbed speeds c1 and c2 and without spatial dispersion, so the width of the bands remains fixed as u and v tend to their steady states. Now suppose that when the predators overtake the prey, the prey try to evade the predators by moving away from them with an extra burst of speed proportional to the predator gradient. In other words, if the overlap is as in Figure 1.4(b), the prey try to move away from the increasing number of predators. By the same token the predators try to move further into the prey and so move in the direction of increasing prey. At a basic, but nontrivial, level we can model this situation by writing the conservation equations (see Chapter 11, Volume I) to include convective effects as ut -[(C1 + h1Vx)u]x = f(u, v), (1.12)

where f and g represent the population dynamics and h1 and h2 are the positive parameters associated with the retreat and pursuit of the prey and predator as a consequence of the interaction. These are conservation laws for u and v so the terms on the left-hand sides of the equations must be in divergence form. We now motivate the various terms in the equations.

The interaction terms f and g are whatever predator-prey situation we are considering. Typically f (u, 0) represents the prey dynamics where the population simply grows or decays to a nonzero steady state. The effect of the predators is to reduce the size of the prey's steady state, so f (u, 0) > f (u, v > 0). By the same token the steady state generated by g(v, u = 0) is larger than that produced by g(v, 0).

To see what is going on physically with the convective terms, suppose, in (1.12), hi = 0. Then ut - ciux = f(u, v), which simply represents the prey dynamics in a travelling frame moving with speed c1. We see this if we use z = x + c1t and t as the independent variables in which case the equation simply becomes ut = f (u,v). If c2 = ci, the predator equation, with h2 = 0, becomes vt = g(v, u). Thus we have travelling waves of changing populations until they have reached their steady states as in Figure 1.4(a), after which they become travelling (top hat) waves of constant shape.

Consider now the more complex case where h1 and h2 are positive and c1 = c2. Referring to the overlap region in Figure 1.4(b), the effect in (1.12) of the h1vx term, positive because vx > 0, is to increase locally the speed of the wave of the prey to the left. The effect of -h2ux, positive because ux < 0, is to increase the local convection of the predator. The intricate nature of interaction depends on the form of the solutions, specifically ux and vx, the relative size of the parameters c1, c2, h1 and h2 and the interaction dynamics. Because the equations are nonlinear through the convection terms (as well as the dynamics) the possibility exists of shock solutions in which u and v undergo discontinuous jumps; see, for example, Murray (1968, 1970, 1973) and, for a reaction diffusion example, Section 13.5 in Chapter 13 (Volume 1).

Before leaving this topic it is interesting to write the model system (1.12), (1.13) in a different form. Carrying out the differentiation of the left-hand sides, the equation system becomes ut -[(c1 + h1 vx)]ux = f(u,v) + h1 uvxx,

In this form we see that the h1 and h2 terms on the right-hand sides represent cross diffusion, one positive and the other negative. Cross diffusion, which, of course, is only of relevance in multi-species models was defined in Section 11.2, Volume I: it occurs when the diffusion matrix is not strictly diagonal. It is a diffusion-type term in the equation for one species which involves another species. For example, in the u-equation, h1uvxx is like a diffusion term in v, with 'diffusion' coefficient h1u. Typically a cross diffusion would be a term d (Dvx)/dx in the u-equation. The above is an example where cross diffusion arises in a practical modelling problem—it is not common.

The mathematical analysis of systems like (1.12)-(1.14) is a challenging one which is largely undeveloped. Some analytical work has been done by Hasimoto (1974), Yoshi-kawa and Yamaguti (1974), who investigated the situation in which h1 = h2 = 0 and

Murray and Cohen (1983), who studied the system with h1 and h2 nonzero. Hasimoto (1974) obtained analytical solutions to the system (1.12) and (1.13), where hi = h2 = 0 and with the special forms f (u,v) = l1uv, g(u,v) = l2uv, where l1 and l2 are constants. He showed how blow-up can occur in certain circumstances. Interesting new solution behaviour is likely for general systems of the type (1.12)-(1.14).

Two-dimensional problems involving convective pursuit and evasion are of ecological significance and are particularly challenging; they have not been investigated. For example, in the first edition of this book, it was hypothesized that it would be very interesting to try and model a predator-prey situation in which species territory is specifically involved. With the wolf-moose predator-prey situation in Canada we suggested that it should be possible to build into a model the effect of wolf territory boundaries to see if the territorial 'no man's land' provides a partial safe haven for the prey. The intuitive reasoning for this speculation is that there is less tendency for the wolves to stray into the neighbouring territory. There seems to be some evidence that moose do travel along wolf territory boundaries. A study along these lines has been done and will be discussed in detail in Chapter 14.

A related class of wave phenomena occurs when convection is coupled with kinetics, such as occurs in biochemical ion exchange in fixed columns. The case of a single-reaction kinetics equation coupled to the convection process, was investigated in detail by Goldstein and Murray (1959). Interesting shock wave solutions evolve from smooth initial data. The mathematical techniques developed there are of direct relevance to the above problems. When several ion exchanges are occurring at the same time in this convective situation we then have chromatography, a powerful analytical technique in biochemistry.

1.3 Competition Model for the Spatial Spread of the Grey Squirrel in Britain

### Introduction and Some Facts

About the beginning of the 20th century North American grey squirrels (Sciurus caro-linensis) were released from various sites in Britain, the most important of which was in the southeast. Since then the grey squirrel has successfully spread through much of Britain as far north as the Scottish Lowlands and at the same time the indigenous red squirrel Sciurus vulgaris has disappeared from these localities.

Lloyd (1983) noted that the influx of the grey squirrel into areas previously occupied by the red squirrel usually coincided with a decline and subsequent disappearance of the red squirrel after only a few years of overlap in distribution.

The squirrel distribution records in Britain seem to indicate a definite negative effect of the greys on the reds (Williamson 1996). MacKinnon (1978) gave some reasons why competition would be the most likely among three hypotheses which had been made (Reynolds 1985), namely, competition with the grey squirrel, environmental changes that reduced red squirrel populations independent of the grey squirrel and diseases, such as 'squirrel flu' passed on to the red squirrels. These are not mutually exclusive of course.

Prior to the introduction of the grey, the red squirrel had evolved without any interspecific competition and so selection favoured modest levels of reproduction with low numerical wastage. The grey squirrel, on the other hand, evolved within the context of strong interspecific competition with the American red squirrel and fox squirrel and so selection favoured overbreeding. Both red and grey squirrels can breed twice a year but the smaller red squirrels rarely have more than two or three offspring per litter, whereas grey squirrels frequently have litters of four or five (Barkalow 1967).

In North America the red and grey squirrels occupy separate niches that rarely overlap: the grey favour mixed hardwood forests while the red favour northern conifer forests. On the other hand, in Britain the native red squirrel must have evolved, in the absence of the grey squirrel, in such a way that it adapted to live in hardwood forests as well as coniferous forests. Work by Holm (1987) also tends to support the hypothesis that grey squirrels may be at a competitive advantage in deciduous woodland areas where the native red squirrel has mostly been replaced by the grey. Also the North American grey squirrel is a large robust squirrel, with roughly twice the body weight of the red squirrel. In separate habitats the two squirrel species show similar social organisation, feeding and ranging ecology but within the same habitat we would expect even greater similarity in their exploitation of resources, and so it seems inevitable that two species of such close similarity could not coexist in sharing the same resources.

In summary it seems reasonable to assume that an interaction between the two species, probably largely through indirect competition for resources, but also with some direct interaction, for example, chasing, has acted in favour of the grey squirrel to drive off the red squirrel mostly from deciduous forests in Britain. Okubo et al. (1989) investigated this displacement of the red squirrel by the grey squirrel and, based on the above, proposed and studied a competiton model. It is their work we follow in this section. They also used the model to simulate the random introduction of grey squirrels into red squirrel areas to show how colonisation might spread. They compared the results of the modelling with the available data.

### Competition Model System

Denote by Si(X, T) and S2(X, T) the population densities at position X and time T of grey and red squirrels respectively. Assuming that they compete for the same food resources, a possible model is the modified competition Lotka-Volterra system with diffusion, (cf. Chapter 5, Volume I), namely,

d S2

where, for i = 1, 2, a, are net birth rates, 1 /b, are carrying capacities, ci are competition coefficients and Di are diffusion coefficients, all non-negative. The interaction (kinetics) terms simply represent logistic growth with competition. For the reasons discussed above we assume that the greys outcompete the reds so

We now want to investigate the possibility of travelling waves of invasion of grey squirrels which drive out the reds. We first nondimensionalise the model system by setting

and (1.15) becomes

d02 2

In the absence of diffusion we analysed this specific competition model system (1.18) in detail in Chapter 5, Volume I. It has three homogeneous steady states which, in the absence of diffusion, by a standard phase plane analysis, are (0, 0) an unstable node, (1,0) a stable node and (0,1) a saddle point. So, with the inclusion of diffusion, by the now usual procedure, there is the possibility of a solution trajectory from (0,1) to (1, 0) and a travelling wave joining these critical points. This corresponds to the ecological situation where the grey squirrels (01) outcompete the reds (02) to extinction: it comes into the category of competitive exclusion (cf. Chapter 5, Volume I).

In one space dimension x = x we look for travelling wave solutions to (1.18) of the form

0i = 0i (z), i = 1,2, z = x — ct, c > 0, (1.20)

where c is the wavespeed. 01(z) and 02 (z) represent wave solutions of constant shape travelling with velocity c in the positive x-direction. With this, equations (1.18) become d201 d01

—2 + c—1 + 01(1 — 01 — Y102) = 0, dz2 dz 1 1 1 2 (1.21)

d202 d02 K—-2 + c—2 + a02(1 — 02 — Y201) = 0, dz2 dz subject to the boundary conditions

01 = 1, 02 = 0, at z =— c», 01 = 0,02 = 1, at z = (1.22)

That is, asymptotically the grey (01) squirrels drive out the red (02) squirrels as the wave propagates with speed c, which we still have to determine.

Hosono (1988) investigated the existence of travelling waves for the system (1.18) with (1.19) and (1.22) under certain conditions on the values of the parameters. In general, the system of ordinary differential equations (1.18) cannot be solved analytically. However, in the special case where k = a = 1,y\ + y2 = 2 we can get some analytical results. We add the two equations in (1.21) to get d2e de

dz2 dz which is the well-known Fisher-Kolmogoroff equation discussed in depth in Chapter 13, Volume I which we know has travelling wave solutions with appropriate boundary conditions at However, the boundary conditions here are different to those for the classical Fisher-Kolmogoroff equation: they are, from (1.22), e = 1, at z = (1.24)

which suggest that for all z, e = 1 ^ ex + e2 = 1. (1.25)

Substituting this into the first of (1.21) we get d2ei dei

dz2 dz which is again the Fisher-Kolmogoroff equation for e1 with boundary conditions (1.22). From the results on the wave speed we deduce that the wavefront speed for the grey squirrels will be greater than or equal to the minimum Fisher-Kolmogoroff wave speed for (1.26); that is, c > cmin = 2(1 — Y1)1/2, Y1 < 1. (1.27)

Similarly, from the second of (1.21) with (1.25) the equation for e2 is d2e2 de2

dz2 dz with boundary conditions (1.22). This gives the result that for the red squirrels c > cmin = 2(Y2 — 1)1/2, Y2 > 1. (1.29)

Since y1 + Y2 = 2 (and remember too that k = a = 1) these two minimum wavespeeds are equal. In terms of dimensional quantities, we thus get the dimensional minimum wavespeed, Cmjn, as

Parameter Estimation

We must now relate the analysis to the real world competition situation that obtains in Britain. The travelling wavespeeds depend upon the parameters in the model system (1.15) so we need estimates for the parameters in order to compare the theoretical wave-speed with available data. As reiterated many times, this is a crucial aspect of realistic modelling.

Let us first consider the intrinsic net growth rates ai and a2. Okubo et al. (1989) used a modified Leslie matrix described in detail by Williamson and Brown (1986). In principle the estimates should be those at zero population density, but demographic data usually refer to populations near their equilibrium density. Three components are considered in estimating the intrinsic net growth rate, specifically the sex ratio, the birth rate and the death rate. The sex ratio is taken to be one to one. Determining the birth and death rates, however, is not easy. It depends on such things as litter size and frequency and their dependence on age, age distribution, where and when the data are collected, food source levels and life expectancy; the paper by Okubo et al. (1989) shows what is involved. After a careful analysis of the numerous, sometimes conflicting, sources they estimated the intrinsic birth rate for the grey squirrels as a1 = 0.82/year with a stable age distribution of nearly three young to one adult and for the red squirrels, a2 = 0.61/year with an age distribution of just over two young for each adult.

Determining estimates for the carrying capacities 1/b1 and 1/b2 involves a similar detailed examination of the available literature which Okubo et al. (1989) also did. They suggested values for the carrying capacities of 1/b1 = 10/hectare and 1/b2 = 0.75/hectare respectively for the grey and red squirrels.

Unfortunately there is no quantitative information on the competition coefficients c1 and c2. In the model, however, only the ratios c1/b2 = y1 and c2/b1 = y2 are needed to estimate the minimum speed of the travelling waves. As far as the speed of propagation of the grey squirrel is concerned, we only need an estimate of y1: recall that 0 < Y1 < 1. Since y1 appears in the expression of the minimum wavespeed in the term (1 — y1)1/2, the speed is not very sensitive to the value of y1 if it is small, in fact unless it is larger than around 0.6. We expect that the competition coefficient c1, that is, red against grey, should have a small value. So, this, together with the smallness of the carrying capacity b—1, it is reasonable to assume that the value of y1 is close to zero, so that the minimum speed of the travelling wave of the grey squirrel, Cmin, is approximately given from (1.30) by 2(^1a1)1/2. In the numerical simulations carried out by Okubo et al. (1989) they used several different values for the y s since the analysis we carried out above was for special values which allowed us to do some analysis.

Let us now consider the diffusion coefficients, D1 and D2. These are crucial parameters in wave propagation and notoriously difficult to estimate. (The same problem of diffusion estimation comes up again later in the book when we discuss the spatial spread of rabies in a fox population, bacterial patterns and tumour cells in the brain.) Direct observation of dispersal is difficult and usually short term. The reported values for movement vary widely. There is also the movement between woodlands.

For grey squirrels, a maximum for a one-dimensional diffusion coefficient of 1.25 km2/yr, and for a two-dimensional diffusion coefficient of 0.63 km2/yr, was derived based on individual movement. However, this may not correspond to the squirrels'

 l (km)
0 0