The effect of the spatial configuration of habitat fragmentation on invasive spread

Theoretical Population Biology 78 (2010) 298–308

Contents lists available at ScienceDirect

Theoretical Population Biology journal homepage: www.elsevier.com/locate/tpb

The effect of the spatial configuration of habitat fragmentation on invasive spread Noriko Kinezaki a , Kohkichi Kawasaki b , Nanako Shigesada b,∗ a

Faculty of Informatics, Nara Sangyo University, Sango-cho, Ikoma-gun, Nara 636-8503, Japan

b

Faculty of Culture and Information Science, Doshisha University, Kyotanabe 610-0321, Japan

article

info

Article history: Received 12 May 2010 Available online 25 September 2010 Keywords: Biological invasion Fragmented environment Corridor Species persistence Rate of spread Traveling wave Reaction–diffusion

abstract To address how the spatial configuration of habitat fragmentation influences the persistence and the rate of spread of an invasive species, we consider three simple periodically fragmented environments, a lattice-like corridor environment, an island-like environment and a striped environment. By numerically analyzing Fisher’s equation with a spatially varying diffusion coefficient and the intrinsic growth rate, we find the following. (1) When the scale of fragmentation is sufficiently large, the minimum favorable area needed for successful invasion reduces in the following order: lattice-like corridor, striped and islandlike environments. (2) When the scale of fragmentation and the fraction of favorable area are sufficiently large, the spreading speeds along contiguous favorable habitats in the lattice-like corridor and striped environments are faster than the speeds across isolated favorable habitats in the island-like environment and the striped environment. (3) When the periodicity of fragmentation is relaxed by stochastically shifting the boundaries between favorable and unfavorable habitats, the average speed increases with increases in the irregularity of fragmentation. © 2010 Elsevier Inc. All rights reserved.

1. Introduction Since the pioneering works of Fisher (1937) and Skellam (1951), diffusion–reaction equations have been widely used to describe the spatial spread of invading species including plants, animals, insects, epidemic agents and so on (Okubo, 1980; Andow et al., 1990; Shigesada and Kawasaki, 1997). Although most of those early theoretical studies assumed environments to be homogeneous, real environments are often fragmented in diverse ways by natural or artificial habitat destructions. Thus recent theoretical developments have largely been directed to more realistic situations involving environmental heterogeneity, temporal variability or interactions with other species (Chesson, 2000; Neubert et al., 2000; Hastings et al., 2005; Facon and David, 2006; Melbourne et al., 2007; Berestycki et al., 2009). In this article, we focus on how the spatial structure of habitat fragmentation influences spatio-temporal patterns of invasive species. There have been several theoretical approaches to integrating spatial heterogeneity into models (Hastings et al., 2005). They may roughly be classified into simulation approaches and theoretically oriented approaches.

∗

Corresponding author. E-mail addresses: [email protected] (N. Kinezaki), [email protected] (K. Kawasaki), [email protected] (N. Shigesada). 0040-5809/$ – see front matter © 2010 Elsevier Inc. All rights reserved. doi:10.1016/j.tpb.2010.09.002

Theoretically oriented approaches have mostly been performed in the framework of reaction–diffusion equations (RD) or integrodifference equations (IDE) (Kot et al., 1996). Shigesada et al. (1986) first introduced a periodic patchy environment that consists of favorable and unfavorable patches alternately arranged in one-dimensional space, and examined spatio-temporal patterns for a single species using the RD model. They showed that there could be a certain threshold of the fraction of favorable area above which the population expands its range in the form of a traveling periodic wave. They further presented an analytical formula for the frontal speed of the traveling periodic wave as a function of the spatially dependent diffusivity and growth rate. Subsequently, this approach was extended to semi-periodically fragmented environments (Shigesada et al., 1987), the spread of genetically modified organisms into a natural population (Cruywagen et al., 1996), directional movement in riverine habitats (Lutscher et al., 2006) and sinusoidally varying environments (Kinezaki et al., 2006). The model for one-dimensional space (Shigesada et al., 1986) was further extended to a two-dimensional case in which belt-like favorable and unfavorable habitats were alternately arranged in a striped pattern (Kinezaki et al., 2003). Recently, rigorous mathematical analyses have been done for more general periodic environments in multiple dimensions to show the existence of a traveling periodic waves and its rate of spread (Weinberger, 2002; Berestycki et al., 2005b; Nadin, 2010), and optimal habitat configurations for species persistence (Berestycki et al., 2005a; Roques and Hamel, 2007). The effect of periodic

N. Kinezaki et al. / Theoretical Population Biology 78 (2010) 298–308

a

c

b l1

l2

Striped environment

299

l1

l2

l1

l1

l1

l2

l2

Island-like environment

l2

Lattice-like corridor environment

Fig. 1. Three periodically fragmented environments. (a) Striped environment. (b) Island-like environment. (c) Lattice-like corridor environment. Areas in white and grey indicate favorable and unfavorable habitats, respectively. l1 and l2 denote the widths of the favorable and unfavorable habitats, respectively. The square area framed by thick line with side length of l (=l1 + l2 ) is designated a unit structure that represents the scale of fragmentation.

fragmentation on the spread rate of a single population in one dimension was also examined in the framework of the IDE model (Van Kirk and Lewis, 1997; Kawasaki and Shigesada, 2007; Weinberger et al., 2008). The IDE model has further been developed to incorporate an Allee effect (Dewhirst and Lutscher, 2009) or interspecies competition (Samia and Lutscher, 2010). On the other hand, individual-based simulation models have been popular because they can easily incorporate two-dimensional spatial structures of fragmentations including corridors or barriers (Hastings et al., 2005). Most of them focused on the analyses of the fragmentation threshold for successful invasion or extinction, but paid little attention to the explicit relationship between the spreading speed and the fragmented structure. In contrast, the RD and IDE models have so far dealt mostly with one-dimensional periodically or semi-periodically fragmented environments except for the two-dimensional striped environment (Kinezaki et al., 2003) or more general periodic environments that are studied in mathematically rigorous manners (Weinberger, 2002; Berestycki et al., 2005b; Nadin, 2010). Thus there has been a criticism that the RD and IDE models tend to neglect the effect of landscape structure on invasions and may not predict invasive spread adequately, even though they could be analytically tractable (With, 2002). With a view to addressing such a criticism, this work tries to extend the RD model to two-dimensional fragmented environments involving corridors or barriers. Specifically, we consider three simple but typical fragmented structures, in which favorable habitats constitute either lattice-like corridors, isolated islandlike patches or alternating stripes as interrupted by unfavorable habitats. By incorporating these fragmented environments to the Fisher equation, we examine how the spatial structure and the scale of fragmentation as well as the relative proportion of favorable habitats influence species persistence, invasion patterns in two dimensions and rates of spread. To further explore how irregularity of the fragmentation affects the rate of spread, the model is extended to randomized fragmented environments that are generated by stochastically varying the boundaries of the favorable and unfavorable habitats in the periodically fragmented environments. 2. An extended Fisher model We consider invasions of single species into three specific types of periodically fragmented environments in two dimensions: the striped environment, island-like environment and lattice-like corridor environment as illustrated in Fig. 1:

(a) A striped environment consists of favorable and unfavorable belt-like habitats with the respective widths of l1 and l2 , which are alternately arranged. (b) An island-like environment consists of square-shaped favorable patches with side l1 that are regularly separated by unfavorable belts with width l2 along both horizontal and vertical axes. (c) A lattice-like corridor environment has an opposite structure to the island-like environment, in that lattice-like favorable belts with width l1 surround square-shaped unfavorable patches with side l2 . These three fragmented environments are all composed of a periodic square unit sized l by l (l = l1 + l2 ) as indicated by a solid square frame. We call it a unit structure, which represents the scale of fragmentation. Hereafter, the area of the unit structure and the favorable area within the unit structure are denoted by U and F , respectively. Note that U and F are given as functions of l1 and l as follows: U = l2 , F =

 l1 l

l2 1 2l1 l − l21

(in the striped environment), (in the island-like environment), (in the lattice-like corridor environment).

(1)

To deal with range expansion in these fragmented environments, we employ an extended Fisher equation, in which the intrinsic growth rate and the diffusion coefficient vary depending on the habitats as follows (Kinezaki et al., 2003):

    ∂n ∂ ∂n ∂ ∂n = D(x, y) + D(x, y) ∂t ∂x ∂x ∂y ∂y + (ε(x, y) − n)n D(x, y) = D1 , D(x, y) = D2 ,

ε(x, y) = ε1 , ε(x, y) = ε2 ,

(2)

(in the favorable habitat), (in the unfavorable habitat),

where the diffusion coefficient D(x, y) is set to D1 and D2 in the favorable and unfavorable habitats, respectively. The intrinsic growth rate ε(x, y) is ε1 and ε2 in the favorable and unfavorable habitats, respectively, where ε1 is positive and greater than ε2 , while ε2 can be negative. Obviously, in the case of D1 = D2 and ε1 = ε2 , (2) is reduced to the Fisher equation for a homogeneous environment (Fisher, 1937; Shigesada and Kawasaki, 1997). As an initial condition, we assume that N0 organisms arrive at the center of the environment at t = 0, i.e., n(x, y, 0) = N0 δ(x)δ(y),

300

N. Kinezaki et al. / Theoretical Population Biology 78 (2010) 298–308

where δ(x) is the Dirac delta function. For convenience of analysis, we introduce the following scale transformation (Shigesada et al., 1986; Kinezaki et al., 2003): x′ =

D1 

x,

ε1

l′1 = d=

ε1



D2

D1

,



y′ = l1 , e=

D1 µ N0′ = N0 ,

l′2 =

ε1 D 1

y,

ε1 D1

ε2 , ε1

l2 ,

n′ =

µ n, ε1

(3) t ′ = ε1 t ,

ε1

√

where x , y′ , l′1 and l′2 are distances measured in D1 /ε1 as a unit length, e and d are the intrinsic growth rate and the diffusion coefficient in the unfavorable habitat relative to those in the favorable habitat, respectively, and n′ is the population density as a fraction of the carrying capacity, ε1 /µ. Substituting (3) into (2), we have an equation for the scaled variables: ′

  ′ ′ ∂ ∂ ∂ n′ ′ ′ ∂n ′ ′ ∂n , y ) , y ) = D ( x + D ( x ∂t′ ∂ x′ ∂ x′ ∂ y′ ∂ y′ + ε(x′ , y′ ) − n′ n′

(4)

D(x , y ) = 1, ε(x , y ) = 1, (in the favorable habitat), D(x′ , y′ ) = d, ε(x′ , y′ ) = e, (in the unfavorable habitat), Initial condition: n′ (x′ , y′ , 0) = N0′ δ(x′ )δ(y′ ), ′

′

′

′

where the intrinsic growth rate and the diffusion coefficient are both transformed to 1 in the favorable habitat, while those in the unfavorable habitat are d (>0) and e (<1), respectively. Thus our model is reduced to an equation involving four parameters: d, e, l′1 and l′2 . Hereafter, we use (4) from which the primes are omitted for notational simplicity. We will first solve (4) numerically to see how the spatiotemporal patterns of range expansions change among the three specific environments. Incidentally, for the striped environment, we have previously obtained an analytical formula for the rate of spread and found that the size of the unit structure profoundly affects the rate of spread even if the fraction of the favorable area within the unit structure is kept the same (Kinezaki et al., 2003; see also Nadin, 2010). Therefore, we will compare the rates of spread under the conditions that both the scale of fragmentation (i.e., U) and the fraction of the favorable area within the unit structure, F /U, are kept the same among the three environments. Hereafter F /U is referred to as the ‘‘favorable area fraction’’, and is denoted by α , which is obviously equivalent to the favorable area fraction in the whole landscape as well. It should be noted that for a given set of U and α , the widths of the favorable and unfavorable habitats are uniquely determined as below:

(l1 , l2 )  √ √   α U , (1 − α) U      (in the striped environment) ,  √ √ √ √   α U , (1 − α) U =  (in the island-like environment) ,    √ √ √  √      (1 − 1 − α) U , 1 − α U 

(5)

(in the lattice-like corridor environment).

In the following, we carry out numerical computations for varying values of U and α , while d and e are fixed to d = 0.5 and e = −10. This means that the death rate exceeds the birth rate in the unfavorable habitat, and that the diffusion coefficient is smaller in the unfavorable habitats than in the favorable habitat. However, we also present some results for the cases when d and e are varied from those values in Section 5.

3. Results 3.1. The pattern of range expansion We numerically solve (4) by using an alternating-direction implicit method on a two-dimensional rectangular grid under the condition that, at the boundaries between the favorable and unfavorable habitats, the population density, n, and the flux, D(x, y)∂ n/∂ x and D(x, y)∂ n/∂ y, are continuous. Fig. 2(a) shows snapshots of solutions for the striped environment, island-like environment and lattice-like corridor environment, when the area of unit structure U is set to 102 and the favorable area fraction α to 0.64. In all cases, populations expand their range, showing higher densities in the favorable habitats and lower densities in the unfavorable habitats. Fig. 2(b) illustrates the contour maps of the population density at t = 50 for the unit structure with U = 132 , while the favorable area fractions α are chosen as 0.852, 0.64 and 0.36 in each type of environment. The outermost contour in each panel corresponds to the population density of 0.001. We refer to it as the front line of invasion (note: the following results remain qualitatively valid irrespective of the exact density assigned to the front line of invasion). We also carry out numerical computations for various other sets of α and U as shown in Sections 3.2 and 3.3, and find some general features in the expansion patterns. The envelope of the front line in the striped environment roughly shows an oval-like shape elongated along the y axis (for more detail, see Kinezaki et al., 2003). On the other hand, the envelopes of front lines in the island-like environment and the lattice-like corridor environment have bulging square shapes whose corners are lying on the x and y axes. When we focus on the range expansions along an arbitrary radial line from the origin, say in direction θ , the advancing speed of the frontal line oscillates in accordance with the environmental change, but its time average converges to a constant value (Weinberger, 2002; Berestycki et al., 2005b), which is denoted by C (θ ). Thus, using polar coordinate representation, the front in direction θ at sufficiently long time t is described by r = C (θ )t. In other words, the region enclosed by the front line expands at a constant rate on average, keeping a similar shape with time. In the next section, we will compare the speeds along the x and y axes among the three environments. The rates of spread along the x axis and the y axis in the striped environment are denoted by Csx and Csy , respectively, while the rates of spread along the x axis in the island-like and lattice-like corridor environments are denoted by Ci and Cc , respectively (note that the speeds along the x axis and y axis are the same in these environments). 3.2. Rate of spread We first compare the rates of spread in the three environments, when the favorable area fraction, α , is changed from 0 to 1, while the area of the unit structure, U, is fixed at 64 and 25, as shown in the two panels of Fig. 3(a). In both panels, Csx and Csy which are obtained from an analytical formula (Kinezaki et al., 2003) are plotted in solid lines, while Ci and Cc numerically derived in this work are indicated by dotted and dashed lines, respectively. In general, each speed remains zero when the favorable area fraction is below a certain threshold value, and then increases monotonically to reach the maximum value of 2 at α = 1, which corresponds to the speed in the homogeneous environment with both the intrinsic growth rate and the diffusion coefficient being 1 (Fisher, 1937; Skellam, 1951). The shapes of the four curves, in either panel, are largely classified into two groups. Csy and Cc rapidly rise from the beginning to show convex upward curves, which turn to being convex downward shortly before α = 1. In contrast, Ci and Csx increase more gradually than the first group

N. Kinezaki et al. / Theoretical Population Biology 78 (2010) 298–308

a

t= 60

t= 75

1 0

2

2 40 1

20

0

20

0

0 0

40 0 Lattice-like corridor environment

100

100

100

80

80

80

60

60

60

40

40

40

20

20

20

0

0

20

40

60

80 100

0

0

20

40

60

80 100

0 0

100

100

100

80

80

80

60

60

60

40

40

40

20

20

20

0 0

20

40

60

80 100

0 0

20

40

60

80 100

0 0

100

100

100

80

80

80

60

60

60

40

40

40

20

20

20

0 20

40

60

80 100

Striped environment

20

40

60

80 100

20

40

60

80 100

0

0 0

20

20

40 0 Island-like environment

0 40 Striped environment

Favorable area flaction

40

20

20

b

t= 30

401

2

0

301

0

20 40 60 80 100 Island-like environment

0

20 40 60 80 100 Lattice-like corridor environment

Fig. 2. Spatial patterns of range expansion in the three periodically fragmented environments, when a few propagules are initially released at the origin. (a) Snapshots of the numerical solution of Eq. (4) derived by an alternating-direction implicit method with step sizes ∆t = 0.01 and ∆x = ∆y = 0.25. For all environments, U = 102 and α = 0.64. (b) Contour maps of population densities at t = 50: the left column, the striped environment; the middle column, the island-like environment; the right column, the lattice-like corridor environment. For all environments, U = 132 . α is chosen as 0.852, 0.64 and 0.36 for the upper, middle and lower rows. The population density on the most exterior contour is 0.001 and the difference in population densities between adjacent contours is 0.2. In both (a) and (b), d = 0.5 and e = −10.

to show slightly convex upward curves at first, which quickly turn convex downward until α = 1. This outcome may be reasoned about as follows. In the latter group, organisms have to cross unfavorable habitats in their forward direction, so their speed should be retarded. On the other hand, Csy and Cc have continuous routes of favorable habitats, so organisms can easily advance forward. Overall, the order of the speeds is Csx < Ci < Cc < Csy , when α is sufficiently large but not too close to α = 1. However, the order of the speeds is not retained when α is small. We also investigate how the rate of spread is affected when the scale of fragmentation is altered without changing the relative spatial pattern within the unit structure, in other words, when the area of the unit structure U is changed with the favorable area

fraction α kept constant. Fig. 3(b) shows the rates of spread, Csx , Csy , Ci and Cc , as functions of U, when α is fixed at 0.64 (thick lines) and 0.36 (dashed lines). In every case, the rate of spread monotonically increases from zero, when U is beyond a certain threshold. Thus it is concluded that as the scale of fragmentation is enlarged, the rate of spread is always accelerated, for any type of fragmentation. Furthermore, the order of the rates of spread is Csx < Ci < Cc < Csy for each α when the scale of fragmentation (i.e., the value of U) is sufficiently large (U > 27 for α = 0.64 and U > 88 for α = 0.36 in the present case). The same order is seen in Fig. 3(a) at sufficiently large α but not close to α = 1. Note that when U is sufficiently large, Cc is close to Csy , and Ci is close to Csx , so we can roughly estimate the upper bound of Cc by Csy and the

302

N. Kinezaki et al. / Theoretical Population Biology 78 (2010) 298–308

a

U = 64

2

Csy

Rate of spread

Csy Cc 1

Ci

Csx

Csx Striped environment

0

0.2

0.4 0.6 0.8 Favorable area fraction, α

1.0

Ci U = 25

2

Island-like

Rate of spread

environment Csy 1

Cc

Ci Csx

Cc Lattice-like corridor environment

0

b

0.2 0.4 0.6 Favorable area fraction, α

0.8

1.0

2

Csy Cc

Rate of spread

Csy Cc 1 Ci Csx Ci Csx 0

50

100 Area of unit structure, U

150

Fig. 3. (a) The rates of spread as functions of the favorable area fraction α , when the unit structure area is fixed to U = 64 and 25. Csx and Csy (both solid) are respectively the speeds along the x axis and the y axis in the striped environment, which are analytically derived. Ci (dotted) and Cc (dashed) are respectively the speeds along the x axis in the island-like and lattice-like corridor environments, which are numerically calculated. (b) The rates of spread as functions of the unit structure area U, when the favorable area fraction α is fixed to 0.64 (solid curves) and 0.36 (dashed curves). In both (a) and (b), d = 0.5 and e = −10.

lower limit of Ci by Csx , where Csx and Csy are analytically obtainable (Kinezaki et al., 2003). 3.3. The minimum favorable area for successful invasion As seen in Fig. 3(a), when the unit structure area U is given, there is a minimum value for α , referred to as α min , above which the spreading speed becomes positive. Obviously, multiplying α min

by U gives the minimum favorable area for successful invasion. Let Fsmin , Fimin and Fcmin denote the minimum favorable areas for the striped, island-like, and lattice-like corridor environments, respectively. Generally the values of the F min ’s differ among the three environments, even when the unit structure area is the same. As seen in Fig. 3(a), when U = 64, F min increases in the order of the island-like, striped, and lattice-like corridor environments (i.e., Fimin < Fsmin < Fcmin ). On the other hand, when U is

N. Kinezaki et al. / Theoretical Population Biology 78 (2010) 298–308

a

303

Fcmin

Minimum favorable area

30

Fsmin 20 min

Fi 10

0

0

40

20

60

16

36

64

2.3

2.3

2.3

3.44

3.58

3.58

80 Unit strucure area, U

b Striped environment

Island-like environment

Lattice-like corridor environment 1.65

1.78

1.82

Fig. 4. The minimum favorable area (or width) for successful invasion as a function of the unit structure area. (a) Fsmin , Fimin and Fcmin represent the minimum favorable areas in the unit structure U for successful invasion for the striped (dotted), island-like (solid), and lattice-like corridor (dashed) environments, respectively. The √ other parameters are d = 0.5, e = −10. At around the origin, F min ’s increase linearly from zero in all environments. At sufficiently large U, Fsmin is approximated by 2.3 U (thin solid), Fimin

√

converges to 12.8, and Fcmin is approximated by 2 × 1.8 U − 1.842 (thin solid). (b) The unit structures for the three kinds of environments at U = 16, 36 and 64, in which the region corresponding to the minimum favorable area F min is illustrated in white. The numerals between facing arrows indicate the minimum width of l1 for successful invasion.

decreased to 25, the order of F min is changed to the striped, islandlike and lattice-like corridor environments (i.e., Fsmin < Fimin < Fcmin ). To see whether the order of F min shows any general pattern with changes in the scale of fragmentation, we plot the F min ’s as functions of the total unit area U, as shown in Fig. 4(a), where Fimin and Fcmin are numerically calculated and Fsmin is analytically obtained (see Appendix(i)). When U is small, the F min ’s in all the three environments increase linearly with U at apparently the same slope. In fact, these slopes are mathematically shown to converge to −e/(1 − e) = 10/11 at the limit of U → 0 (see Appendix(ii)). As U further increases, Fsmin for the striped environment (dotted line) and Fcmin for the lattice-like corridor environment (dashed line) increase monotonically with gradually decelerating rates, while Fimin for the island-like environment (solid line) rapidly rises at first and quickly reaches a plateau (after passing around U = 20 in the present case). Consequently, Fimin crosses with the other two curves, first Fcmin and then Fsmin . Generally speaking, when U is sufficiently large (U > 30 in the present case), F min gets larger in the order of the island-like, striped and lattice-like corridor environments (i.e., Fimin < Fsmin < Fcmin ), and hence, invasion is more difficult in this order. To explain why the F min ’s in the three environments have such different dependencies on the scale of fragmentation, we diagrammatically illustrate the unit structures for each type of

environment at U = 16, 36 and 64 in Fig. 4(b), in which the minimum favorable areas are shown in white. The width l1 of the minimum favorable area for each fragmented environment is calculated from (5) by substituting F min /U for α . Its value seems to monotonically increase, tending to a certain limit as U increases. In fact, in the striped environment, Fsmin = 9.2, 13.8 and 18.4 for U = 16, 36 and 64, and √ the corresponding widths of the favorable stripe (i.e., l1 = Fsmin / U) are 2.3, 2.3 and 2.3, respectively. Therefore the minimum width l1 of the favorable stripe seems to already converge to 2.3 at U = 16. This implies that, in the striped environment, invasion always succeeds when l1 exceeds 2.3, no matter how large l2 (i.e., the unfavorable area) is. As a matter of fact, Kinezaki et al. (2003) analytically derived the invasion condition for the striped environment, from which the minimum width of the favorable stripe at U = ∞ is calculated to be 2.3005 (see√ Appendix(iii)). Thus by using (1), Fsmin is approximated by 2.3 U for sufficiently large U. This approximated function is plotted for 20 < U < 85 as indicated by a thin curve, which closely matches with the exact value. On the other hand, in the island-like environment, Fimin is 11.8, 12.8, 12.8 for U = 16, 36 and 64, and the corresponding widths of the favorable island are 3.44, 3.58, 3.58, respectively. Accordingly, invasion in the island-like environment may always succeed as long as the favorable area, F , exceeds 12.8. Similarly, in the

304

N. Kinezaki et al. / Theoretical Population Biology 78 (2010) 298–308

a

b

σ/l=0.1

σ/l=0.3

50

50

40

40

40

30

30

30

20

20

20

10

10

10

0

c

σ/l=0.0

50

0 0

10

20

30

40

σ/l=0.0

0 0

50

10

20

30

40

50

σ/l=0.1

0

50

50

40

40

40

30

30

30

20

20

20

10

10

10

0

10

20

30

40

50

0

10

20

30

40

50

0 0

20

30

40

50

40

50

σ/l=0.3

50

0

0

10

10

20

30

Fig. 5. (a) How to generate a randomized island-like environment. For each side of the favorable square in every unit structure, draw a parallel line shifted by a random variable chosen from a uniform distribution on the interval [−σ , σ ] as indicated by arrows. The rectangle bounded by the resultant four lines makes a new favorable habitat as indicated by the thick line. The dotted lines represent part of the favorable habitats newly generated from the adjacent unit structures. When the favorable habitat overlaps with adjacent one(s) as shown on the lower right, the overlapping area is regarded as a favorable area. When the opposing sides of the favorable square cross over, and so unfavorable areas overlap with each other (not shown), the overlapping area is taken to be unfavorable. (b) Randomized island-like environments that are generated from a unit structure of U = 25 and α = 0.64 (i.e., l = 5 and l1 = 4). The randomization index is set to σ /l = 0, 0.1 and 0.3. (c) Contour maps of population densities at t = 70 in the respective randomized island-like environments as indicated in (b), which are numerically calculated from Eq. (4) with an initial distribution localized at the origin. The population range is apparently extended more broadly as σ /l increases.

lattice-like corridor environment, we found that the minimum favorable width l1 converges to 1.84 as U increases further, so Fcmin is

√

approximated by 2 × 1.84 U − 1.842 for sufficiently large U, as indicated by a thin curve. Finally, it should be noted that the minimum fraction of the favorable habitat in the unit structure α min as given by F min /U is a monotonically decreasing function of U for all the three environments, because all F min ’s increase in a convex upward manner (see Fig. 4(a)). In particular, √ for sufficiently large U, they are approximated by αsmin ≈ 2.3/ U, αimin ≈ 12.8/U and αcmin ≈

√

3.68/ U for the striped, island-like and lattice-like corridor environments, respectively. Translating these results into the context of environmental conservation, we may say that within each type of fragmented environment, the minimum fraction of favorable area for successful invasion decreases and hence species persistence becomes easier, as the scale of fragmentation is larger. Among the three fragmented environments, the islandlike environment is most favorable for species persistence if the

area of the unit structure, U, exceeds a certain value (U > 30 in the present case). 4. Range expansion in semi-periodic environments So far, we have examined invasions into periodically fragmented environments. However, natural environments are not generally periodic, but rather constitute a complex mosaic of favorable and unfavorable patches with various shapes and sizes. Here, we will consider randomized environments that are generated from the three periodic fragmented environments (striped, island-like and lattice-like corridor) by stochastically shifting the boundaries between favorable and unfavorable habitats. We respectively call them the randomized striped, randomized islandlike and randomized lattice-like corridor environments. Fig. 5(a) shows how to generate the randomized island-like environment in detail. In a unit structure, each side of the favorable square is shifted in parallel by a random variable, which is chosen from a

a

b

305

1.5 1.0 0.5



0 Average favorable area fraction, <α>

uniform distribution on interval [−σ , σ ]. Since σ is the maximum shift length, we call σ /l the randomization index, which represents the maximum shift length relative to the side length of the unit structure, l. When the randomization index is large, there are possibilities that a part of the favorable patch overlaps with an adjacent one(s), or that the opposing boundaries of a favorable patch cross over and so the unfavorable areas overlap each other. In such cases, we unite the overlapping area into one. Similar procedures are applied to the striped and lattice-like corridor environments. Fig. 5(b) shows the randomized island-like environments that are generated from a unit structure of U = 25 and α = 0.64 (i.e., l = 5 and l1 = 4) with the randomization index set to σ /l = 0, 0.1 and 0.3. As the mathematical model for randomized environments, we again adopt the extended Fisher equation (4). Fig. 5(c) illustrates the snapshots of contour maps of the population densities at t = 70 in the respective environments as specified in Fig. 5(b). The population range is apparently extended more broadly as σ /l increases. In other words, the rate of spread seems to be faster on average when the randomization index is larger. Such a tendency is also observed for the randomized striped environment and the randomized lattice-like corridor environment. To examine the effects of the randomization on the rate of spread along the x axis or the y axis more quantitatively, we perform numerical computations of the extended Fisher model in a sufficiently large rectangular domain with width W and height L, where L > W , under the conditions that N0 organisms are initially distributed uniformly on the bottom side, and the four boundaries are subject to the Neumann condition. If A(t ) denotes the total range area at time t, in which the population density is larger than a certain level (0.001 in the present case), then C rand (t ) = A(t )/W /t represents the average rate of spread during t. As t → ∞, C rand (t ) would converge to a certain limit, which is defined as the rate of spread in the randomized environments. In fact, C rand (t ) fluctuates initially but its amplitude tends to decrease with time. In real simulations, we set the height and the bottom width of the rectangular domain to be W = 6 × l and L = 60 × l, where l is the side length of the unit structure. We take a sample average of C rand (t ) over 100 runs. This quantity tends to a constant at a time well before the front line reaches the upper boundary. Thus we use it as the average rate of spread in the randomized environment, which is denoted as ⟨C ⟩. Fig. 6(a) shows the average rates of spread as functions of randomization index σ /l in the three randomized environments, which are generated from the corresponding periodic environments with U = 25 and α = 0.64. Note that the points at σ /l = 0 represent the speeds in the periodically fragmented environments, which are equivalent to the speeds previously obtained at U = 25 in Fig. 3(b). As the randomization index increases, the average rates of spread increase monotonically in all the randomized environments. In particular, the speed rises most prominently with the randomized island-like environment, becoming nearly 3 times higher at σ /l = 0.5 than at σ = 0. In the other environments, the increase in the speed is still at least 1.4 times. It should be noted here that although the favorable area fraction α is fixed at 0.64 for the three periodic environments, the average value of α in the randomized environments as denoted by ⟨α⟩ varies with the randomization index, as shown in Fig. 6(b). This is because the overlapping areas generated upon randomization are united as described before. In particular, in the randomized island environment, the average favorable area fraction is reduced from 0.64 to 0.55 as σ /l is increased from 0 to 0.5. Thus it is obvious that accelerated speeds with randomization are not caused by changes in the favorable area fraction. Instead, as we have seen in Fig. 5(b), for σ /l = 0.3 the randomization tends to produce large lumped favorable patches, and hence the resulting environment

Average rate of spread,

N. Kinezaki et al. / Theoretical Population Biology 78 (2010) 298–308

0.1

0.2

0.3

0.4

0.5

0.7 <αc > <α s >

0.6

<α i > 0.5 0

0.1

0.2 0.3 0.4 Randamization index, σ/l

0.5

Fig. 6. The effect of randomization on the average rate of spread and the average favorable area fraction. (a) The rates of spread in the three environments as functions of the randomization index σ /l. ⟨Csx ⟩ and ⟨Csy ⟩ are the average rates of spread along the x axis and the y axis in the randomized striped environment, respectively, and ⟨Ci ⟩ and ⟨Cc ⟩ are the average rates of spread along the x axis in the randomized island-like and randomized lattice-like corridor environments, respectively. The unit structure area and the favorable area fraction in the starting periodic environment are set to U = 25 and α = 0.64, respectively. The average rates of spread increase with increases in σ /l in all the environments. (b) Average favorable area fractions as functions of the randomization index σ /l. ⟨αs ⟩, ⟨αi ⟩ and ⟨αc ⟩ are the average favorable area fractions in the striped, island-like, and the lattice-like corridor environments, respectively. ⟨αs ⟩ and ⟨αi ⟩ decrease with increases in σ /l, while ⟨αc ⟩ is almost constant. Nevertheless, the average speeds monotonically increase with increases in σ /l in all the environments.

qualitatively resembles the situation when the unit structure is enlarged in scale. Incidentally, we have shown in Section 3.2 that the speed increases as the scale of fragmentation is enlarged. A similar scale effect could apply to the increased speed in the randomized island-like environment, and plausibly also in the other randomized environments. Another feature specific to the randomized island-like environment is that stripes or corridor-like favorable areas could be newly generated in local regions, so the average speed would increase all the more. We further carried out numerical simulations for varying values of U and α , and found that randomization generally increases the average speed, though the degree of increment tends to decrease as U increases. 5. Discussion We have previously demonstrated that in the two-dimensional striped environment, the critical favorable area fraction for species persistence decreases and the spreading speed of invasion increases, as the scale of fragmentation increases (Kinezaki et al., 2003). More recently, Berestycki et al. (2005a) and Nadin (2010) gave rigorous mathematical proofs for these results for general periodic environments in any dimension. On the other hand, it has remained largely unexplored how the exact spatial configuration of fragmented environments influences such properties in the framework of the RD model. Thus we take up three simple but typical periodically fragmented environments in two dimensions, striped, island-like, and lattice-like corridor environments. We use the extended Fisher reaction–diffusion equation that includes essentially four parameters, the unit structure area U, the fraction of the favorable habitat within the unit structure α , and the

N. Kinezaki et al. / Theoretical Population Biology 78 (2010) 298–308

a

2

Rate of spread

d = 0.1

Csy

Cc 1 Ci Csx

0

0.2

0.4

0.6

0.8

1.0

Favorable area fraction,α

b

2 d = 0.5 Csy Rate of spread

diffusion coefficient d and the intrinsic growth rate e in the unfavorable habitat relative to those in the favorable habitat. Analyses are done for various values of U and α , while the other parameters are fixed as d = 0.5 and e = −10. Randomized fragmented environments are also examined to learn how irregularity in fragmentation influences the spreading speed. The major results are: (1) In each periodically fragmented environment, there is a minimum fraction of favorable area α min for population persistence, which monotonically decreases as U increase, in accord with Berestycki et al. (2005a). This means that invasion becomes easier as the scale of fragmentation is increased without changing the relative structure of fragmentation. On the other hand, when U is fixed, the value of α min differs among the three environments. In particular, when U is sufficiently large, α min is larger in the order of the island-like, striped and lattice-like corridor environments (i.e., αimin < αsmin < αcmin ), and hence, invasion is more difficult in this order. (2) When α exceeds the invasion threshold, the population expands its range in the form of a propagating wave. The spreading speed differs among the three fragmented environments, even when both U and α are kept the same in all the three environments. When U and α are sufficiently large, the speed along contiguous favorable habitats in the lattice-like corridor environment, Cc , or the striped environment, Csy , is faster than the speed across isolated favorable habitats separated by unfavorable habitats as seen in the island-like environment, Ci , or in the striped environment, Csx . This tendency becomes more remarkable as the scale of fragmentation is enlarged. Overall the order of the speed becomes Csx < Ci < Cc < Csy . One might expect here that there could be a positive correlation between this order of spreading speeds and the order of invasibility. However, this is not the case, because the island-like environment is most favorable for successful invasion among the three environments (i.e., αimin < αsmin < αcmin ) as mentioned in (1). (3) When periodicity of fragmentation is relaxed by stochastically shifting the boundaries between favorable and unfavorable habitats, the average spreading speed increases with increases in irregularity of fragmentation. Since the above results are derived in the special case of d = 0.5 and e = −10, we further investigate whether these characteristic features apply to a wider range of parameter values. By performing similar analyses for various values of d > 0 and e < 0, we confirm that the above results hold qualitatively when d < 1 and e < 0 except for in the vicinity of d = 1 and e = 0. The speed is particularly sensitive to d when this parameter is close to 0. As d decreases to near zero, the curves of Csx and Ci become more prominently convex downward, whereas Csy and Cc globally increase (see Fig. 7(a) and (b)). Consequently, for α sufficiently larger than α min , the differences between two groups, Csy and Cc versus Ci and Csx , are more enhanced as d approaches 0. The growth rate in the unfavorable habitat, e, has somewhat different effects on the rate of spread. As e increases (decreases) while d < 1, all the speeds, Csy , Cc , Ci and Csx , increase (decrease). However, their order of Csx < Ci < Cc < Csy is kept unchanged when α and U are sufficiently large. When d exceeds 1, namely, when organisms diffuse faster in the unfavorable habitat than in the favorable habitats to find better places, the order of the four speeds could change. In particular the speeds in the island-like environment, Ci , and the speed along the x axis in the striped environment, Csx , increase remarkably, because the higher dispersal in the unfavorable habitats facilitates access between adjacent isolated favorable habitats. As a result, Ci and Csx show convex upward curves for wider ranges of α , so all the four speeds become closer to each other (see Fig. 7(c)). On the other hand, as d increases, the chance of population persistence

Cc 1

Ci Csx

0

c

0.2

0.4 0.6 0.8 Favorable area fraction,α

1.0

2 d = 3.0 Csy Rate of spread

306

Cc 1

Ci Csx

0

0.2

0.4

0.6

0.8

1.0

Favorable area fraction,α Fig. 7. The effect of diffusivity in the unfavorable habitat on the rate of spread. The rates of spread as functions of the favorable area fraction α are illustrated for (a) d = 0.1, (b) d = 0.5, (c) d = 3, while U = 64 and e = −10 for all cases. As d decreases from 0.5 to 0.1, the curves of Csx and Ci become more convex downward, whereas Csy and Cc globally increase. Consequently, the difference between two groups, Csy and Cc versus Ci and Csx , is more amplified for sufficiently large α . On the other hand, as d increases from 0.5 to 3, both Csy and Cc decrease globally, while Ci and Csx increase except for α close to α min , so all the curves become closer to each other. Note that the panel (b) is the same as the upper panel in Fig. 3(a).

decreases for all the three environments, since the minimum fraction of favorable area increases in the order of (a), (b) and (c) as seen in Fig. 7. Similar results for more general periodic environments have previously been suggested by Shigesada et al. (1986) and Berestycki et al. (2005a). Combining the overall results on invasion thresholds and spreading speeds, we may draw the following conclusion. When U is sufficiently large, e < 0 and d < 1, the island-like environment is easiest to invade in the sense that the minimum fraction of

N. Kinezaki et al. / Theoretical Population Biology 78 (2010) 298–308

favorable area is smallest, whereas the speed in the island-like environment is relatively slow compared with those in the other environments. The lattice-like corridor environment is hardest to invade, but the speed is relatively high. The striped environment is intermediate as regards ease of invasion, while the speed along stripes is highest and the speed across the stripes is slowest. On the other hand, when U ≫ 1, e < 0 and d > 1, the order of invasibility remain unchanged as above, while the spreading speeds are close to each other among the three environments. Many of the previous studies have focused on how the minimum fraction of favorable area is influenced by the configuration of fragmentation in percolation-based models, neutral landscape models or individual-based models (Gardner et al., 1987; Fahrig, 1991, 2001; Higgins et al., 1996; Collingham and Huntley, 2000; With, 2002, 2004; Holt et al., 2005; Facon and David, 2006) and mathematically oriented models such as the RD equations in two dimensions (Cantrell and Cosner, 2003; Kinezaki et al., 2003; Berestycki et al., 2005a; Fagan et al., 1999). Recently, Roques and Stoica (2007) presented a hybrid model that combines an RD equation with stochastically fragmented environments. A general conclusion from their model is that more fine-grained fragmentation makes invasion more difficult for a given fraction of favorable habitat, while aggregated configurations have an opposite effect. As another interesting approach, using an RD in which only ε(x, y) is spatially heterogeneous, Roques and Hamel (2007) examined through numerical computations what spatial configuration of habitat maximizes the chances of survival for given α and U. They found that the optimal habitat configuration of favorable habitat in a unit structure changes from a ball shape, to a stripe shape, to the complement of a ball shape as α increases from 0 to 1. The focus of their analyses appears to be on which pattern is less vulnerable to disturbance under the condition that α exceeds the critical threshold for species persistence. On the other hand, we have compared the critical thresholds for species persistence among the three chosen fragmented environments. It remains an interesting question whether and how the results of these approaches are related to each other. Finally, although the present work focused on single-species reaction–diffusion models in two dimensions, real ecosystems are obviously more complicated. We would have to incorporate various additional factors such as long-distance dispersal, the Allee effect, population structure, demographic stochasticity, environmental variability, and interactions between invasive species and the local populations (Lewis and Kareiva, 1993; Shigesada et al., 1995; Neubert and Caswell, 2000; Keitt et al., 2001; Shigesada and Kawasaki, 2002; Petrovskii and Malchow, 2001; Neubert and Parker, 2004; Kawasaki et al., 2006; Fagan et al., 2009; Schreiber and Lloyd-Smith, 2009; Samia and Lutscher, 2010; Kolpasa and Nisbet, 2010). We hope that the present analysis will provide some insights for pursuit in that direction.

307

Kinezaki et al. (2003) have demonstrated that, in a striped environment, the minimum favorable area required for successful invasion is given in terms of l1 , l2 , d and e (<0) as tan

l1 2



√ −ed tanh

=

−e l2 d 2

 .

(A.1)

√

√

Substituting l1 = Fsmin / U and l2 = (U − Fsmin )/ U as given by (5) in the above equation, we have Fsmin

tan √ 2 U

   √ 1 −e min = −ed tanh (U − Fs ) . 2

dU

(A.2)

From (A.2), Fsmin is obtained as a function of U, which is illustrated by the dotted line in Fig. 4(a). (ii) The minimum favorable area Fsmin , when U is sufficiently small. When U ≪ 1, (A.2) is approximately given as Fsmin = −eU /(1 − e), since tan z ≈ z and tanh z ≈ z for small z. In particular when e = −10, the gradient, −e/(1 − e), is equal to 10/11. The slopes for the other two environments also can be shown to be −e/(1 − e), by using the formula for speed at U → 0 as presented by Nadin (2010). For Eq. (4) in the text, the spreading speed is given as lim C = 2 ⟨ε⟩A ⟨D⟩H ,



U →0

where

⟨ε⟩A =

1

∫

U

ε(x)dx =

F + e(U − F ) U

U

and

⟨D⟩H =

1 1 U

1 dx U D(x)



=

1 1 U

F 1

+

U −F d

.

Since ⟨D⟩H > 0, the minimum favorable area corresponds to the value of F at which ⟨ε⟩A = 0. Thus we have F min =

−e 10 U = U, 1−e 11

which fits well with the initial slope of F min for all the environments as shown in Fig. 4(a). (iii) The minimum width of the favorable stripe, when U is sufficiently large. When U ≫ 1, (A.2) is approximated, by using tanh z ≈ 1 for z ≫ 1, as Fsmin

√

U

= 2 arctan

√ −e d .

(A.3)

√

Fsmin / U corresponds to the minimum width of the favorable stripe, which gives 2.3005 when e = −10 and d = 1. References

Acknowledgments We thank F. Takasu and S. Takahashi for providing valuable comments and suggestions on the manuscript. We also acknowledge the reviewers for their insightful comments. This work was supported in part by a Grant-in-Aid for Scientific Research from the Japan Ministry of Education, Culture, Sports, Science, and Technology No. 18570029 to N.S. and K.K., and No. 22101004 to K.K. Appendix (i) The minimum favorable area Fsmin in the striped environment as a function of U.

Andow, D.A., Kareiva, P.M., Levin, S.A., Okubo, A., 1990. Spread of invading organisms. Landscape Ecol. 4, 177–188. Berestycki, H., Hamel, F., Roques, L., 2005a. Analysis of the periodically fragmented environment model: I—species persistence. J. Math. Biol. 51, 75–113. Berestycki, H., Hamel, F., Roques, L., 2005b. Analysis of the periodically fragmented environment model: II—biological invasions and pulsating travelling fronts. J. Math. Pures Appl. 84, 1101–1146. Berestycki, H., Diekmann, O., Nagelkerke, C.J., Zegeling, P.A., 2009. Can a species keep pace with a shifting climate?. B. Math. Biol. 71, 399–429. Cantrell, R.S., Cosner, C., 2003. Spatial Ecology via Reaction–Diffusion Equations. John Wiley and Sons, Chichester, Sussex, UK. Chesson, P., 2000. General theory of competitive coexistence in spatially-varying environments. Theor. Popul. Biol. 58, 211–237. Collingham, Y.C., Huntley, B., 2000. Impacts of habitat fragmentation and patch size upon migration rates. Ecol. Appl. 10, 131–144. Cruywagen, G.C., Kareiva, P., Lewis, M.A., Murray, J.D., 1996. Competition in a spatially heterogeneous environment: modeling the risk of spread of a genetically engineered population. Theor. Popul. Biol. 49, 1–38.

308

N. Kinezaki et al. / Theoretical Population Biology 78 (2010) 298–308

Dewhirst, C., Lutscher, F., 2009. Dispersal in heterogeneous habitats: thresholds, spatial scales, and approximate rates of spread. Ecology 90, 1338–1345. Facon, B., David, P., 2006. Metapopulation dynamics and biological invasions: a spatially explicit model applied to a freshwater snail. Am. Nat. 168, 769–783. Fagan, W.F., Cantrell, R.S., Cosner, C., 1999. How habitat edges change species interactions. Am. Nat. 153, 165–182. Fagan, W.F., Cantrell, R.S., Cosner, C., Ramakrishnan, S., 2009. Interspecific variation in critical patch size and gap-crossing ability as determinants of geographic range size distributions. Am. Nat. 173, 363–375. Fahrig, L., 1991. Simulation methods for developing general landscape-level hypotheses of single-species dynamics. In: Turner, M.G., Gardner, R.H. (Eds.), Quantitative Methods in Landscape Ecology. Springer-Verlag, USA, pp. 417–442. Fahrig, L., 2001. How much habitat is enough? Biol. Conserv. 100, 65–74. Fisher, R.A., 1937. The wave of advance of advantageous genes. Ann. Eugenic. 7, 355–369. Gardner, R.H., Milne, B.T., Turnei, M.G., O’Neill, R.V., 1987. Neutral models for the analysis of broad-scale landscape pattern. Landscape Ecol. 1, 19–28. Hastings, A., Cuddington, K., Davies, K.F., Dugaw, C.J., Elmendorf, S., Freestone, A., Harrison, S., Holland, M., Lambrinos, J., Malvadkar, U., Melbourne, B.A., Moore, K., Taylor, C., Thomson, D., 2005. The spatial spread of invasions: new developments in theory and evidence. Ecol. Lett. 8, 91–101. Higgins, S.I., Richardson, D.M., Cowling, R.M., 1996. Modeling invasive plant spread: the role of plant–environment interactions and model structure. Ecology 77, 2043–2054. Holt, R.D., Keitt, T.H., Lewis, M.A., Maurer, B.A., Taper, M.L., 2005. Species’ borders: single species approaches. Oikos 108, 18–27. Kawasaki, K., Shigesada, N., 2007. An integrodifference model for biological invasions in a periodically fragmented environment. Japan. J. Indust. Appl. Math. 24, 3–15. Kawasaki, K., Takasu, F., Caswell, H., Shigesada, N., 2006. How does stochasticity in colonization accelerate speeds of invasion in a cellular automaton model? Ecol. Res. 21, 334–345. Keitt, T.H., Lewis, M.A., Holt, R.D., 2001. Allee effects, invasion pinning, and species’ borders. Am. Nat. 157, 203–216. Kinezaki, N., Kawasaki, K., Takasu, F., Shigesada, N., 2003. Modeling biological invasions into periodically fragmented environments. Theor. Popul. Biol. 64, 291–302. Kinezaki, N., Kawasaki, K., Shigesada, N., 2006. Spatial dynamics of invasion in sinusoidally varying environments. Popul. Ecol. 48, 263–270. Kolpasa, A., Nisbet, R.M., 2010. Effects of demographic stochasticity on population persistence in advective media. B. Math. Biol. doi:10.1007/s11538-009-9489-4. Kot, M., Lewis, M.A., van den Driessche, P., 1996. Dispersal data and the spread of invading organisms. Ecology 77, 2027–2042. Lewis, M.A., Kareiva, P., 1993. Allee dynamics and the spread of invading organisms. Theor. Popul. Biol. 43, 141–247. Lutscher, F., Lewis, M.A., McCauley, E., 2006. Effects of heterogeneity on spread and persistence in rivers. B. Math. Biol. 68, 2129–2160. Melbourne, B.A., Cornell, H.V., Davies, K.F., Dugaw, C.J., Elmendorf, S., Freestone, A.L., Hall, R.J., Harrison, S., Hastings, A., Holland, M., Holyoak, M., Lambrinos, J., Moore, K., Yokomizo, H., 2007. Invasion in a heterogeneous world: resistance, coexistence or hostile takeover?. Ecol. Lett. 10, 77–94.

Nadin, G., 2010. The effect of the Schwarz rearrangement on the periodic principal eigenvalue of a nonsymmetric operator. SIAM J. Math. Anal. 41, 2388–2406. Neubert, M.G., Caswell, H., 2000. Demography and dispersal: calculation and sensitivity analysis of invasion speed for structured population. Ecology 81, 1613–1628. Neubert, M.G., Kot, M., Lewis, M.A., 2000. Invasion speeds in fluctuating environments. Proc. R. Soc. Lond. B 267, 1603–1610. Neubert, M.G., Parker, I.M., 2004. Projecting rates of spread for invasive species. Risk Anal. 24, 817–831. Okubo, A., 1980. Diffusion and Ecological Problems: Mathematical Models. Springer-Verlag, Berlin. Petrovskii, S.V., Malchow, H., 2001. Wave of chaos: new mechanism of pattern formation in spatio-temporal population dynamics. Theor. Popul. Biol. 59, 157–174. Roques, L., Hamel, F., 2007. Mathematical analysis of the optimal habitat configurations for species persistence. Math. Biosci. 210, 34–59. Roques, L., Stoica, R.S., 2007. Species persistence decreases with habitat fragmentation: an analysis in periodic stochastic environments. J. Math. Biol. 55, 189–205. Samia, Y., Lutscher, F., 2010. Coexistence and Spread of Competitors in Heterogeneous Landscapes. B. Math. Biol. doi:10.1007/s11538-010-9529-0. Schreiber, S.J., Lloyd-Smith, J.O., 2009. Invasion dynamics in spatially heterogeneous environments. Am. Nat. 174, 490–505. Shigesada, N., Kawasaki, K., 1997. Biological Invasions: Theory and Practice. Oxford University Press, Oxford, UK. Shigesada, N., Kawasaki, K., 2002. Invasion and the range expansion of species: effects of long-distance dispersal. In: Bullock, J., Kenward, R., Hails, R. (Eds.), Dispersal. Blackwell Science, Oxford, pp. 350–373. Shigesada, N., Kawasaki, K., Teramoto, E., 1986. Traveling periodic waves in heterogeneous environments. Theor. Popul. Biol. 30, 143–160. Shigesada, N., Kawasaki, K., Teramoto, E., 1987. The speeds of traveling frontal waves in Heterogeneous environments. In: Teramoto, E., Yamaguti, M. (Eds.), Mathematical Topics in Population Biology, Morphogenesis and Neurosciences. In: Lecture Notes in Biomathematics, vol. 71. Springer, Berlin, pp. 87–97. Shigesada, N., Kawasaki, K., Takeda, Y., 1995. Modeling stratified diffusion in biological invasions. Am. Nat. 146, 229–251. Skellam, J.G., 1951. Random dispersal in theoretical populations. Biometrika 38, 196–218. Van Kirk, R.W., Lewis, M.A., 1997. Integrodifference models for persistence in fragmented habitats. B. Math. Biol. 59, 107–137. Weinberger, H.F., 2002. On spreading speeds and traveling waves for growth and migration models in a periodic habitat. J. Math. Biol. 45, 511–548. Weinberger, H.F., Kawasaki, K., Shigesada, N., 2008. Spreading speeds of spatially periodic integro-difference models for populations with nonmonotone recruitment functions. J. Math. Biol. 57, 387–411. With, K.A., 2002. The landscape ecology of invasive spread. Conserv. Biol. 16, 1192–1203. With, K.A., 2004. Assessing the risk of invasive spread in fragmented landscape. Risk Anal. 24, 803–815.