Climate Change and Integrodifference Equations in a Stochastic Environment

Abstract

Climate change impacts population distributions, forcing some species to migrate poleward if they are to survive and keep up with the suitable habitat that is shifting with the temperature isoclines. Previous studies have analysed whether populations have the capacity to keep up with shifting temperature isoclines, and have mathematically determined the combination of growth and dispersal that is needed to achieve this. However, the rate of isocline movement can be highly variable, with much uncertainty associated with yearly shifts. The same is true for population growth rates. Growth rates can be variable and uncertain, even within suitable habitats for growth. In this paper, we reanalyse the question of population persistence in the context of the uncertainty and variability in isocline shifts and rates of growth. Specifically, we employ a stochastic integrodifference equation model on a patch of suitable habitat that shifts poleward at a random rate. We derive a metric describing the asymptotic growth rate of the linearised operator of the stochastic model. This metric yields a threshold criterion for population persistence. We demonstrate that the variability in the yearly shift and in the growth rate has a significant negative effect on the persistence in the sense that it decreases the threshold criterion for population persistence. Mathematically, we show how the persistence metric can be connected to the principal eigenvalue problem for a related integral operator, at least for the case where isocline shifting speed is deterministic. Analysis of dynamics for the case where the dispersal kernel is Gaussian leads to the existence of a critical shifting speed, above which the population will go extinct, and below which the population will persist. This leads to clear bounds on rate of environmental change if the population is to persist. Finally, we illustrate our different results for butterfly population using numerical simulations and demonstrate how increased variances in isocline shifts and growth rates translate into decreased likelihoods of persistence.

Appendices

Appendix 1: Proof of Theorem 1 and 2

The proof of Theorems 1 and 2 mostly follows from Hardin et al. (1988)[Theorem 4.2 and Theorem 5.3]. Indeed Hardin et al. (1988) study the general stochastic model

$$\begin{aligned} X_{t+1}=F_{\alpha _t}(X_t), \end{aligned}$$

(69)

with $(X_t)_t$ a random process that takes values in the set on nonnegative continuous function in $\varOmega$, $(\alpha _t)_t$, independent identically distributed random variable taking values in the set $\mathcal {S}$. And they have the following assumption for $F_\alpha$

• (H1) For each $\alpha \in \mathcal {S}$, $F_\alpha$ is a continuous map of $C_+(\varOmega )$, into itself such that $F_\alpha (u)=0\in C_+(\varOmega )$ if and only if $u=0\in C_+(\varOmega )$.

• (H2) If $u,\, v\in C_+(\varOmega )$ and $u\ge v$ then $F_\alpha (u)\ge F_\alpha (v)$.

• (H3) There exists some $b>0$ such that for $u\in C_+(\varOmega )$

  • (a) $||F_\alpha (u)||_\infty \le b$ for all $\alpha \in \mathcal {S}$ and $||u||_\infty \le b$,

  • (b) there exists some time t (depending on $u_0$) such that

    $$\begin{aligned} ||F_{\alpha _t}\circ \cdots \circ F_{\alpha _0}u_0||_\infty <b, \end{aligned}$$

    (70)

    for all $\alpha _0,\dots ,\alpha _t\in \mathcal {S}$,

  • (c) there exists some $d>0$ such that

    $$\begin{aligned} F_\alpha (b)\ge d \end{aligned}$$

    (71)

    for all $\alpha \in \mathcal {S}$.

• (H4) Let $B_b:=\{u\in C_+(\varOmega ):\, ||u||_\infty \le b\}$, then there is some compact set $D\subset C_+(\varOmega )$ such that $F_\alpha (B_b)\subset D$ for all $\alpha \in \mathcal {S}$.

• (H5) There exists some $h>0$ such that

  $$\begin{aligned} ||F_\alpha (u)||_\infty \le h||u||_\infty , \end{aligned}$$

  (72)

  for all $\alpha \in \mathcal {S}$ and $u\in C_+(\varOmega )$.

• (H6) There exists some $\xi >0$ such that

  $$\begin{aligned} F_\alpha (B_b)\subset K_\xi , \end{aligned}$$

  (73)

  where $K_\xi :=\{u\in C_+(\varOmega ):\,u\ge \xi ||u||_\infty \}$.

• (H7) For each $a>0$, there exists a continuous function $\tau :(0,1]\rightarrow (0,1]$ such that $\tau (s)>s$ for all $s\in (0,1)$ and such that

  $$\begin{aligned} \tau (s)F_\alpha (u)\le F_\alpha (su) \end{aligned}$$

  (74)

  for all $\alpha \in \mathcal {S}$ and $u\in C_+(\varOmega )$ such that $a\le u\le b$.

• (H8) $F_\alpha$ is Fréchet differentiable (with respect to $C_+(\varOmega )$) at $0\in C_+(\varOmega )$. We denote by $\mathcal {L}_\alpha$ the operator $F'_\alpha (0)$

• (H9) There exists a function $\mathcal {N}:{\mathbb R}_+\rightarrow [0,1]$ such that

  $$\begin{aligned} \underset{u\rightarrow 0,\,u>0}{\lim } \mathcal {N}(u)=1 \text { and } \mathcal {N}(||u||_\infty )\mathcal {L}_\alpha u\le F_\alpha (u)\le \mathcal {L}_\alpha u, \end{aligned}$$

  (75)

  for all $u\in C_+(\varOmega )$.

We want to prove that the previous hypotheses are satisfied in our framework and then apply Theorem 4.2, Theorem 5.3 by Hardin et al. (1988) and Theorem 2 by Jacobsen et al. (2015). One can check, as it is done by Jacobsen et al. (2015, Section 5.3), that under Hypotheses 1–5, the previous hypotheses (H1)–(H9) are satisfied and one can then apply Theorem 4.2, Theorem 5.3 from Hardin et al. (1988) and Theorem 2 from Jacobsen et al. (2015). For completeness, we will write the main steps referring most of the time to Jacobsen et al. (2015, Section 5.3). First, we denote by,

$$\begin{aligned} \overline{K}:=\sup \{K(x), x\in {\mathbb R}\}<+\infty , \end{aligned}$$

(76)

and

$$\begin{aligned} \underline{K}:=\inf \left\{ K(x), \, x\in ( \inf \varOmega -\sup \varOmega +c,\sup \varOmega -\inf \varOmega +c)\right\} >0. \end{aligned}$$

(77)

These two constants will be used several times in the proof of (H1)–(H9) below.

• (H1) The continuity of $F_\alpha$ follows from the continuity of $f_r$, and the boundedness of K and $\int _\varOmega g_0(y)dy$. K positive in ${\mathbb R}$, $f_r$ and $g_0$ nonnegative yield the second statement.

• (H2) This follows from the monotonicity of $f_r$ for all $\alpha \in \mathcal {S}$.

• (H3) The constant $b>0$ will be defined later in the proof,

  • (a) for all $\alpha \in \mathcal {S}$, $u\in C_+(\varOmega )$,

    $$\begin{aligned} ||F_\alpha (u)||_\infty&=\underset{x\in \varOmega }{\max }\int _ \varOmega K(x-y+c)g_0(y-\sigma )f_r(u(y))\mathrm{d}y\end{aligned}$$

    (78)
    $$\begin{aligned}&< m \cdot \overline{K}\cdot \int _{\varOmega _0} g_0(y)\mathrm{d}y=b. \end{aligned}$$

    (79)

    This proves the statement.

  • (b) Using the same argument as before, for all $u_0\in C_+(\varOmega )$, $t\in {\mathbb N}^*$, $\alpha _t,\dots ,\alpha _0\in \mathcal {S}$,

    $$\begin{aligned} ||F_{\alpha _t}\circ \cdots \circ F_{\alpha _0}(u_0)||_\infty =||F_{\alpha _t}(u)||_\infty <b, \end{aligned}$$

    (80)

    where $u\in C_+(\varOmega )$.

  • (c) For all $x\in \varOmega$,

    $$\begin{aligned} F_\alpha (b)(x)\ge \underset{\alpha \in \mathcal {S}}{\inf }\,f_r(b)\int _\varOmega K(x-y+c)g_0(y-\sigma )dy\ge \underset{\alpha \in \mathcal {S}}{\inf }\,f_r(b)\cdot d_1\nonumber \\ \end{aligned}$$

    (81)

    with

    $$\begin{aligned} d_1:= \underline{K}\cdot \int _{\varOmega _0} g_0(y)\mathrm{d}y>0 \end{aligned}$$

    (82)

    One concludes using the positivity of K and Hypothesis 4(iii)b that there exists $d>0$ such that for all $x\in \varOmega$,

    $$\begin{aligned} F_\alpha (b)(x)\ge d. \end{aligned}$$

    (83)

• (H4) This statement follows from the continuity of K, uniform boundedness of $F_\alpha$ and $f_r$ and Hypothesis 5, details can be found in Jacobsen et al. (2015)[Section 5.3].

• (H5) From assumption 4(ii)c, we have that for all $u>0$, $f'_r(0)u\ge f_r(u)$. Thus, using this inequality , with assumptions 1(ii), 2(ii) and 4 (iii)a, we get for all $u\in C_+(\varOmega )$, for all $\alpha \in \mathcal {S}$,

  $$\begin{aligned} ||F_\alpha (u)||_\infty \le h||u||_\infty , \end{aligned}$$

  (84)

  with $h:=\overline{r}\cdot \overline{K}\cdot \int _\varOmega g_0(y)\mathrm{d}y$.

• (H6) First, notice that for all $\alpha \in \mathcal {S}$, $u\in B_b$

  $$\begin{aligned} ||F_\alpha (u)||_\infty\le & {} \overline{r}\cdot \overline{K}\cdot \int _\varOmega g_0(y-\sigma )u(y)\mathrm{d}y \, \implies \, \int _\varOmega g_0(y-\sigma )u(y)\nonumber \\\ge & {} \frac{||F_\alpha (u)||_\infty }{\overline{r}\cdot \overline{ K}}, \end{aligned}$$

  (85)

  and using Hypothesis 4(ii)c,

  $$\begin{aligned} f_\alpha (u)\ge \frac{f_\alpha (b)}{b}u. \end{aligned}$$

  (86)

  Then for all $x\in \varOmega$

  $$\begin{aligned} F_\alpha (u)(x)&\ge \underline{K}\cdot \frac{f_\alpha (b)}{b} \int _\varOmega g_0(y-\sigma )u(y)dy \end{aligned}$$

  (87)
  $$\begin{aligned}&\ge \frac{\underline{K}}{\overline{K}}\cdot \frac{f_\alpha (b)}{b} \frac{||F_\alpha (u)||_\infty }{\overline{r}} \end{aligned}$$

  (88)

  and the statement is proved.

• (H7) This proof is also derived from Jacobsen et al. (2015). We want to find a continuous function $\tau :\, (0,1]\rightarrow (0,1]$ such that for all $s\in (0,1)$,

  $$\begin{aligned}&\tau (s)F_\alpha (u)\le F_\alpha (su)\\ \Leftrightarrow&\int _\varOmega K(x-y+c)g_0(y-\sigma )\tau (s)f_r(u(y)) dy\le \int _\varOmega K(x-y+c)\\&g_0(y-\sigma )f_r(su(y))\mathrm{d}y \end{aligned}$$

  Thus, it is sufficient to have $\tau (s)f_r(u(y))\le f_r(su(y))$ for all $y\in \varOmega$. From Hypothesis 4(ii)c, as $s\in (0,1)$, we have that for all $\alpha \in \mathcal {S}$

  $$\begin{aligned} s<\frac{f_r(su)}{f_r(u)}, \end{aligned}$$

  (89)

  and letting $\tau (s)=\min \left\{ \frac{f_r(su)}{f_r(u)},\, \alpha \in \mathcal {S},\, a\le u\le b\right\}$ we have, for all $s\in (0,1)$, $\alpha \in \mathcal {S}$ and $u\in [a,b]$,

  $$\begin{aligned} \tau (s)>s \text { and } \frac{f_r(su)}{f_r(u)}\ge \tau (s) \end{aligned}$$

  (90)

  and the statement is proved.

• (H8) We want to prove that

  $$\begin{aligned} \underset{h\rightarrow 0}{\lim }\frac{||F_\alpha (0+h)-F_\alpha (0)-\mathcal {L}_\alpha h||_\infty }{||h||_\infty }=0, \end{aligned}$$

  (91)

  where $\mathcal {L}_\alpha :=F'_\alpha (0)$ is a linear operator. Using the differentiability of $f_r$ at 0, one proves that the limit exists and $\mathcal {L}_\alpha h=r\int _\varOmega K(x-y+c)g_0(y-\sigma )h(y)\mathrm{d}y$.

• (H9) The second part of the inequality follows from assumption 4(ii)c. For each $\alpha \in \mathcal {S}$ let $\mathcal {N}_\alpha :[0,+\infty )\rightarrow [0,1]$ be such that

  $$\begin{aligned} \mathcal {N}(u)= {\left\{ \begin{array}{ll} \frac{f_\alpha (u)}{r\cdot u} &{}\text {if } u>0,\\ 1 &{}\text {if } u=0. \end{array}\right. } \end{aligned}$$

  (92)

  The function $\mathcal {N}_\alpha$ is continuous and defining $\mathcal {N}(u)=\min \left\{ \mathcal {N}_\alpha (u),\, \alpha \in \mathcal {S}\right\}$, it gives the wanted statement (using the fact that $\mathcal {N}_\alpha (u(x))\ge \mathcal {N}_\alpha (||u||_\infty )$).

We can thus apply Theorem 4.2 from Hardin et al. (1988) to prove Theorem 1 and use Theorem 5.3 by Hardin et al. (1988) and Theorem 2 by Jacobsen et al. (2015) to prove Theorem 2. For sake of completeness, we include the three Theorems from Hardin et al. (1988) and Jacobsen et al. (2015) below.

Theorem 3

[Hardin et al. (1988, Theorem 4.2)] Suppose that Hypotheses (H1)–(H7) above and Hypothesis 3(i) in Sect. 2 are satisfied and that $u_0\ne 0\in C_+(\varOmega )$ with probability one. Then $u_t$ solution of (14) converges in distribution to a stationary distribution $\mu ^*$, independent of $u_0$, such that either $\mu ^*(\{0\})= 0$ or $\mu ^*(\{0\}) = 1$.

The above Theorem from Hardin et al. (1988) states that the distribution $\mu ^*$ is stationary in the sense that it satisfies

$$\begin{aligned} P \mu ^*=\mu ^* \end{aligned}$$

(93)

where P is the Markov process given by

$$\begin{aligned} P \mu (B)=\int \mu (F^{-1}_\alpha (B))d\mathcal {P}_\alpha , \end{aligned}$$

(94)

with $B\in \mathcal {B}(C_+(\varOmega )$ the family of borel sets in $C_+(\varOmega )$. This is equivalent to writing that the random variable $u^*$ with distribution $\mu ^*$ satisfies

$$\begin{aligned} u^*=F_\alpha ^*[u^*] \end{aligned}$$

(95)

for some $\alpha ^*$ taking its values in $\mathcal {S}$ with distribution $\mathcal {P}_\alpha$.

Theorem 4

[Hardin et al. (1988, Theorem 5.3)] Suppose that Hypotheses (H1)–(H9) above and Hypothesis 3(i) in Sect. 2 are satisfied and that $u_0\ne 0\in C_+(\varOmega )$ with probability one. Let $\mu ^*$ be as in Theorem 4.2 and define $R=\lim _{t\rightarrow +\infty } ||\mathcal {L}_{\alpha _t}\circ \dots \circ \mathcal {L}_{\alpha _0}||$, with $\mathcal {L}$ the linearised operator around zero ((19))

• (a) If $R< 1$ then $\mu ^*(\{0\}) = 1$ and $u_t\rightarrow 0$ with probability one.

• (b) If $R > 1$ then $\mu ^*(\{0\}) =0$.

Theorem 5

[Jacobsen et al. (2015, Theorem 2)] Let R be defined as in the previous theorem by $R=\lim _{t\rightarrow +\infty }||\mathcal {L}_{\alpha _t}\circ \dots \circ \mathcal {L}_{\alpha _0}||$, then $\varLambda =R$ ($\varLambda$ defined in (18)).

Appendix 2: Persistence Condition and Invasion Speed

In this section, we derive an heuristic criteria for persistence inspired from Neubert et al. (2000) linking the critical patch size and the asymptotic invasion speed to get necessary conditions for persistence. We consider problem (9), where the suitability functions $g_0$ are indicator functions, that is

$$\begin{aligned} g_0(y)=\mathbbm {1}_{\varOmega _0}= {\left\{ \begin{array}{ll} 1 &{}\text {if }y\in \varOmega _0,\\ 0 &{}\text {otherwise.} \end{array}\right. } \end{aligned}$$

In addition to the assumptions made in Sect. 2, we will also assume that K is a thin tail dispersal kernel in the sense that it has exponentially bounded tails, i.e., there exists $s>0$,

$$\begin{aligned} \int _{\mathbb R}e^{s|x|}K(x)\mathrm{d}x<+\infty . \end{aligned}$$

(96)

This last assumption guarantees that the moment generating function (29) exists on some open interval of the form $(0,s^+)$. This assumption was not necessary to derive the persistence condition in Sect. 3 but it will be used to study the speed of the stochastic wave.

1.1 Appendix 2.1: Critical Domain Size in the Constant Environment

Let us first consider the problem in the nonshifted frame and thus assume that $s_t\equiv 0$, where $(s_t)_t$, defined in Sect. 2, is the centre of the suitable habitat at generation t. Thus, we consider the following integrodifference equation

$$\begin{aligned} u^0_{t+1}(x)=\int _{\varOmega _0} K(\xi -\eta )f_{r_t}(u^0_t(\eta ))\mathrm{d}\eta , \end{aligned}$$

(97)

denoting by $u^0$ the solution in this nonshifted framework. Using the theorems in Sect. 3, we get that as time goes to infinity, $u^0$ persists if

$$\begin{aligned} E[\ln (r_0)] >-\ln (\lambda _0), \end{aligned}$$

(98)

where $E[\cdot ]$ is the expectation of a random variable and $\lambda _0$ is the principal eigenvalue of the linear operator $\mathcal {K}_0$:

$$\begin{aligned} \mathcal {K}_0[u](x)=\int _{\varOmega _0}K(x-y)u(y)\mathrm{d}y. \end{aligned}$$

(99)

Now assume that K, $\varOmega _0$ and $(r_t)_t$ are such that (98) is satisfied and study the problem in the nonmoving, homogeneous framework.

1.2 Appendix 2.2: Invasion Speed in a Stochastic Homogeneous Environment

Now we are interested in deriving the asymptotic invasion speed of the population in an homogeneous environment to compare it with the forced shifting speed c. We thus consider the homogeneous problem on ${\mathbb R}$ in the nonmoving frame, i.e., let $(n_t)_t$ be the solution of the equation

$$\begin{aligned} n_{t+1}(\xi )=\int _{\mathbb R}K(\xi -\eta )f_{r_t}(n_t(\eta ))\mathrm{d}\eta . \end{aligned}$$

(100)

As we are considering the initial problem in the nonmoving frame, K does not depend on c and $g_0\equiv 1$ in ${\mathbb R}$ and thus the stochasticity comes only from the growth term. From the analysis of Neubert et al. (2000), there are two different approaches to estimate the invasion speed of the stochastic process $(n_t)_t$. One can either consider the invasion speed of the expected wave or the asymptotic speed of the stochastic wave. We will only consider the latter approach and assume that the speed is governed by the linearisation at 0. Denote by $(\tilde{n}_t)_t$, the solution of the linearised operator at 0, i.e., for all $t\in {\mathbb N}$,

$$\begin{aligned} \tilde{n}_{t+1}(\xi )=\int _{\mathbb R}K(\xi -\eta )r_t\tilde{n}_t(\eta )\mathrm{d}\eta . \end{aligned}$$

(101)

We define the random variable $\varXi _t$ as the most rightward position such that $\tilde{n}_t$ is greater that some threshold, i.e.,

$$\begin{aligned} \varXi _t=\sup \left\{ \xi \in {\mathbb R},\, \tilde{n}_t>\overline{n}\right\} , \end{aligned}$$

(102)

where $\overline{n}\in (0,1)$ is a fixed critical threshold. Assume that $\forall \xi \in {\mathbb R}$, $n_0(\xi )=\alpha e^{-s\xi }$, for some $s>0$, i.e., the initial condition has a wave shape, then for all $t\in {\mathbb N}$, $\xi \in {\mathbb R}$,

$$\begin{aligned} \tilde{n}_{t+1}(\xi )=\alpha \prod _{i=0}^t(r_iM(s))e^{-s\xi }, \end{aligned}$$

(103)

where M is the moment generating function of K ((29)). This function exists in some interval $(0,s^+)$ because of assumption (96). Moreover, $\overline{n}=n_0(\varXi _0)=\tilde{n}_{t+1}(\varXi _{t+1})$, thus denoting by $\overline{c}_t(s)$ the invasion speed of $(\varXi _t)_t$ starting with $n_0(\xi )=\alpha e^{-s\xi }$ for all $\xi \in {\mathbb R}^+$, we have

$$\begin{aligned} \overline{c}_{t+1}(s)&=\frac{\varXi _{t+1}-\varXi _0}{t+1}\\&=\frac{1}{t+1}\sum _{i=0}^t \frac{1}{s}\ln (r_iM(s))\\&=\left( \frac{1}{s}\ln (M(s))+\frac{1}{t+1}\sum _{i=0}^t \frac{1}{s}\ln (r_i)\right) . \end{aligned}$$

Thus, $(\overline{c}_t(s))$ is the sum of independent identically distributed variables and thus converges in distribution to a random variable that is normally distributed with mean $\mu (s)$ and variance $\sigma ^2(s)$ such that

$$\begin{aligned} \mu (s)=E[\frac{1}{s}\ln (r_0M(s))] \end{aligned}$$

(104)

and

$$\begin{aligned} \sigma ^2(s)=\underset{t\rightarrow +\infty }{\lim }\frac{1}{t} V\left[ \frac{1}{s}\ln (r_0M(s))\right] =0. \end{aligned}$$

(105)

As $\sigma ^2(s)\equiv 0$, this implies that $\overline{c}_t(s)$ converges in probability to the constant $\frac{1}{s}E[\ln (r_0M(s))]$. This is true for all s such that M(s) exists. Now if we want to consider the more general cases when $n_0$ is a compactly supported function, the minimal speed over all the s will be the relevant one and we have that the invasion speed of the stochastic wave at time t, $\overline{c}_t$, has mean $\mu =\underset{s>0}{\inf }\mu (s)$ and variance $(\sigma ^*_t)^2=\sigma ^2_t(s^*)$, where $s^*$ is such that $\mu (s^*)=\mu$, and thus converges in probability to

$$\begin{aligned} \overline{c}^*=\underset{s>0}{\inf }\frac{1}{s}E[\ln (r_0M(s))]. \end{aligned}$$

(106)