Printer Friendly


Abstract. A fundamental characteristic of any biological invasion is the speed at which the geographic range of the population expands. This invasion speed is determined by both population growth and dispersal. We construct a discrete-time model for biological invasions that couples matrix population models (for population growth) with integrodifference equations (for dispersal). This model captures the important facts that individuals differ both in their vital rates and in their dispersal abilities, and that these differences are often determined by age, size, or developmental stage. For an important class of these equations, we demonstrate how to calculate the population's asymptotic invasion speed. We also derive formulas for the sensitivity and elasticity of the invasion speed to changes in demographic and dispersal parameters. These results are directly comparable to the familiar sensitivity and elasticity of population growth rate. We present illustrative examples, using published data on two plants: teasel (Dipsacus sylvestris) and Calathea ovandensis. Sensitivity and elasticity of invasion speed is highly correlated with the sensitivity and elasticity of population growth rate in both populations. We also find that, when dispersal contains both long- and short-distance components, it is the long-distance component that governs the invasion speed--even when long-distance dispersal is rare.

Key words: dispersal; integrodifference equations; matrix population models; sensitivity; speed of invasion; stage-structure; traveling waves.


When a species invades new habitat, its population size grows and its range expands. The rate at which the range expands--the invasion speed--is the basic descriptive statistic for invasion dynamics. It plays the same role as does the intrinsic population growth rate ([lambda] or r = In [lambda]) in demography. Like [lambda], it is determined by the environment and by the life cycle of the species. Like [lambda], it characterizes potential dynamics in a simplified situation, and provides a reference point for discussing effects of more complicated scenarios. And, like [lambda], it provides a starting point for management, either to slow the spread of pests or increase the spread of reintroduced species.

Like most useful theories, the theory of biological invasions began by oversimplifying everything, in this case, both demography and dispersal. The earliest studies of invasion speed used reaction--diffusion models of the form

[frac{[partial]n}{[partial]t}] = f(n)n + D[frac{[[partial].sup.2]n}{[partial][x.sup.2]}] (1)

(Fisher 1937, Skellam 1951), where n(x, t) is population density at location x and time t. These models neglect demographic structure, attempting to capture population dynamics in the density-dependent per capita growth rate f(n). They also neglect possible complexities in the dispersal process that are incompatible with the diffusion formulation.

In the absence of Allee effects, single-species reaction--diffusion models like Eq. 1 predict an asymptotically constant rate of spread equal to 2[sqrt{f(0)D}] (Kolmogorov et al. 1937, Aronson and Weinberger 1975, Kametaka 1976, Fife 1979, Bramson 1983). With parameters estimated from field data, they have been applied to a wide variety of organisms, from muskrats (Skellam 1951, Andow et al. 1990 a, b), sea otters (Lubina and Levin 1988), and ceareal leaf beetles (Andow et al. 1990 a, b), to rabies (Kallen et al. 1985) and the plague (Nobel 1974). (See Hengeveld [1989] for other examples.) Reaction--diffusion models have also been extended to include the effects of advection, spatial heterogeneity, Allee effects, and species interactions (Shigesada and Kawasaki 1997, and references therein).

An alternative to diffusion models like Eq. 1 is the use of a probability density function to describe dispersal. This function, often called a "dispersal kernel," gives the probability distribution of the location of an individual as a function of its starting location. This approach allows for a wider range of dispersal mechanisms than diffusion. The diffusion model arises from a Gaussian dispersal kernel, but field measurements of dispersal rarely yield normal distributions (e.g., Wolfenbarger 1946, 1959, 1975, Taylor 1978, Fitt et al. 1987). Instead, leptokurtic distributions seem to be the rule (Okubo 1980, Howe and Westley 1986, Okubo and Levin 1989); such kernels also arise from biologically reasonable models of dispersal and settlement (Broad-bent and Kendall 1953, Beer and Swaine 1977, Okubo and Levin 1989, Neubert et al. 1995). Kot et al. (1996) have emphasized the importance of the shape of the dispersal kernel in the determination of invasion speeds.

The model we explore here is based on a discrete-time dispersal kernel model, an integrodifference equation. There is a growing literature on integrodifference equation models for unstructured populations (Kot and Schaffer 1986). These models were first used in biology to describe changes in gene frequency (Slatkin 1973, Weinberger 1978, 1982, Lui 1982a, b, 1983, 1985, 1986). Since then, they have also been used to model ecological invasions (Kot 1992, Allen et al. 1996, Kot et al. 1996, Veit and Lewis 1996, Sherratt et al. 1997, Clark et al. 1998).

In this paper, we will show how to develop and analyze models that couple the integrodifference equation format for dispersal with matrix models for structured populations. The importance of population structure for demography is well known (e.g., Caswell 1989, 2000, Tuljapurkar and Caswell 1997). Individuals differ in their vital rates and response to the environment, and much of that difference is determined by age, size, or developmental stage. The combination of these individual responses determines population growth rate, the stable stage distribution, reproductive value, sensitivity and elasticity of population growth rate, and many other demographic statistics.

In studies of invasion, all these demographic factors are important, and in addition we recognize that individuals also differ in their dispersal characteristics. Some species have special dispersal stages (seeds, larvae), while the rest of the life cycle is sedentary. Others may be mobile throughout their life cycle, but disperse farther when they are young (e.g., birds) or older (e.g., mammals). Different stages may have different adaptations for dispersal (e.g., seed polymorphism). Even identical propagules produced by different stages may have different dispersal characteristics (e.g., seeds from a tall tree might disperse farther than identical seeds from a short tree).

Our results provide not only the invasion wave speed, but also its perturbation analysis. We will show how to compute the sensitivity and elasticity of wave speed to both demographic are dispersal parameters. This allows us to answer questions like: "By how much must we decrease second-year survivorship in order to stop an invasion?" or "How important is wind vs. bird-mediated dispersal of seeds to the rate of spread?"

In the next section we describe how to assemble the model from its component parts. In the section Analysis, we show how to calculate the invasion wave speed, and how to analyze the sensitivity of wave speed to changes in demographic and dispersal parameters. In Examples, we present two examples. We end with a brief discussion.


Unstructured integrodifference equation models for population growth and spread have been described in detail in a number of places (e.g., Kot et al. 1996). Briefly, the model is composed of two parts. The first is a difference equation that describes population growth at each spatial location:

n(y, t + 1) = b[n(y, t), y]n(y, t) (2)

where n(y, t) is population density at location y and time t, and b[[cdotp]] is the per capita population growth rate, which may depend on both density and location. The second part is an integral operator that accounts for the movement of individuals. We define k(x, y) to be the probability density function for the location x to which an individual at y disperses. To find the population density at x at time t + 1 we simply sum the contributions from all locations y to obtain

n(x, t + 1) = [[[integral of].sup.[infty]].sub.-[infty]] k(x, y)b[n(y, t), y]n(y, t) dy. (3)

The probability density function k is variously known as the "redistribution kernel" or "dispersal kernel." In plant population biology, where the dispersing individuals are often seeds, the kernel is called the "seed shadow."

To add stage structure, we must expand Eq. 2 to a system of m difference equations (one equation for each of the m stages). Designating [n.sub.i] as the population density in the ith stage and [b.sub.ij] as the per capita production of stage i individuals at time t + 1 by stage j individuals at time t, we have

[n.sub.i](y, t + 1) = [[[sum].sup.m].sub.j=1] [b.sub.ij]([n.sub.l](y, t),...,[n.sub.m](y, t), y)[n.sub.j](y, t)

i = 1,..., m (4)

or, in the more familiar matrix notation,

n(y, t + 1) = [B.sub.n](y)n(y, t) (5)

where [B.sub.n](y) is the density-dependent population projection matrix at location y.

To allow for stage-specific dispersal, we must specify a dispersal kernel for each of the [m.sup.2] possible transitions between stages. We define [k.sub.ij](x, y) to be the probability that an individual making the transition from stage j to stage i moves from location y to location x. If there is no dispersal during a given transition, the associated kernel is the Dirac delta function [delta](x - y) (Bracewell 1978), a function which (very roughly speaking) is zero if x [neq] y, is infinite when x = y, and integrates to 1. That is, with probability 1, such an individual stays where it is.

The stage-structured analog of Eq. 3 is then given as follows:

[n.sub.i](x, t + 1)

= [[[integral of].sup.[infty]].sub.-[infty]] [[[sum].sup.m].sub.j=1] [k.sub.ij](x,y)[b.sub.ij]([n.sub.1](y,t),...,[n.sub.m](y,t),y)[n.sub. j](y,t)dy

i = 1,...,m. (6)

Just as Eq. 4 is simplified by Eq. 5, if we create the matrix K(x, y) from the kernels [k.sub.ij], the complex notation of Eq. 6 can be greatly simplified:

n(x, t + 1) = [[[integral of].sup.[infty]].sub.-[infty]][K(x,y)o[B.sub.n](y)]n(y,t) dy. (7)

The symbol "o" stands for the Hadamard product (Horn and Johnson 1985), wherein multiplication is computed element by element. The element in the ith row and jth column of K(x, y) o [B.sub.n](y) is [k.sub.ij](x, y) [b.sub.ij](y).

We now adopt a major simplifying assumption: that the environment is spatially homogeneous. This implies that dispersal from x to y depends only on the relative locations of the two points. As a result, dispersal kernels depend only on the single variable x - y. It also implies that the vital rates depend only on local population density and not explicitly on spatial location. With this extra assumption, Eq. 7 becomes

n(x, t + 1) = [[[integral of].sup.[infty]].sub.-[infty]][K(x-y)o[B.sub.n]]n(y,t) dy. (8)

This assumption is required for much of the analysis of invasion speed that follows below, but we emphasize that Eq. 7 can be used to describe more complicated, heterogeneous situations, and to study them by numerical simulation.

An example

As a simple example, consider a model with two stages, juveniles ([n.sub.1](x, t)) and adults ([n.sub.2](x, t)), as drawn in Fig. 1. Adult fertility is an exponentially decreasing function of local density (measured by total population size [n.sub.1](x, t) + [n.sub.2](x, t)). The other vital rates--juvenile survival ([[sigma].sub.1]), maturation ([gamma]), and adult survival ([[sigma].sub.2])--are density independent. The demographic part of this model is


To describe dispersal we need to decide which stages disperse, and pick a dispersal kernel to describe their movement. Suppose that dispersal occurs only as juveniles become adults, and that the dispersal kernel is a Laplace distribution (back-to-back exponential). The dispersal kernel matrix is then


where [delta](x) is the delta function described above.

In Fig. 2, we show the result of iterating this model forward from an initially small population concentrated at the origin. As the invasion develops, the speed with which it proceeds converges to a constant value (Fig. 3). This is the invasion speed that we want to calculate from the matrices [B.sub.n] and K.


Unstructured models

The unstructured model (Eq. 3) has solutions called traveling waves--solutions that maintain a fixed shape in space and move at a constant speed. Weinberger (1978, 1982), Lui (1982a, b, 1983) and Kot (Kot 1992, Kot et al. 1996) have studied the mathematical properties of these invasion waves. Provided that (1) the growth function b(n)n increases monotonically with population density; that (2) the growth function does not exhibit an Allee effect (i.e., b(n) [leq] b(0) for n [greater than] 0); and that (3) the dispersal kernel possesses a moment-generating function

m(s) = [[[integral of].sup.[infty]].sub.-[infty]] k(x)[]dx (10)

then traveling wave solutions have a minimum speed

[c.sup.*] = [min.sub.s[greater than]0] [frac{1}{s}]ln[b(0)m(s)]. (11)

A population initially concentrated in a finite region of space will never spread faster than [c.sup.*] and asymptotically will spread at a rate of exactly [c.sup.*] (Weinberger 1978, 1982). The consequences of violating some of these conditions have been explored by Weinberger (1982), Lui (1983), Kot et al. (1996) and, in continuous-time models, by Mollison (1972) and Lewis and Kareiva (1993).

Under the assumptions given above, [c.sup.*] depends on the dispersal kernel (via m(s)) but only on the linear part (b(0)) of the demographic model. That is, the wave speed [c.sup.*] can be computed from the linear approximation to Eq. 3 at low densities, i.e., near n = 0. More generally, the linear conjecture (Van den Bosch et al. 1990, Mollison 1991) states that the speed of invasion for any nonlinear model is governed by its linearization at low population densities, as long as there are no Allee effects and no long-distance density dependence. The linear conjecture is supported by many numerical studies. There is also evidence that it holds in higher-dimensional models (Lui 1989a, b; but see Hosono 1995, 1998, for a case in which insidious Allee effects arise in a Lotka-Voterra competition model).

The invasion wave speed [c.sup.*] thus gives the asymptotic invasion speed, for a large and important class of population growth models and dispersal kernels, and for initial conditions corresponding to biologically reasonable invasion scenarios.

Stage-structured models

Speed of invasion.--Relying on the linear conjecture, we now calculate the invasion wave speed for models with stage-structured demography. We make the following assumptions about the demographic and dispersal components of the model:

1) The matrix [B.sub.n] is positive, or nonnegative and primitive (Caswell 1989). In practice most population projection matrices have this property.

2) Let A = [B.sub.0] (i.e., [B.sub.n] evaluated at n = 0). This matrix represents the set of vital rates at low population densities. We will require the largest eigenvalue [lambda] of A to be larger than one. With this assumption we are guaranteed that the population will grow when small.

3) Increased population density has a negative effect on the organism's vital rates:

[B.sub.n]n [leq] An for all n [geq] (12)

where the inequalities are evaluated elementwise. A sufficient condition for the inequalities in Eq. 12 to hold is that all of the elements of [B.sub.n] are nonincreasing functions of the elements of n.

4) All of the kernels have a moment-generating function (cf. Eq. 10); i.e., they have exponentially bounded tails. Models with dispersal kernels that do not satisfy this requirement can produce accelerating invasions (Mollison 1972, Kot et al. 1996).

According to the linear conjecture, the rate of spread of the population modeled by Eq. 8 is governed by its linearization near n = 0:

n(x, t + 1) = [[[integral of].sup.[infty]].sub.-[infty]] [K(x - y)oA]n(y, t) dy. (13)

We look for a rightward traveling wave solution to Eq. 13, i.e., a solution of the form n(x, t) = u([xi]), with [xi] = x - ct and c [greater than] 0. Substitution into Eq. 13 yields

U([xi] - c) = [[[integral of].sup.[infty]].sub.-[infty]] [K([xi] - [eta])*A]u([eta]) d[eta] (14)

where [eta] = y - ct.

Because Eq. 14 is a linear equation, we expect to find solutions that are exponential functions of the new variable [xi], i.e., solutions of the form

u([xi]) = w[e.sup.-s[xi]]. (15)

In this solution, w is a vector that gives the relative abundance of the different stages in the traveling wave, and s determines the shape of the advancing edge of the wave. The problem is to calculate s and c.

We begin by substituting Eq. 15 into Eq. 14 to obtain

w[e.sup.-s[xi]] [] = [A*[[[integral of].sup.[infty]].sub.-[infty]] K([xi] - [eta])[e.sup.-s[eta]] d[eta]]w. (16)

Since we are looking for a rightward traveling wave, we will restrict our attention to s [greater than] 0, so that the solution decays with x.

With the change of variables

[zeta] = [xi] - [eta] (17)

Eq. 16 becomes

w[] = [A*[[[integral of].sup.[infty]].sub.-[infty]] K([zeta])[e.sup.s[zeta]] d[zeta]]w. (18)

The right-hand side of this equation contains a matrix of the moment-generating functions of the dispersal kernels (cf. Eq. 10). We will call this matrix M(s),

M(s) = [[[integral of].sup.[infty]].sub.-[infty]] K([zeta])[e.sup.s[zeta]] d[zeta] (19)

With elements [m.sub.ij](s). Typically, the moment-generating function of a kernel only exists for some finite interval around s = 0. We will therefore assume that there exists some [hat{s}] such that all of the elements of M(s) exist for all 0 [leq] s [less than] [hat{s}], and we will restrict s to this range.

In terms of M(s) we can rewrite Eq. 18 as

[]w = [A*M(s)]w (20)

= H(s)w. (21)

Thus [] is an eigenvalue and w is an eigenvector of H(s). Since H(s) has dimension m X m, there will be m of these eigenvalues, which we denote [[rho].sub.1](s), ..., [[rho].sub.m](s), in decreasing order of magnitude. (Note that H(s) includes both the stage-structured demography and dispersal.)

We show in Appendix A that the asymptotic wave speed of an invasion with initial condition

n(x, 0) = [n.sub.0][e.sup.-sx] (22)

is determined by the largest of the eigenvalues of H(s):

c(s) = [frac{1}{s}]1n [[rho].sub.1](s). (23)

Real invasions, of course, do not begin with initial conditions like Eq. 22. Instead, they have compact support--they are finite in size and restricted to a finite range in space. In Appendix A, we show that the minimum value of c(s) from Eq. 23,

[c.sup.*] = [min.sub.0[less than]s[less than][hat{s}]] [[frac{1}{s}]ln [[rho].sub.1](s)] (24)

is an upper bound on the invasion wave speed for a model governed by Eq. 8 with initial conditions with compact support. We will use the symbol [s.sup.*] for the value of the spatial decay rate s that corresponds to the minimum speed, [c.sup.*].

The linear conjecture and our numerical simulations strongly suggest that [c.sup.*] is not only the upper bound on the wave speed, but also the actual asymptotic wave speed. For example, consider the model (Eq. 9a, b) that we assembled in the Introduction. For this model,



(The moment-generating function of a delta function is the constant 1, and the moment-generating function of the Laplace distribution is 1/(1 - [[alpha].sup.2][s.sup.2]).) Fig. 4 is a plot of c(s) vs. s as calculated from the dominant eigenvalue of H(s) (using Eq. 23). The minimum wave speed obtained from this calculation ([c.sup.*] [approx] 0.208) agrees with the speed we calculated from our simulations (cf. Fig. 3). Fig. 5 shows that [s.sup.*] is the rate at which the population density decays with x as t becomes large.

Sensitivity of invasion speed.--Demographic analysis has been greatly extended by the application of perturbation theory, which permits calculation of the sensitivity and elasticity of [lambda] to changes in the vital rates (Caswell 1978, 1989). In this section we present formulas for the sensitivity and elasticity of [c.sup.*], the asymptotic invasion wave speed, to changes in the entries of the matrix A and to changes in the dispersal parameters. In the context of the models we examine, the parameter changes represent changes in the life history. But there are other, equally compelling reasons for a sensitivity analysis of this sort. For example, errors in estimating the parameters result in errors in [c.sup.*]. The most important errors will be those in the parameters to which [c.sup.*] is most sensitive. Information on the sensitivity of [c.sup.*] to changes in the parameters can thus be used to design sampling procedures that maximize the precision of the estimates of the most critical parameters. Sensi tivity analyses of invasion wave speeds may also be valuable in the evaluation of management strategies for the control of invasive species (cf. Sharov and Liebhold 1998).

The results of our analysis (presented in Appendix B) are formulas for the sensitivity of [c.sup.*] to changes in the demographic parameters [a.sub.ij],

[frac{[partial][c.sup.*]}{[partial][a.sub.ij]}] = [frac{[m.sub.ij]}{[s.sup.*][[rho].sub.1]}] [frac{[partial][[rho].sub.1]}{[partial][h.sub.ij]}] (26)

and for the elasticity of [c.sup.*] to changes in the demographic parameters,

[frac{[a.sub.ij]}{[c.sup.*]}] [frac{[partial][c.sup.*]}{[partial][a.sub.ij]}] = [frac{1}{ln [[rho].sub.1]}] [[frac{[h.sub.ij]}{[[rho].sub.1]}] [frac{[partial][[rho].sub.1]}{[partial][h.sub.ij]}]]. (27)

The quantities [m.sub.ij], [[rho].sub.1], [h.sub.ij] and their derivatives are all functions of s; in Eqs 26 and 27 they are all evaluated at s = [s.sup.*].

We have also calculated the sensitivities and elasticities of [c.sup.*] to changes in the dispersal parameters. Assume that the kernel [k.sub.ij] depends upon a parameter [[alpha].sub.ij]. The sensitivity of [c.sup.*] to [[alpha].sub.ij] is

[frac{[partial][c.sup.*]}{[partial][[alpha].sub.ij]}] = [frac{[a.sub.ij]}{[s.sup.*][[rho].sub.1]}] [frac{[partial][m.sub.ij]}{[partial][[alpha].sub.ij]}] [frac{[partial][[rho].sub.1]}{[partial][h.sub.ij]}] (28)

and the elasticity is

[frac{[[alpha].sub.ij]}{[c.sup.*]}] [frac{[partial][c.sup.*]}{[partial][[alpha].sub.ij]}] = [frac{[a.sub.ij][[alpha].sub.ij]}{[[rho].sub.1]ln[[rho].sub.1]}] [frac{[partial][[rho].sub.1]}{[partial][h.sub.ij]}] (29)

where [m.sub.ij], [[rho].sub.1], [h.sub.ij] and their derivatives are evaluated at s = [s.sup.*].

In each of these four formulas we use the well-known eigenvalue sensitivity formula to evaluate [partial][[rho].sub.1]/[partial][h.sub.ij]:

[frac{[partial][[rho].sub.1]}{[partial][h.sub.ij]}] = [frac{[v.sub.i][w.sub.j]}{(v, w)}]. (30)

Here w and v are the right and left eigenvectors of H([s.sup.*]) corresponding to [[rho].sub.1]([s.sup.*]) and [langle][cdotp][rangle] denotes the scalar product (Caswell 1978).

A step-by-step procedure to implement these calculations is outlined in Table 1.


In this section, we analyze invasion wave speeds for two plant species. Ideally, the analyses would be based on stage-specific data on both demography and dispersal. Such data are, however, seldom presented together (possibly because until now it has been difficult to combine them in a single analysis). The examples we present here are not ideal, but they will demonstrate the analyses, show the kinds of conclusions that are possible, and, we hope, inspire studies with more complete data.

Dipsacus sylvestris

Teasel (Dipsacus sylvestris) is an herbaceous mono-carpic perennial found in old fields and disturbed areas in northeastern North America, where it was introduced from Europe in the late 19th century (Werner 1975a).

Stage-structured matrix population models for teasel were developed by Werner and Caswell (1977) and Caswell and Werner (1978), using data from populations experimentally introduced into eight old fields in southern Michigan. As modified by Caswell (1989), individuals were classified as dormant seeds, small, medium, or large rosettes, and flowering plants. Only flowering plants produce dispersing propagules; these propagules may remain dormant or produce various-sized rosettes in the next year. Fig. 6 shows the life cycle. Annual population growth rates ([lambda], measured as the dominant eigenvalue of the projection matrix A) in the eight fields ranged from 0.275 to 3.792.

Teasel seeds are relatively large, and have no obvious adaptations for dispersal. Werner (1975b) measured the probability distribution of dispersal distances surrounding individual flowering plants using seed traps. She found a Laplace distribution of dispersal distances, with a mean dispersal distance [alpha] [approx] 0.257 m. More than 99% of the seeds land within a 1.5 m radius of the parent plant.

Conjecturing that the population with the highest population growth rate (Field M, [lambda] = 3.792) is representative of ideal, low-density conditions, we calculated the invasion wave speed for this population, using Eq. 24. The population projection matrix is


The dispersal matrix K(x) contains delta functions everywhere except for the last column, describing the fates of seeds. Each element of this column is equal to

[k.sub.i6](x) = [frac{1}{2[alpha]}] [e.sup.-[frac{\x\}{[alpha]}]] i = 1,..., 6 (32)

where [alpha] = 0.257 is the mean dispersal distance. The moment-generating function matrix M(s) contains ones everywhere except for the last column. Each element of the last column is

[m.sub.i6](s) = [frac{1}{1 - [([alpha]s).sup.2]}] i = 1,..., 6 (33)

where \s\ [less than] 1/[alpha].

Fig. 7 shows c(s) calculated from A and M according to Eq. 24. The minimum, [c.sup.*] = 0.5639 (m/yr) occurs for [s.sup.*] = 3.1973. Under these conditions, a teasel invasion wave would move each year roughly twice the mean seed dispersal distance.

Fig. 8 shows the wave speed [c.sup.*] as a function of [lambda] for the eight fields. In the absence of any information to the contrary, we have assumed that the same dispersal kernel applies to all fields. This might not, of course, be true.

The sensitivities and elasticities of [lambda] and of [c.sup.*] for Field M are shown in Fig. 9. To facilitate comparison with the elasticities of [lambda], which sum to 1, we have normalized the elasticities of [c.sup.*] to do the same. The patterns of both sensitivity and elasticity are the same for [lambda] and [c.sup.*] (the correlation of the logs of the two sensitivities is 0.991; the correlation of the two elasticities is 0.980). In other words, if [lambda] has a high sensitivity or elasticity to a change in a particular element of A, invasion wave speed does too. The patterns shown in Fig. 9 are typical. Both [lambda] and [c.sup.*] are most sensitive to changes in the lower left corner of A: transitions from seeds to large, reproducing individuals.

In this example, there is only one dispersal parameter, [alpha]; it appears in every entry in the last column of M. To calculate the sensitivity of [c.sup.*] to [alpha], we use the formula Eq. 28 six times, remembering that every [[alpha].sub.i6] = [alpha], and sum the results over i. The result is d[c.sup.*]/d[alpha] 2.17, which is approximately half of the largest sensitivity of [c.sup.*] to changes in the elements of A.

Calathea ovandensis

Calathea ovandensis is an herbaceous perennial plant found in the understory of neotropical successional forests. Horvitz and Schemske (1995) constructed stage-classified matrix models for populations at four sites in Mexico, over four years, yielding a total of 16 projection matrices. Individuals were classified as seeds, seedlings, juveniles, pre-reproductives, and small, medium, large, and extra-large reproductive plants. Fig. 10 shows the life cycle. In contrast to teasel, in which dispersal occurred during transitions from one stage (flowering plants) to several (dormant seeds or rosettes), in Calathea dispersal involves transitions from several stages (different-sized flowering plants) to a single stage (seeds).

Seeds of C. ovandensis are dispersed by ants. Horvitz and Schemske (1986) carried out dispersal experiments, and reported the mean distance that seeds were carried by each of several species of ants, and the proportion [p.sub.i] of seeds dispersed by each ant species. One ant (Pheidole sp.) did not move seeds at all, but seeds that it found were later dispersed by other ant species. Eliminating Pheidole and rescaling the proportions of seeds dispersed by each species yields the data in Table 2.

Horvitz and Schemske (1986) did not measure dispersal kernels, but they did report the mean dispersal distance. For the sake of this example, we assume a Laplace distribution of dispersal distances. We also assume that the fate of seeds is independent of the species of ant by which they are dispersed. Under these assumptions, the dispersal kernel matrix is

K (x - y) = [[sum].sub.i] [p.sub.i][K.sub.i](x - y) (34)

Where [K.sub.i]/(x - y) is a matrix whose first row is composed of Laplace dispersal kernels with parameter [[alpha].sub.i]. All other entries of [K.sub.i] are delta functions. Similarly, the moment-generating function matrix is

M(s) = [[sum].sub.i] [p.sub.i][M.sub.i](s). (35)

The overall average dispersal distance, weighted by the probabilities [p.sub.i], is 1.1779 m.

Population growth rates for the 16 matrices ranged from [lambda] = 0.7356 to [lambda] = 1.2477. Horvitz et al. (1997) analyze the causes of this variation. The matrix with the largest population growth rate (C3 in the notation of Horvitz and Schemske 1995) represents the best conditions for Calathea, so we examine its wave speed in detail. The invasion wave speed is [c.sup.*] =212.0 mI yr; this is more than two mean dispersal distances per year. The invasion wave speed increases with population growth rate [lambda], as shown in Fig. 11. The sensitivities and elasticities of[lambda] and [c.sup.*] to changes in the vital rates are shown in Fig. 12. Once again, there is a close correlation between the sensitivity and elasticity of [c.sup.*] and of [lambda] to changes in the vital rates.

Short- and long-distance dispersal

Each of the four ant species that disperse the seeds of Calathea carries them a different distance. The species with the longest mean dispersal distance (925 cm) disperses only 7% of the seeds. The species that disperses most (66%) of the seeds carries them only 76 cm. We can ask how important is rare long-distance dispersal by Pachycondyla apicalis compared to the much more common short-distance dispersal by P. harpax. Similar problems have been considered in other frameworks (Goldwasser et al. 1994, Lewis and Schmitz 1996, Clark 1998, Clark et al. 1998).

We explored this problem in Calathea ovandensis by setting the dispersal distances for all ant species other than Pachycondyla apicalis to zero and recalculating [c.sup.*] for all plots in which [lambda] = 1. The results are shown in Table 3. The resulting wave speeds are all within 1% of the values calculated from the entire suite of ant species. Thus in Calathea ovandensis, the long-distance dispersal events, rare though they may be, determine the invasion wave speed almost completely.

Clark et al. (1998) suggested that rare long-distance dispersal might be responsible for observed rates of spread of plant populations over geological time, which are sometimes much faster than would be predicted by mean seed-dispersal distances. Teasel may be an example. It was introduced into North America from eastern Europe, and was first collected in Canada in 1877, in Niagara Falls, Ontario. By 1900, it had reached the east coast of North America, an east-west distance of roughly 650 km (Werner 1975a). This corresponds to an average speed of [sim]27 km/yr which is faster than the calculated speed of [c.sup.*] = 0.5639 m/yr by a factor of 4.8 x [10.sup.4]. Teasel obviously did not travel this far this fast by the mechanisms included in our analysis.

Teasel seeds float, and the species has dispersed in North America mainly along waterways (Werner 1975a). Seed trap experiments, like those of Werner (1975b), could not, of course, detect or quantify this mode of dispersal. Suppose that some small fraction of teasel seeds land in water, and that when they do they are dispersed according to a Laplace distribution with some larger mean distance. Fig. 13 shows the invasion wave speed as a function of the mean distance for water-dispersed seeds, for two proportions (1.0 X [10.sup.-3] and 1.0 X [10.sup.-6]) of seeds dispersed by water. Regardless of the probability of dispersal by water, it requires a mean waterborne dispersal distance of [sim]20 km to reach a wave speed of 27 km/yr. While this seems far, it is certainly not impossible. (A seed afloat for 12 h in a river flowing at 0.5 m/s will cover over 21 km.) Another possible explanation is that multiple introductions of teasel have contributed to its distribution in North America.

An interesting aspect of Fig. 13 is the very small sensitivity of [c.sup.*] to the probability of long-distance dispersal. The wave speed resulting from a one-in-a-million event is almost identical to that resulting from a one-in-a-thousand event. Thus, undetectably rare long-distance dispersal events can completely determine the wave speed of an invading species. Fig. 14 demonstrates this by plotting the invasion wave speed as a function of the probability of long-distance dispersal.


The main result of this paper is a framework for calculating asymptotic invasion speed from stage-specific data on demography and dispersal. It combines two previously separate bodies of theory: dispersal models (which have been applied to such important issues as habitat fragmentation, epidemic spread, spatial pattern formation, and competitive coexistence) and matrix population models (which have been used to study population dynamics, life history evolution, ecotoxicology, and management strategies).

Perhaps the most useful part of this framework is the recipe for calculating the sensitivity and elasticity of wave speed to changes in both the demographic vital rates and the dispersal parameters (Table 1). Such perturbation analyses have become an important part of demographic analysis, and they should become standard practice in analysis of wave speed as well. In two relatively simple examples, we found that the sensitivity and elasticity of wave speed are highly correlated with the corresponding sensitivity and elasticity of population growth rate [lambda]. If generally true, this implies that the perturbation analysis of matrix population models can be interpreted as providing approximate information on effects of life history changes on wave speed. This possibility warrants further investigation, especially in cases including more stage-specific dispersal information than the examples we analyzed here.

Other approaches

Ours is not the first attempt to include population structure in dispersal models. Continuous-time models including age (but not stage) structure have been developed by Thieme (1977, 1979), Diekmann (1978, 1979), Van den Bosch et al. (1990), and Mollison (1991); also see Rotenberg (1972), Gurtin (1973), and Gopalsamy (1976, 1977). These models are essentially spatial extensions of Lotka's integral equation for the intrinsic rate of increase r. They have been applied to the spread of birds (Van den Bosch et al. 1992, Lensink 1997), earthworms (Marinissen and Van den Bosch 1992), and potato blight (Van den

Bosch et al. 1994). Hengeveld (1994) provides a nontechnical review. In a challenging theoretical paper, Diekmann et al. (1998) have generalized these models to include some kinds of stage structure.

The situation is analogous to that in demography, were several different mathematical frameworks (matrix models, delay differential equations (Nisbet 1997), partial differential equations (de Roos 1997), and integral equations (Diekmann et al. 1994)) coexist peacefully. Each method has its strengths and weaknesses (Caswell et al. 1997), and which method is preferred depends upon the situation--the life cycle of the organism, the form of the data, the biological question to be addressed--and the tastes and training of the analyst.

The matrix-based approach described in this paper has several strengths:

1) For many organisms, classifying individuals by stage is more accurate than classification by age, and stage-classified models are particularly easy to develop in matrix form (Caswell 1989). They are also relatively simple to analyze: the computation of invasion speed (Eq. 24) only involves finding the eigenvalue of a matrix.

2) Although we did not examine questions of life history evolution here, that theory is currently more influenced by matrix models than by other stage-structured models (Roff 1992, Stearns 1992), and we hope to visit these questions at another time.

3) The parameters in matrix population model can be estimated in a variety of ways (cohort studies, transition frequencies, mark--recapture data, inverse methods; see Caswell 2000). As a result, there is a rich and rapidly growing literature of empirical studies using these models. Dispersal kernels also lend themselves to empirical measurement (Silverman 1986, H[ddot{a}]rdle, 1990, de Wit and Floriani 1998). Without belittling the challenges involved in collecting either demographic or dispersal data, we believe that the models we describe are empirically practical.

Just as the various mathematical frameworks for incorporating demographic structure complement one another, we expect that the study of stage-structured integrodifference equation models will illuminate work on other types of spatial models, and vice versa.

The importance of population structure

It is reasonable to ask whether it is worth the trouble to include the complications of stage-structure in predicting invasion speeds. One way to address this is to compare the stage-structured result with the invasion speed computed as if the population were unstructured, using Eq. 11 and replacing b(0) with [lambda] Fig. 15 shows the result of doing just that for the teasel and Calathea populations we modeled in the Examples section. In both cases, the unstructured model overestimates the invasion speed. It is a straightforward exercise to show that, ceteris paribus, this procedure always overestimates the invasion speed because it ignores the fact that some transitions do not involve dispersal. It should also be noted that Eq. 11 implies that invasion speed increases monotonically with intrinsic growth rate. But this is not the case in a structured model; parameter changes that increase [lambda] can produce a slower invasion. (This happens in the teasel model; cf. Fig. 15a.)

However, at a conceptual level, we believe that this comparison misses the point. Including population structure is not only, or even primarily, a way to increase the accuracy of calculating wave speed. It is a way to link a result ([c.sup.*] in this case) with processes occurring within the life cycle of the individual. That is, it is a way to include the mechanisms responsible for determining the end results. As such, its value is that it provides much more insight, as well as somewhat more accuracy.

The importance of long-distance dispersal

In the examples analyzed here, when dispersal is composed of long- and short-distance dispersal, it is the long-distance component that governs the invasion speed, even when long-distance dispersal is rare. This is not a new result (e.g., Goldwasser et al. 1994, Kot et al. 1996, Lewis and Schmitz 1996, Clark 1998, Clark et al. 1998), and it makes the estimation of dispersal kernels all the more challenging because it is particularly difficult to distinguish kernels with different tails (Kot et al. 1996). However, our results suggest that most attention should be paid to those modes of dispersal with the greatest potential for long-distance movement. For example, it might be nearly impossible to estimate the proportion of teasel seeds that disperse by water, but much more feasible to conduct an experiment that would measure the distance traveled by seeds introduced into water. Figs. 13 and 14 show that knowing this distance is much more important than knowing the proportion of seeds dispersed in this way, as long as that proportion is small.

Data requirements

Our analysis clarifies the data required to analyze invasion speed. The data requirements for matrix population models are well known: a set of stages, a set of possible transitions, and estimates of transition probabilities and fertilities (e.g., Caswell 2000). The underlying idea is that the demographic properties of individuals differ among stages, and the model must incorporate those differences. Since dispersal may also differ among stages, we require a kernel describing the movement that accompanies each transition in the demographic model. Transitions that involve only negligible movement (e.g., among size classes of trees) are easy. But in principle a complete model would include an independently estimated dispersal kernel for each transition involving mobile stages.

This is not a trivial task. But it is not, we submit, hopeless. The model suggests ways to focus the effort, by concentrating on stages that disperse, by concentrating on modes of dispersal that are markedly stage-specific (something that none of the data we have encountered have included), by concentrating on stages and modes of dispersal likely to involve long-distance movement, and by concentrating on stages to whose movement [c.sup.*] is particularly sensitive.

Calculated and observed wave speeds: what if they differ?

The wave speed calculated for teasel disagrees, by a large factor, with the apparent rate of spread of the species in North America. How should such disagreements by interpreted? Like calculations of population growth rate, calculations of wave speed are projections from estimates of current conditions. If the vital rates have these values and if the dispersal kernel has this shape and if those values and that shape are maintained indefinitely, then an invasion wave will spread at the speed [c.sup.*]. Disagreement of [c.sup.*] with an observed rate of spread is not a failure of the analysis; it is evidence that the mechanisms included in the model are not those, or all of those, operating in the population, and an invitation to develop alternative hypotheses.

In the case of teasel, field biologists were already aware of the discrepancy between the short distances traveled by most seeds, and the speed with which the plant spread after its introduction to North America. This kind of discrepancy is well known in plant ecology (Clark 1998, Clark et al. 1998), as is the fact that wave speeds are highly influenced by the tails of the dispersal kernel (Kot et al. 1996). Teasel's habit of spreading along river systems was known (Werner 1975a), and suggests a possible mechanism for rapid spread. The combination of models for demography and dispersal, together with sensitivity analysis, promises to be an even more powerful tool for exploring alternative hypotheses. Our examples here only scratch the surface of such explorations.

Generalizations and directions for future research

The model we have analyzed can be extended in a number of directions. For example, we have assumed away temporal variation and spatial heterogeneity. The effects of stochastic temporal variation (environmental stochasticity) on wave speed can be explored using ergodic theory methods familiar from stochastic demography (e.g., Tuljapurkar 1990). The effects of spatial heterogeneity on wave speed have been studied in reaction-diffusion models by Shigesada et al. (1986, 1987); similar approaches can be applied to integro-difference equation models (Van Kirk and Lewis 1997; M. G. Neubert and A. LaBonte, unpublished results).

Demographic stochasticity, the result of probabilistic application of vital rates to individuals (rather than to populations), is another source of variation that we have ignored. The effects of demographic stochasticity are most evident when population size is small, and so are particularly important in the establishment and spread of an invading species. Lewis and Pacala (1999) studied an unstructured invasion model that incorporates demographic stochasticity. Their method produces deterministic integrodifference equations like Eq. 3 for the spatial moments of the system. Thus the deterministic theory is necessary for the analysis of the stochastic model. We expect the same will be true for structured models with demographic stochasticity; i.e., their analysis should require the results we have presented here.

Finally, our analyses also assume that the nonlinearities in population growth are well behaved. Allee effects, which would produce positive density dependence at the low densities in front of an advancing wave, can have major effects. In particular, models with Allee effects typically have both a critical population density and a critical population extent (Kot et al. 1996; also see Aronson and Weinberger 1975, Lewis and Kareiva 1993). If the population is initially too small, or its range is initially too restricted, the invasion will fail and the population will become extinct. How these effects might manifest themselves in a structured integrodifference equation model is an open question.

These examples, all subjects of current research, certainly do not encompass all of the complexities and vagaries of real invasions. The combination of the full richness of life cycle structure, density dependence, and dispersal in realistically complicated environments is a challenge for future research.


The manuscript was improved as a result of discussions with Mark Kot and Mark Lewis and suggestions from Roger Nisbet and two anonymous reviewers. We would like to acknowledge the National Science Foundation (DEB-9527400) who supported our work on this project. M.G. Neubert was also supported by a Texaco Research Award in Science, Technology, and Public Policy. WHOI Contribution Number 9859.


Allen, L. J. S., E. J. Allen, and S. Ponweera. 1996. A mathematical model for weed dispersal and control. Bulletin of Mathematical Biology 58:815-834.

Andow, D. A., P. M. Kareiva, S. A. Levin, and A. Okubo. 1990a. Spread of invading organisms: patterns of spread. Pages 219-242 in K. C. Kim, editor. Evolution of insect pests: the pattern of variations. Wiley, New York, New York, USA.

Andow, D. A., P. M. Kareiva, S. A. Levin, and A. Okubo. 1990b. Spread of invading organisms. Landscape Ecology 4: 177-188.

Aronson, D. G., and H. F Weinberger. 1975. Nonlinear diffusion in population genetics, combustion, and nerve propagation. Pages 5-49 in J. Goldstein, editor. Partial differential equations and related topics. Lecture Notes in Mathematics 446.

Beer, T., and M. D. Swaine. 1977. On the theory of explosively dispersed seeds. New Phytologist 78:681-694.

Bracewell, R. N. 1978. The Fourier transform and its applications. McGraw-Hill, New York, New York, USA.

Bramson, M. 1983. Convergence of solutions of the Kolmogorov equation to travelling waves. American Mathematical Society Memoirs Volume 44.

Broadbent, S. R., and D. G. Kendall. 1953. The random walk of Trichostrongylus retortaeformis. Biometrika 9:460-465.

Caswell, H. 1978. A general formula for the sensitivity of population growth rate to changes in life history parameters. Theoretical Population Biology 14:215-230.

Caswell, H. 1989. Matrix population models. Sinauer, Sunderland, Massachusetts, USA.

Caswell H. 2000. Matrix population models, 2nd Edition. Sinauer, Sunderland, Massachusetts, USA.

Catwell, H., R. M. Nisbet, A. M. de Roos, and S. Tuljapurkar. 1997. Structured-population models: many methods, a few basic concepts. Pages 3-17 in S. Tuljapurkar and H. Caswell, editors. Structured population models in marine, terrestrial and freshwater systems. Chapman & Hall, New York, New York, USA.

Caswell, H., and P. A. Werner. 1978. Transient behavior and life history analysis of teasel (Dipsacus sylvestris Huds.). Ecology 59:53-66.

Chandrasekhar, S. 1943. Stochastic problems in physics and astronomy. Reviews of Modern Physics 15: 1-91.

Clark, J. S. 1998. Why trees migrate so fast: confronting theory with dispersal biology and the paleorecord. American Naturalist 152:204-224.

Clark, J. S., C. Fastie, G. Hurtt, S. T. Jackson, C. Johnson, G. A. King, M. Lewis, J. Lynch, S. Pacala, C. Prentice, E. W. Schupp, T. Webb III, and P. Wyckoff. 1998. Reid's paradox of rapid plant migration: dispersal theory and interpretation of paleoecological records. BioScience 48:13-24.

Crosby, A. W. 1986. Ecological imperialism: the biological expansion of Europe, 900-1900. Cambridge University Press, Cambridge, UK.

de Roos, A. M. 1997. A gentle introduction to physiologically structured population models. Pages 119-204 in S. Tuljapurkar and H. Caswell, editors. Structured population models in marine, terrestrial and freshwater systems. Chapman & Hall, New York, New York, USA.

de Wit, T. D., and E. Floriani. 1998. Estimating probability densities from short samples: a parametric maximum likelihood approach. Physical Review E 58:5115-5122.

Diekmann, O. 1978. Thresholds and travelling waves for the geographical spread of infection. Journal of Mathematical Biology 6:109-130.

Diekmann, O. 1979. Run for your life. A note on the asymptotic speed of propagation of an epidemic. Journal of Differential Equations 33:58-73.

Diekmann, O., M. Gyllenberg, J. A. J. Metz, and H. R. Thieme. 1998. On the formulation and analysis of general deterministic structured population models. I. Linear theory. Journal of Mathematical Biology 36:349-388.

Diekmann, O., and J. A. J. Metz. 1994. On the reciprocal relationship between life histories and population dynamics. Pages 263-279 in S. A. Levin, editor. Frontiers in mathematical biology. Lecture Notes in Biomathematics 100.

Fisher, R. A. 1937. The wave of advance of advantageous genes. Annals of Eugenics 7:353-369.

Fitt, B. D. L., P. H. Gregory, A. D. Todd, H. A. McCartney, and O. C. MacDonald. 1987. Spore dispersal and plant disease gradients; a comparison between two empirical models. Journal of Phytopathology 118:227-242.

Goldwasser, L., Cook, and E. D. Silverman. 1994. The effects of variability on metapopulation dynamics and rates of invasion. Ecology 75:40-47.

Gopalsamy, K. 1976. On the asymptotic age distribution in dispersive populations. Mathematical Biosciences 31:191-205.

Gopalsamy, K. 1977. Age dependent population dispersion in linear habitats. Ecological Modelling 3:119-132.

Gurtin, M. E. 1973. A system of equations for age-dependent population diffusion. Journal of Theoretical Biology 40: 389-392.

H[ddot{a}]rdle, W. 1990. Applied nonparametric regression. Cambridge University Press, Cambridge, UK.

Hengeveld, R. 1989. Dynamics of Biological Invasions. Chapman & Hall, London, UK.

Hengeveld, R. 1994. Small-step invasion research. Trends in Ecology and Evolution 9:339-342.

Horn, R. A., and C. R. Johnson. 1985. Matrix analysis. Cambridge University Press, Cambridge, UK.

Horvitz, C. C., and D. W. Schemske. 1986. Seed dispersal of a neotropical myrmecochore: variation in removal rates and dispersal distance. Biotropica 18:319-323.

Horvitz, C. C., and D. W. Schemske. 1995. Spatiotemporal variation in demographic transitions of a tropical under-story herb: projection matrix analysis. Ecological Monographs 65:155-192.

Horvitz, C., D. Schemske, and H. Caswell. 1997. The "importance" of life history stages to population growth: prospective and retrospective analyses. Pages 247-272 in S. Tuljapurkar and H. Caswell, editors. Structured population models in marine, terrestrial and freshwater systems. Chapman and Hall, New York, New York, USA.

Hosono, Y. 1995. Travelling waves for a diffusive Lotka-Volterra competition model II: a geometric approach. Forma 10:235-257.

Hosono, Y. 1998. The minimal speed of travelling fronts for a diffusive Lotka-Volterra competition model. Bulletin of Mathematical Biology 60:435-448.

Howe, H. E, and L. C. Westley. 1986. Ecology of pollination and seed dispersal. Pages 185-216 in M. J. Crawley, editor. Plant ecology. Blackwell Scientific, Oxford, UK.

K[ddot{a}]llen, A., P. Arari, and J. D. Murray. 1985. A simple model for the spatial spread and control of rabies. Journal of Theoretical Biology 116:377-393.

Kametaka, Y. 1976. On the nonlinear diffusion equations of Kolmogorov-Petrovsky-Piskunov type. Osaka Journal of Mathematics 13:11-66.

Keyfitz, N. Applied mathematical demography. Springer-Verlag, New York, New York, USA.

Kolmogorov, A., N. Petrovsky, and N. S. Piscounov. 1937. [acute{E}]tude de l'[acute{e}]quations de la diffusion avec croissance de la quantit[acute{e}] de mati[grave{e}]re et son application a un probl[grave{e}]me biologique. Moscow University Bulletin of Mathematics 1:1-25.

Kot, M. 1992. Discrete-time travelling waves: ecological examples. Journal of Mathematical Biology 30:413-436.

Kot, M., M. A. Lewis, and P. van den Driessche. 1996. Dispersal data and the spread of invading organisms. Ecology 77:2027-2042.

Kot, M., and W. M. Schaffer. 1986. Discrete-time growth-dispersal models. Mathematical Biosciences 80:109-136.

Lensink, R. 1997. Range expansion of raptors in Britain and the Netherlands since the 1960s: testing an individual-based diffusion model. Journal of Animal Ecology 66:811-826.

Lewis, M. A., and P. Kareiva. 1993. Allee dynamics and the spread of invading organisms. Theoretical Population Biology 116:221-247.

Lewis, M. A., and S. Pacala. 1999. Modelling and analysis of stochastic invasion processes. Journal of Mathematical Biology, in press.

Lewis, M. A., and G. Schmitz. 1996. Biological invasion of an organism with separate mobile and stationary states: modeling and analysis. Forma 11:1-25.

Lotka, A. J. 1939. A contribution to the theory of self-renewing aggregates, with special reference to industrial replacement. Annals of Mathematical Statistics 10:1-25.

Lubina, J., and S. A. Levin. 1988. The spread rate of a reinvading organism: range expansion of the California sea Otter. American Naturalist 131:526-543.

Lui, R. 1982a. A nonlinear integral operator arising from a model in population genetics. I. Monotone initial data. SIAM Journal of Mathematical Analysis 13:913-937.

Lui, R. 1982b. A nonlinear integral operator arising from a model in population genetics. II. Initial data with compact support. SIAM Journal of Mathematical Analysis 13:93 8-953.

Lui, R. 1983. Existence and stability of travelling wave solutions of a nonlinear integral operator. Journal of Mathematical Biology 16:199-220.

Lui, R. 1985. A nonlinear integral operator arising from a model in population genetics. III. Heterozygote inferior case. SIAM Journal of Mathematical Analysis 16:1180-1206.

Lui, R. 1986. A nonlinear integral operator arising from a model in population genetics. IV. Clines. SIAM Journal of Mathematical Analysis 17:152-168.

Lui, R. 1989a. Biological growth and spread modeled by systems of recursions. I. Mathematical theory. Mathematical Biosciences 93:269-295.

Lui, R. 1989b. Biological growth and spread modeled by systems of recursions. II. Biological theory. Mathematical Biosciences 93:297-312.

Marinissen, J. C. Y., and F. van den Bosch. 1992. Colonization of new habitats by earthworms. Oecologia 91:371-376.

Markoff, A. A. 1912. Wahrscheinlichkeitsrechnung. Teubner, Leipzig, Germany.

Mollison, D. 1972. The rate of spatial propagation of simple epidemics. Proceedings of the Sixth Berkeley Symposium on Mathematical Statistics and Probability 3:579-614.

Mollison, D. 1991. Dependence of epidemic and population velocities on basic parameters. Mathematical Biosciences 107:255-287.

Mollison, D., and S. A. Levin. 1995. Spatial dynamics of parasitism. Pages 384-389 in B. T. Grenfell and A. P. Dobson, editors. Ecology of infectious diseases in natural populations. Cambridge University Press, Cambridge, UK.

Nisbet, R. M. 1997. Delay-differential equations for structured populations. Pages 89-118 in S. Tuljapurkar and H. Caswell, editors. Structured population models in marine, terrestrial and freshwater systems. Chapman & Hall, New York, New York, USA.

Neubert, M. G., M. Kot, and M. Lewis. 1995. Dispersal and pattern formation in a discrete-time predator--prey model. Theoretical Population Biology 48:7-43.

Nobel, J. V. 1974. Geographic and temporal development of plagues. Nature 250:726-728.

Okubo, A. 1980. Diffusion and ecological problems: mathematical models. Springer-Verlag, Berlin, Germany.

Okubo, A., and S. A. Levin. 1989. A theoretical framework for data analysis of wind dispersal of seeds and pollen. Ecology 70:329-338.

Pearson, K. 1905. Nature 72:294.

Roff, D. A. 1992. The evolution of life histories. Chapman & Hall, New York, New York, USA.

Rotenberg, M. 1972. Theory of population transport. Journal of Theoretical Biology 37:291-305.

Sharov, A. A., and A. M. Liehbold. 1998. Bioeconomics of managing the spread of exotic pest species with barrier zones. Ecological Applications 8:833-845.

Sharpe, F. R., and A. J. Lotka. 1911. A problem in age-distribution. Philosophical Magazine 21:435-438.

Shigesada, N., and K. Kawasaki. 1997. Biological invasions: theory and practice. Oxford University Press, Oxford, England.

Shigesada, N., K. Kawasaki, and E. Teramoto. 1986. Traveling periodic waves in heterogeneous environments. Theoretical Population Biology 30:143-160.

Shigesada, N., K. Kawasaki, and E. Teramoto. 1987. The speeds of traveling frontal waves in heterogeneous environments. Pages 88-97 in E. Teramoto and M. Yamaguti, editors. Mathematical topics in population biology, morphogenesis and neurosciences. Springer-Verlag, Berlin, Germany.

Skellam, J. G. 1951. Random disperal in theoretical populations. Biometrika 38:196-218.

Slatkin, M. 1973. Gene flow and selection in a cline. Genetics 75:733-756.

Stearns, S. C. 1992. The evolution of life histories. Oxford University Press, Oxford, UK.

Taylor, R. A. J. 1978. The relationship between density and distance of dispersing insects. Ecological Entomology 3: 63-70.

Thieme, H. R. 1977. A model for the spatial spread of an epidemic. Journal of Mathematical Biology 4:337-351.

Thieme, H. R. 1979. Density-dependent regulation of spatially distributed populations and their asymptotic speed of spread. Journal of Mathematical Biology 8: 173-187.

Tuljapurkar, S. 1990. Population dynamics in variable environments. Springer-Verlag, New York, New York, USA.

Van den Bosch, F., R. Hengeveld, and J. A. J. Metz. 1992. Analysing the velocity of animal range expansion. Journal of Biogeography 19:135-150.

Van den Bosch, F., J. A. J. Metz, and 0. Diekmann. 1990. The velocity of spatial population expansion. Journal of Mathematical Biology 28:529-565.

Van den Bosch, F., J. C. Zadoks, and J. A. J. Metz. 1994. Continental expansion of plant disease: a survey of some recent results. Pages 274-281 in J. Grasman and G. van Straten, editors. Predictability and nonlinear modelling in natural sciences and economics. Kluwer Academic, Dordrecht, The Netherlands.

Van Kirk, R. W., and M. A. Lewis. 1997. Integrodifference models for persistence in fragmented habitats. Bulletin of Mathematical Biology 59:107-138.

Weinberger, H. F. 1978. Asymptotic behavior of a model of population genetics. Pages 47-96 in J. Chadam, editor. Nonlinear partial differential equations and applications. Lecture Notes in Mathematics 684.

Weinberger, H. F. 1982. Long-time behavior of a class of biological models. SIAM Journal of Mathematical Analysis 13: 353-396.

Werner, P. A. 1975a. The biology of Canadian weeds. 12. Dipsacus sylvestris Huds. Canadian Journal of Plant Science 55:783-794.

Werner, P. A. 1975b. A seed trap for determining patterns of seed deposition in terrestrial plants. Canadian Journal of Botany 53:810-813.

Werner, P. A., and H. Caswell. 1977. Population growth rates and age versus stage-distribution models for teasel (Dipsacus sylvestris Huds.). Ecology 58:1103-1111.

Wolfenbarger, D. 0. 1946. Dispersion of small organisms, distance dispersion rates of bacteria, spores, seeds, pollen and insects: incidence rates of disease and injuries. American Midland Naturalist 35:1-152.

Wolfenbarger, D. 0. 1959. Dispersion of small organisms, incidence of viruses and pollen; dispersion of fungus, spores and insects. Lloydia 22:1-106.

Wolfenbarger, D. 0. 1975. Factors affecting dispersal distances of small organisms. Exposition, Hicksville, New York, USA.

Calculation of invasion speed and its sensitivity.

1) Generate H(s) A * M(s) using Eq. 19.

2) Determine [[rho].sub.1](s), the largest eigenvalue of H(s), for a range of s values.

3) Find [c.sup.*] and [s.sup.*] by finding the value of s that minimizes

c(s) = [frac{1}{s}] ln [[rho].sub.1] (s)

(see, for example, Fig. 4).

4) Form the matrices M([s.sup.*]) and H([s.sup.*]).

5) Calculate the sensitivity matrix for [[rho].sub.1]([s.sup.*]) using Eq. 30.

6) Compute the sensitivity and elasticity of [c.sup.*] using Eqs. 26 and 27.


We have relegated to Appendices A and B the derivation of certain results that, while mathematically important, are not essential for the understanding or use of our results.


The role of the dominant eigenvalue of H(s)

The matrix H(s) = A*M(s) is an m x m matrix, with eigenvalues and eigenvectors [[rho].sub.i](s) and [w.sub.i], i = 1,..., m. We have assumed that A is nonnegative and primitive, and by definition (Eq. 10), every element of the matrix M(s) is real and greater than or equal to one. Thus H(s) is also nonnegative and primitive, and the Perron-Frobenius theorem (Horn and Johnson 1985) guarantees that H(s) has a real positive eigenvalue [[rho].sub.1](s) that is strictly greater in magnitude than any other eigenvalue. In addition, [[rho].sub.1](s) is the only eigenvalue with a real and strictly positive eigenvector. Any other real eigenvalue will have an eigenvector with at least one negative entry.

Any linear combination of solutions of the form of Eq. 15 is also a solution of Eq. 13. We can write any such linear combination as

n(x, t) = [[[beta].sub.1][[[rho].sup.t].sub.1](s)[w.sub.1] + [[beta].sub.2][[[rho].sup.t].sub.2](s)[w.sub.2] + ... + [[beta].sub.m][[[rho].sup.t].sub.m](s)[w.sub.m]][e.sup.-sx] (A.1)

for some set of coefficients [[beta].sub.i]. This is the general solution to Eq. 13 for any initial condition of the form

n(x, 0) = [n.sub.0][e.sup.-sx]. (A.2)

Dividing each side of Eq. A.l by [[[rho].sup.t].sub.1](s) yields

[frac{n(x, t)}{[[[rho].sup.t].sub.1](s)}] = [[[beta].sub.1][w.sub.1] + [[beta].sub.2] [[lgroup][frac{[[rho].sub.2](s)}{[[rho].sub.1](s)}][rgroup].sup.t] [w.sub.2] + ... + [[beta].sub.m] [[lgroup][frac{[[rho].sub.m](s)}{[[rho].sub.1](S)][rgroup].sup.t] [w.sub.m]][e.sup.-sx]. (A.3)

Since [[rho].sub.1](s) is the largest eigenvalue, we have

[lim.sub.t[rightarrow]x] [frac{n(x, t)}{[[[rho].sup.t].sub.1](S)}] = [[beta].sub.1][w.sub.1][e.sup.-sx]. (A.4)

Thus for each s between 0 and [hat{s}], a population that begins with a spatial distribution given by (A.2) and any stage distribution will eventually achieve a stage distribution proportional to [w.sub.1]. Identifying [[rho].sub.1](s) with [], Eq. A.4 also shows that the population will spread with speed c(s) given by

c(s) = [frac{1}{s}] ln [[rho].sub.1](s). (A.5)

Since we assumed that the dominant eigenvalue of A is greater than I, and since [h.sub.ij](s) [geq] [a.sub.ij], then [[rho].sub.1](s) [greater than] 1 and c(s) [grater than] 0.

An upper bound on invasion wave speed

Following in the footsteps of Weinberger (1978, 1982), Lui (19890), and Kot et al. (1996), we can demonstrate that the minimum value of c(s) is an upper bound on the invasion speed of any population governed by Eq. 8 that begins with compact support.

First note that by Eq. 12, the nonlinear model (Eq. 8) cannot have a faster rate of spread than its linearization (Eq. 13). Next note that if

n(x, 0) [leq] [[beta].sub.1][w.sub.1][e.sup.-sx] [[beta].sub.1] [greater than] 0 (A.6)

then by Eq. 13

n(x, 1) [leq] [[beta].sub.1] [[[integral of].sup.[infty]].sub.-[infty]] [K(x - y)*A][w.sub.1][e.sup.-sy] dy (A.7)

which, after using the change of variables z = x - y, may be rewritten as follows:

n(x, 1) [leq] [[beta].sub.1]H(s)[w.sub.1][e.sup.-sx] (A.8)

But, since [w.sub.1] is the eigenvector of H(s) with eigenvalue [[rho].sub.1](s), we have

[[beta].sub.1]H(s)[w.sub.1][e.sup.-sx] = [[beta].sub.1][w.sub.1][[rho].sub.1](s)[e.sup.-sx] (A.9)

= [[beta].sub.1][w.sub.1][e.sup.-s(x-c(s))] (A.10)


n(x, 1) [leq] [[beta].sub.1][w.sub.1][e.sup.-s(x-c(s))] (A.11)

Continuing in this fashion we find that

n(x, t) [leq] [[beta].sub.1][w.sub.1][e.sup.-s(x-c(s)t)]. (A.12)

For initial conditions with compact support we can choose any s between 0 and [hat{s}]. (We can always adjust [[beta].sub.1] to satisfy Eq. A.6.) Minimizing c(s) with respect to s then gives the minimum wave speed

[c.sup.*] = [min.sub.0[less than]s[less than][hat{s}]] [[frac{1}{s}]ln [[rho].sub.1](s)]. (A.13)


Derivation of sensitivity and elasticity formulas

Our analysis begins with Eq. 24. To minimize c(s) with respect to s we differentiate Eq. 23 and set dc/ds = 0. This gives

c(s) = [frac{1}{[[rho].sub.1](s)}] [frac{d[[rho].sub.1](s)}{ds}] (B.1)

Thus the wave speed [c.sup.*] and wave shape [s.sup.*] are obtained by the simultaneous solution of the system

[c.sup.*] = [frac{1}{[s.sup.*]}] ln [[rho].sub.1]([s.sup.*]) (B.2)

[c.sup.*] = [frac{1}{[[rho].sub.1]([s.sup.*])}] [frac{d[[rho].sub.1]([s.sup.*])}{ds}] (B.3)

Both [s.sup.*] and [c.sup.*] are functions of the demographic parameters [a.sub.ij] (as are, implicitly, H(s) and [[rho].sub.1](s)); we want to know the sensitivity of the invasion speed [c.sup.*] to changes in these paramters (i.e., d[c.sup.*]/d[a.sub.ij]).

Differentiating Eq. B.2 with respect to [a.sub.ij], while holding the rest of A fixed, yields

[frac{[partial][c.sup.*]}{[partial][a.sub.ij]}] + [lgroup][frac{[c.sup.*]}{[s.sup.*]}] - [frac{1}{[s.sup.*][[rho].sub.1]([s.sup.*])}] [frac{[partial][[rho].sub.1]([s.sup.*])}{[partial][s.sup.*]}][rgroup] [frac{[partial][s.sup.*]}{[partial][a.sub.ij]}] = [frac{1}{[s.sup.*][[rho].sub.1]([s.sup.*])}] [frac{[partial][[rho].sub.1]([s.sup.*])}{[partial][a.sub.ij]}] (B.4)

By Eq. B.3 the second term vanishes and we have

[frac{[partial][c.sup.*]}{[partial][a.sub.ij]}] = [frac{1}{[s.sup.*][[rho].sub.1]([s.sup.*])}] [frac{[partial][[rho].sub.1]([s.sup.*])}{[partial][a.sub.ij]}] (B.5)

= [frac{1}{[s.sup.*][[rho].sub.1]([s.sup.*])}] [frac{[partial][[rho].sub.1]([s.sup.*])}{[partial][h.sub.ij]}] [frac{[partial][h.sub.ij]}{[partial][a.sub.ij]}] (B.6)

= [frac{[m.sub.ij]([s.sup.*])}{[s.sup.*][[rho].sub.1]([s.sup.*])}] [frac{[partial][[rho].sub.1]([s.sup.*])}{[partial][h.sub.ij]}] (B.7)

where the [h.sub.ij] are the elements of H. This expresses [partial][c.sup.*]/[partial][a.sub.ij] in terms of the sensitivity of [[rho].sub.1]([s.sup.*]) to changes in [a.sub.ij]. But since [[rho].sub.1]([s.sup.*]) is an eigenvalue of H(s), its sensitivity is given by the standard formula

[frac{[partial][[rho].sub.1]([s.sup.*])}{[partial][h.sub.ij]}] = [frac{[v.sub.i][w.sub.j]}{[langle]v, w[rangle]}], (B.8)

where w and v are the right and left elgenvectors of H([s.sup.*]) corresponding to [[rho].sub.1]([s.sup.*]) and [langle][cdotp][rangle] denotes the scalar product (Caswell 1978). Thus

[frac{[partial][c.sup.*]}{[partial][a.sub.ij]}] = [frac{[m.sub.ij]([s.sup.*])}{[s.sup.*][[rho].sub.1]([s.sup.*])}] [frac{[v.sub.i][w.sub.j]}{[langle]v, w[rangle]}]. (B.9)

The elasticity, or proportional sensitivity, of [c.sup.*] follows directly from Eqs. B.7 and B.2:

[frac{[a.sub.ij]}{[c.sup.*]}] [frac{[partial][c.sup.*]}{[partial][a.sub.ij]}] = [frac{1}{ln[[rho].sub.1]([s.sup.*])}] [[frac{[h.sub.ij]([s.sup.*])}{[[rho].sub.1]([s.sup.*])}] [frac{[partial][[rho].sub.1]([s.sup.*])}{[partial][h.sub.ij]}]]. (B.10)

In other words, the elasticity of [c.sup.*] to changes in [a.sub.ij] is proportional to the elasticity of the dominant eigenvalue of the matrix H(s). Eq. B.8 can again be used for actual calculations, i.e.,

[frac{[a.sub.ij]}{[c.sup.*]}] [frac{[partial][c.sup.*]}{[partial][a.sub.ij]}] = [frac{[a.sub.ij][m.sub.ij]([s.sup.*])}{[[rho].sub.1]([s.sup.*])ln[[rh o].sub.1]([s.sup.*])}] [frac{[v.sub.i][w.sub.j]}{[langle]v, w[rangle]}]. (B.11)
COPYRIGHT 2000 Ecological Society of America
No portion of this article can be reproduced without the express written permission from the copyright holder.
Copyright 2000 Gale, Cengage Learning. All rights reserved.

Article Details
Printer friendly Cite/link Email Feedback
Geographic Code:1USA
Date:Jun 1, 2000

Terms of use | Privacy policy | Copyright © 2021 Farlex, Inc. | Feedback | For webmasters