Stochastic Nicholson's blowflies delay differential equation with regime switching

TL;DR

This paper studies the existence and behavior of positive solutions to a stochastic Nicholson's blowflies delay differential equation with regime switching, extending previous results in the field.

Contribution

It provides new theoretical results on the global existence and path estimation of solutions for stochastic delay differential equations with regime switching.

Findings

Proved the global existence of positive solutions.

Derived estimates for the solution paths.

Extended previous results to more general stochastic models.

Abstract

In this paper, we investigate the global existence of almost surely positive solution to a stochastic Nicholson's blowflies delay differential equation with regime switching, and give the estimation of the path. The results presented in this paper extend some corresponding results in Wang et al. Stochastic Nicholson's Blowflies Delayed Differential Equations, Appl. Math. Lett. 87 (2019) 20-26.

Equations

$$X^{\prime }(t)=-\delta X(t)+pX(t-\tau )e^{-aX(t-\tau )}$$

$$dX(t)=[-\delta X(t)+pX(t-\tau )e^{-aX(t-\tau )}]dt+\sigma X(t)dB(t),$$

$$\mathbb{P}(r_{t+\vartriangle}=j|r_{t}=i)=\left\{\begin{array}[]{ll}&q_{ij}\vartriangle+o(\vartriangle),\;\mbox{for}\;i\not=j,\\ &1+q_{ij}\vartriangle+o(\vartriangle),\;\mbox{for}\;i=j,\end{array}\right.as\vartriangle\rightarrow 0^{+},$$

$$dX(t)=$$

$$X(s)=\varphi (s),\text{mbox}fors\in [-\tau ,0],\varphi \in C([-\tau ,0],{\mathbb{R}}^{+}),$$

$$\hat{g}=\underset{i\in \mathcal{S}}{\min }\{g_i\},\check{g}=\underset{i\in \mathcal{S}}{\max }\{g_i\},{\mathbb{R}}_0^{+}=[0,\infty ),{\beta }_{r_t}=\theta \left[{\delta }_{r_t}+(1-\theta ){\sigma }_{r_t}^2/2\right],{\gamma }_{r_t}=p_{r_t}(a_{r_t}e{)}^{-1},{\rho }_{r_t}={\delta }_{r_t}/{\sigma }_{r_t}^2,$$

$$M_{r_t}=\underset{x\in {\mathbb{R}}_0^{+}}{\max }\{\theta {\gamma }_{r_t}{x}^{\theta -1}-{\beta }_{r_t}{x}^{\theta }\},W_{r_t}=\underset{x\in {\mathbb{R}}_0^{+}}{\max }\{(\alpha -{\beta }_{r_t})x^{\theta }+\theta {\gamma }_{r_t}{x}^{\theta -1}\},0<\alpha <\hat{\beta },d_{r_t}={\delta }_{r_t}+{\sigma }_{r_t}^2/2,$$

$$\displaystyle\mu_{r_{t}}=2\delta_{r_{t}}-p_{r_{t}}-\sigma^{2}_{r_{t}},\,\;\mbox{and}\;\;C_{i}=\left\{\begin{array}[]{ll}\frac{\sigma^{2}_{i}-2\delta_{i}+\sqrt{(\sigma^{2}_{i}-2\delta_{i})^{2}+4\gamma^{2}_{i}}}{2},&\mbox{if}\;\;\delta_{i}\leqslant\sigma^{2}_{i}/2,\\ \frac{\gamma^{2}_{i}}{2\delta_{i}-\sigma^{2}_{i}},&\mbox{if}\;\;\delta_{i}>\sigma^{2}_{i}/2,\end{array}\right.\;\;i\in\mathcal{S}.$$

$$dX(t)=$$

$$...$$

$$dX(t)=$$

$$X(t)=\Phi_{r_{1}j_{1}}(t)\biggl{[}X(t_{0}+j_{1}\tau_{r_{1}})+\int_{t_{0}+j_{1}\tau_{r_{1}}}^{t}\Phi^{-1}_{r_{1}j_{1}}(s)G(X(s-\tau_{r_{1}}),r_{1})ds\biggr{]}>0\;\;a.s.$$

$$dV(X(t),r_t)=$$

$$\triangleq$$

$${\int }_0^{{\Lambda }_n\wedge t}dV(x(t),r_0)⩽{\int }_0^{{\Lambda }_n\wedge t}{M}_{r_0}dt+{\int }_0^{{\Lambda }_n\wedge t}\theta {\sigma }_{r_0}{X}^{\theta }(s)dB(s),$$

$$\underset{t\to \infty }{\mathrm{limsup}}\mathbb{E}[X^{\theta }(t)]⩽\hat{M},\underset{t\to \infty }{\mathrm{limsup}}\frac{1}{t}{\int }_0^t\mathbb{E}[X^{\theta }(s)]ds⩽\hat{W}/\alpha ,\text{mbox}where\theta \in [1,1+2\hat{\rho });\text{mbox}and$$

$$-\check{d}⩽\underset{t\to \infty }{\mathrm{liminf}}\frac{1}{t}\log X(t)⩽\underset{t\to \infty }{\mathrm{limsup}}\frac{1}{t}\log X(t)⩽\check{C}/2a.s.$$

$$\alpha {\int }_0^t\mathbb{E}[X^{\theta }(s)]ds⩽{\int }_0^t{W}_{r_s}ds+V(\varphi (0)),$$

$$\log (1+X^2(t))=\log (1+{\varphi }^2(0))+{\int }_0^{t}F(X(s),r_s)ds-2{\int }_0^t\frac{{\sigma }_{r_s}^2{X}^4(s)}{(1+X^2(s){)}^2}ds+M(t),$$

$$F(X(t),r_t)⩽[({\sigma }_{r_t}^2-2{\delta }_{r_t})X^2(t)+2{\gamma }_{r_t}X(t)]/(1+X^2(t))⩽C_{r_t}.$$

$$\log (1+X^2(t))⩽\log (1+{\varphi }^2(0))+{\int }_0^t{C}_{r_s}ds+\log \sqrt{k}.$$

$$\frac{1}{t}\log (1+X^2(t))⩽\frac{1}{k-1}\left[\log (1+{\varphi }^2(0))+\log \sqrt{k}\right]+\frac{1}{t}{\int }_0^t{C}_{r_s}ds.$$

$$\underset{t\to \infty }{\mathrm{limsup}}\frac{1}{t}\log (X^2(t))⩽\underset{t\to \infty }{\mathrm{limsup}}\frac{1}{t}\log (1+X^2(t))⩽\underset{t\to \infty }{\mathrm{limsup}}\frac{1}{t}{\int }_0^t{C}_{r_s}ds⩽\check{C}a.s.$$

$$\log X(t)=\log \varphi (0)-{\int }_0^t\left({\delta }_{r_s}+\frac{1}{2}{\sigma }_{r_s}^2\right)ds+{\int }_0^t{p}_{r_s}\frac{X(s-{\tau }_{r_s})e^{-a_{r_s}X(s-{\tau }_{r_s})}}{X(s)}ds+{\int }_0^t{\sigma }_{r_s}dB(s),$$

$$\underset{t\to \infty }{\mathrm{limsup}}\frac{1}{t}{\int }_0^t\mathbb{E}[X^{\theta }(s)]ds⩽\frac{1}{\alpha }\sum _{i\in \mathcal{S}}{\pi }_i{W}_i,\text{mbox}where\theta \in [1,1+2\hat{\rho });\text{mbox}and$$

$$-\sum _{i\in \mathcal{S}}{\pi }_i{d}_i⩽\underset{t\to \infty }{\mathrm{liminf}}\frac{1}{t}\log X(t)⩽\underset{t\to \infty }{\mathrm{limsup}}\frac{1}{t}\log X(t)⩽\frac{1}{2}\sum _{i\in \mathcal{S}}{\pi }_i{C}_{i}a.s.$$

$$0⩽\underset{t\to \infty }{\mathrm{liminf}}X(t)⩽\sum _{i\in \mathcal{S}}{\pi }_i{\gamma }_i/\sum _{i\in \mathcal{S}}{\pi }_i{d}_i:=N^{\ast }a.s.$$

$$\log X(t)\leqslant\log\phi(0)+\int_{0}^{t}\biggl{(}-d_{r_{s}}+\frac{\gamma_{r_{s}}}{X(s)}\biggr{)}ds+\int_{0}^{t}\sigma_{r_{s}}dB(s)\;\;a.s.$$

$$\limsup_{t\rightarrow\infty}\frac{1}{t}\log X(t)\leqslant\limsup_{t\rightarrow\infty}\frac{1}{t}\int_{0}^{t}\biggl{(}-d_{r_{s}}+\frac{\gamma_{r_{s}}}{N^{*}+\varepsilon}\biggr{)}ds=\sum_{i\in\mathcal{S}}\pi_{i}\left(-d_{i}+\frac{\gamma_{i}}{N^{*}+\varepsilon}\right)<0\;\;a.s.$$

$$\underset{t\to \infty }{\mathrm{limsup}}\frac{1}{t}\log X(t)⩽-\sum _{i\in \mathcal{S}}{\pi }_i{d}_i<0\text{mbox}and\underset{t\to \infty }{\lim }X(t)=0a.s.;$$

$$\underset{t\to \infty }{\mathrm{limsup}}\frac{1}{t}\log X(t)⩽-\sum _{i\in \mathcal{S}}{\pi }_i(d_i-p_i):=-\lambda ,$$

$$\underset{t\to \infty }{\mathrm{limsup}}\frac{1}{t}\log \mathbb{E}[X^2(t)]⩽-\kappa <0,\underset{t\to \infty }{\mathrm{limsup}}\frac{1}{t}\log X(t)⩽-\kappa /2<0\text{mbox}and\underset{t\to \infty }{\lim }X(t)=0a.s.,$$

$$d[e^{\epsilon t}V(X(t),t)]⩽$$

Full text

Stochastic Nicholson’s blowflies delay differential equation with regime switching∗

** Yanling Zhu†, Kai Wang*†,∗, Yong Ren‡,∗, Yingdong Zhuang†*

†*** School of Statistics and Applied Mathematics, Anhui University of Finance and Economics,

Bengbu 233030, Anhui, China

‡ Department of Mathematics, Anhui Normal University, Wuhu 241000, Anhui, China*****

Abstract  In this paper, we investigate the global existence of almost surely positive solution to a stochastic Nicholson’s blowflies delay differential equation with regime switching, and give the estimation of the path. The results presented in this paper extend some corresponding results in Wang et al. Stochastic Nicholson’s blowflies delayed differential equations, Appl. Math. Lett. 87 (2019) 20-26.

Keywords  Nicholson’s blowflies model, Regime switching, Global solution

AMS(2000)  34K50;   60H10;   65C30

1 Introduction

The following delay differential equation, known as the famous Nicholson’s blowflies model

[TABLE]

was presented by Gurney et al. (1990) [1] to model the population of Lucilia cuprina(the Australian sheep-blowfly). The biological meaning of the parameters in this model is that: $\delta$ denotes the adult death rate of per capita daily, $p$ denotes the maximum egg production rate of per capita daily, $\tau$ denotes the generation time, and $1/a$ is the size that the Lucilia cuprina reproduces at its maximum rate. Since then, the dynamics of this model was investigated by many researchers, see Berezansky et al. (2010) [2] for an overview.

Considering that the per capita daily adult death rate $\delta$ is affected by the white noise of the environment, one often uses $\delta-\sigma dB(t)$ in lieu of $\delta$, where $\sigma$ denotes the intensity of the white noise and $B(t)$ is a one-dimensional Brownian motion defined on $(\Omega,\{\mathfrak{F}_{t}\}_{t\geqslant 0},\mathbb{P})$, which is a complete probability space. Wang et al. (2019)[3] studied the following stochastic Nicholson’s blowflies delay differential equation:

[TABLE]

with initial conditions $X(s)=\phi(s),\;\mbox{for}\;s\in[-\tau,0],\phi\in C([-\tau,0],[0,+\infty)),\phi(0)>0.$ Under the condition $\delta>\sigma^{2}/2$, they prove that equation (1.1) has a unique global solution on $[-\tau,\infty)$, which is positive a.s. on $[0,\infty)$, and give the estimations of $\limsup_{t\rightarrow\infty}\mathbb{E}[X(t)]$ and $\limsup_{t\rightarrow\infty}\frac{1}{t}\int_{0}^{t}\mathbb{E}[X(s)]ds$, respectively.

But we find from numerical simulation that the condition $\delta>\sigma^{2}/2$ is too strict to the nonexplosion of the solution to equation (1.1). For example, let $\delta=1$, $p=5$, $\tau=1$, $a=1$, $\sigma=2$, then the numerical simulation in Figure 1(left) shows that the solution $X(t)$ is globally existent on $[0,+\infty)$, and oscillating about $N^{*}=\log p/\delta$, which is the positive equilibrium of the corresponding deterministic equation of (1.1), but $\delta<\sigma^{2}/2$.

Meanwhile, as one knows that besides the white noises there is a so called telegraph noise in the real ecosystem. A switching process between two or more regimes can be used to illustrate this type of noise. One can model the regime switching by a continuous time Markov chain $(r_{t})_{t\geqslant 0}$, which takes values in a finite state space $\mathcal{S}=\{1,2,...,m\}$ and satisfies

[TABLE]

where $Q=(q_{ij})\in R^{m\times m}$ is the infinitesimal generator, and $q_{ij}\geqslant 0$ denotes the transition rate from $i$ to $j$ for $i\not=j$, and $\sum_{j\in\mathcal{S}}q_{ij}=0$, $\forall\,i\in\mathcal{S}$. In consideration of that the white noise and the telegraph noise are deferent types of noise, we assume that the Markov chain $(r_{t})_{t\geqslant 0}$ is independent of the Brownian motion $B(t)$. For more details on the theorem of regime switching system, one can refer to Mao et al. (2006) [4], Khasminskii et al. (2007) [5] and Yin et al. (2009) [6].

Considering both the effects of white noise and telegraph noise of environment, we present the following stochastic Nicholson’s blowflies delay differential equation with regime switching

[TABLE]

with initial conditions:

[TABLE]

where $\tau=\max_{i\in S}\{\tau_{i}\}$, $\mathbb{R}^{+}=(0,+\infty)$. The contribution of this paper is as follows:

$\bullet$ Model (1.2) considered in this paper includes regime switching, which is more general than model (1.1);

$\bullet$ We show the global existence of almost surely positive solution to equation (1.2) without the restricted condition $\delta_{i}>\sigma_{i}^{2}/2$, $i\in\mathcal{S}$, which extends the one in [3];

$\bullet$ We give the estimations of $\limsup_{t\rightarrow\infty}\mathbb{E}[X^{\theta}(t)]$, $\limsup_{t\rightarrow\infty}\frac{1}{t}\int_{0}^{t}\mathbb{E}[X^{\theta}(s)]ds$, $\liminf_{t\rightarrow\infty}X(t)$ and the sample Lyapunov exponent of $X(t)$, which are more general than the ones in literature.

For simplicity, in the following and next sections we use the following notations:

[TABLE]

Assumption (A): The Markov chain $(r_{t})_{t\geqslant 0}$ is irreducible with an invariant distribution $\pi=(\pi_{i},i\in\mathcal{S})$.

Before the main results, we first give two lemmas and omit the proof here for saving layouts.

Lemma 1.1**.**

Let $f(x)=xe^{-ax},\,a\in\mathbb{R}^{+},x\in\mathbb{R}^{+}_{0}$, then $f(x)\leqslant(ae)^{-1}$.

Lemma 1.2**.**

For $a\in\mathbb{R}$, $b\in\mathbb{R}^{+}$, we have $\frac{ax^{2}+bx}{1+x^{2}}\leqslant C(a)$ for $x\in\mathbb{R},$ where C(a)=\left\{\begin{array}[]{ll}({a+\sqrt{a^{2}+b^{2}}})/{2},&a\geqslant 0;\\ -{b^{2}}/{4a},&a<0.\end{array}\right.

2 Main Results

In this section, we will show the global existence of almost surely positive solution of equation (1.2) with initial value (1.3), and give some estimations of the solution $X(t)$.

Theorem 2.1**.**

For any given initial value (1.3), equation (1.2) has a unique solution $X(t)$ on $[-\tau,\infty)$, which is positive on $\mathbb{R}^{+}_{0}$ a.s.

Proof.

Since all the coefficients of equation (1.2) are locally Lipschitz continuous on $\mathbb{R}^{+}_{0}$, then for any given initial value (1.3), there is a unique maximal local solution $X(t)$ on $[-\tau,\Lambda_{e})$, where $\Lambda_{e}$ is the explosion time. Assume that the subsystems of (1.2) switch along the Markov chain $r_{0},r_{1},r_{2},...,r_{n}$, and the times on each state are $[0,t_{0}],[t_{0},t_{1}],[t_{1},t_{2}],...,$ $[t_{n-1},\Lambda_{e})$, respectively.

Step 1. We first show that $X(t)$ is positive almost surely on $[0,\Lambda_{e})$. System (1.2) switching along the Markov chain $r_{0},r_{1},r_{2},...,r_{n}$ gives

[TABLE]

where $G(X(t-\tau_{r_{i}}),r_{i})=p_{r_{i}}X(t-\tau_{r_{i}})e^{-a_{r_{i}}X(t-\tau_{r_{i}})}$. Obviously, each equation is a linear stochastic differential equation, which can be solved step by step, see Mao [7] (p. 98-99) for more details.

In fact, for $t\in[0,\tau_{r_{0}}]\subset[0,t_{0}]$, the solution of the first equation of (2.4) is

\;\;\;\;X(t)=\Phi_{r_{0}0}(t)\biggl{[}\phi(0)+\int_{0}^{t}\Phi^{-1}_{r_{0}0}(s)G(X(s-\tau_{r_{0}}),r_{0})ds\biggr{]}>0\;\;a.s.,

where $\Phi_{r_{0}i}(t)=e^{-d_{r_{0}}(t-i\tau_{r_{0}})+\sigma_{r_{0}}[B(t)-B(i\tau_{r_{0}})]}$, and for $t\in[\tau_{r_{0}},2\tau_{r_{0}}]\subset[\tau_{r_{0}},t_{0}]$ we have

\;\;\;\;X(t)=\Phi_{r_{0}1}(t)\biggl{[}X(\tau_{r_{0}})+\int_{\tau_{r_{0}}}^{t}\Phi^{-1}_{r_{0}1}(s)G(X(s-\tau_{r_{0}}),r_{0})ds\biggr{]}>0\;\;a.s.,

$\quad\qquad\;\,...$

and for $t\in[j_{0}\tau_{r_{0}},t_{0}]$, where $j_{0}=\lfloor\frac{t_{0}}{\tau_{r_{0}}}\rfloor$, we get

\;\;\;\;X(t)=\Phi_{r_{0}j_{0}}(t)\biggl{[}X(j_{0}\tau_{r_{0}})+\int_{j_{0}\tau_{r_{0}}}^{t}\Phi^{-1}_{r_{0}j_{0}}(s)G(X(s-\tau_{r_{0}}),r_{0})ds\biggr{]}>0\;\;a.s.,

For $t\in[t_{0},t_{0}+\tau_{r_{1}}]\subset[t_{0},t_{1}]$, the solution of the second equation of (2.4) is

\;\;\;\;X(t)=\Phi_{r_{1}0}(t)\biggl{[}X(t_{0})+\int_{t_{0}}^{t}\Phi^{-1}_{r_{1}0}(s)G(X(s-\tau_{r_{1}}),r_{1})ds\biggr{]}>0\;\;a.s.,

where $\Phi_{r_{1}i}(t)=e^{-d_{r_{1}}[t-(t_{0}+i\tau_{r_{1}})]+\sigma_{r_{1}}[B(t)-B(t_{0}+i\tau_{r_{1}})]}$, and, for $t\in[t_{0}+\tau_{r_{1}},t_{0}+2\tau_{r_{1}}]\subset[t_{0}+\tau_{r_{1}},t_{1}]$,

\;\;\;\;X(t)=\Phi_{r_{1}1}(t)\biggl{[}X(t_{0}+\tau_{r_{1}})+\int_{t_{0}+\tau_{r_{1}}}^{t}\Phi^{-1}_{r_{1}1}(s)G(X(s-\tau_{r_{1}}),r_{1})ds\biggr{]}>0\;\;a.s.,

$\quad\qquad\;\,...$

and for $t\in[t_{0}+j_{1}\tau_{r_{1}},t_{1}]$ with $j_{1}=\lfloor\frac{t_{1}-t_{0}}{\tau_{r_{1}}}\rfloor$,

[TABLE]

Repeating this procedure, we can obtain the explicit solution $X(t)$ of (2.4), which is positive on $[0,\Lambda_{e})$ a.s.

Step 2. Now, we prove that this solution $X(t)$ to equation (1.2) with (1.3) is global existence, that is, the solution will not explode in finite time. We only need to show $\Lambda_{e}=\infty$ a.s. Let $k_{0}>0$ be sufficiently large that $\max_{t\in[-\tau,0]}\phi(t)<k_{0}$, and for each integer $k\geqslant k_{0}$ define the stopping time $\Lambda_{k}=\inf\{t\in[0,\Lambda_{e})|\;X(t)\geqslant k\}.$ It is clear that $\Lambda_{k}$ is increasing as $k\rightarrow\infty.$ Let $\Lambda_{\infty}=\lim_{k\rightarrow\infty}\Lambda_{k}$, whence $\Lambda_{\infty}\leqslant\Lambda_{e}$ a.s. In order to show $\Lambda_{e}=\infty$, we only need to prove that $\Lambda_{\infty}=\infty$ a.s. If it is false, then there is a pair of constants $T>0$ and $\epsilon\in(0,1)$ such that $\mathbb{P}\{\Lambda_{\infty}\leqslant T\}>\epsilon,$ which yields that there exists an integer $k_{1}\geqslant k_{0}$ such that $\mathbb{P}\{\Lambda_{k}\leqslant T\}\geqslant\epsilon\;\;\mbox{for all}\;\;k\geqslant k_{1}.$

Let $V(x,i)=x^{\theta}$ with $\theta\in[1,1+2\rho)$, then by Itô’s formula, we have

[TABLE]

It follows from $X(t)>0$ a.s. on $[0,\Lambda_{e})$, $\beta_{r_{t}}>0$ and Lemma 1.1 that $\mathcal{L}V(X(t),r_{t})\leqslant-\beta_{r_{t}}X^{\theta}(t)+\theta\gamma_{r_{t}}X^{\theta-1}(t)\leqslant M_{r_{t}}.$ Thus, we have $dV(X(t),r_{t})\leqslant M_{r_{t}}dt+\theta\sigma_{r_{t}}X^{\theta}(t)dB(t).$ For $t\in[0,\tau_{r_{0}}]\subset[0,t_{0}]$, integrating both sides of above inequality on $[0,\Lambda_{n}\wedge t]$ and taking expectation lead to

[TABLE]

then $\mathbb{E}[V(\Lambda_{n}\wedge t,r_{0})]\leqslant V(\phi(0))+M_{r_{0}}\mathbb{E}[(\Lambda_{n}\wedge t)]\leqslant V(\phi(0))+M_{r_{0}}\tau_{r_{0}},$ which yields $\Lambda_{\infty}\geqslant\tau_{r_{0}}$ a.s. Similarly, we can obtain that $\Lambda_{\infty}\geqslant t_{0}$. Repeating this procedure, we can show $\Lambda_{\infty}\geqslant t_{1}$,…, and $\Lambda_{\infty}\geqslant t_{n-1}+j_{n}\tau_{r_{n}}$ a.s. for any integer $j_{n}$, which implies $\Lambda_{\infty}=\infty$ a.s. ∎

Theorem 2.2**.**

The solution $X(t)$ of equation (1.2) with initial value (1.3) has the following properties:

[TABLE]

Proof.

It follows from It$\hat{o}$ formula that $\mathbb{E}[e^{t}V(X(t))]\leqslant V(\phi(0))+\int_{0}^{t}\hat{M}e^{s}ds=V(\phi(0))+\hat{M}(e^{t}-1),$ which yields $\limsup_{t\rightarrow\infty}\mathbb{E}[X^{\theta}(t)]\leqslant\hat{M}\;\,a.s.$ From equation (2.5) we have $dV(t)+\alpha X^{\theta}(t)dt\leqslant W_{r_{t}}dt+\theta\sigma_{i}X^{\theta}(s)dB(s).$ Integrating this inequality from [math] to $t$ and taking expectation on both sides of it gives

[TABLE]

which implies $\limsup_{t\rightarrow\infty}\frac{1}{t}\int_{0}^{t}\mathbb{E}[X^{\theta}(s)]ds\leqslant\hat{W}/\alpha\;\,a.s.$ Applying It$\hat{o}$ formula we have

[TABLE]

where $F(X(t),r_{t})=\frac{1}{1+X^{2}(t)}[(\sigma^{2}_{r_{t}}-2\delta_{r_{t}})X^{2}(t)+2p_{r_{t}}X(t)G(t-\tau_{r_{t}},r_{t})]$ and $M(t)=2\int_{0}^{t}\frac{\sigma_{r_{s}}X^{2}(s)}{1+X^{2}(s)}dB(s)$. It follows from Lemmas 1.1-1.2 that

[TABLE]

Meanwhile, the exponential martingale inequality yields, for any $k\in\mathbb{N}$, $\mathbb{P}\{\sup_{0\leqslant t\leqslant k}[M(t)-2\langle M(t),M(t)\rangle]\geqslant\log\sqrt{k}\}\leqslant 1/k^{2}.$ It follows from $\sum_{k}k^{-2}<\infty$ and the Borel-Cantelli lemma that there exists an $\Omega_{0}\subset\Omega$ with $\mathbb{P}(\Omega_{0})=1$ such that for every $\omega\in\Omega_{0}$ there is an integer $k_{0}=k_{0}(\omega)$ such that, for all $k\geqslant k_{0}$ and $0\leqslant t\leqslant k$, we get $\sup_{0\leqslant t\leqslant k}\left[M(t)-2\langle M(t),M(t)\rangle\right]\leqslant\log\sqrt{k}.$ Substituting this inequality and (2.7) into equation (2.6), for all $\omega\in\Omega_{0}$ and $0\leqslant t\leqslant k$, we obtain

[TABLE]

Thus, for all $\omega\in\Omega_{0}$ and $t\in[k-1,k]$ we have

[TABLE]

By letting $t\rightarrow\infty$, we get

[TABLE]

This leads to the desired assertion. It follows from It$\hat{o}$ formula and equation (1.3) that

[TABLE]

which together with the large number theorem for martingales yields $\liminf_{t\rightarrow\infty}\frac{1}{t}\log X(t)\geqslant-\check{d}.$ ∎

Theorem 2.3**.**

Assume (A) holds, then the solution $X(t)$ of equation (1.2) with (1.3) has the properties:

[TABLE]

[TABLE]

Proof.

The Markov $r_{t}$ has an invariant distribution $\pi=(\pi_{i},i\in\mathcal{S})$, thus we have $\limsup_{t\rightarrow\infty}\frac{1}{t}\int_{0}^{t}W_{r_{s}}ds=\sum_{i\in\mathcal{S}}\pi_{i}W_{i},$ $\limsup_{t\rightarrow\infty}\frac{1}{t}\int_{0}^{t}C_{r_{s}}ds=\sum_{i\in\mathcal{S}}\pi_{i}C_{i}$ and $\limsup_{t\rightarrow\infty}\frac{1}{t}\int_{0}^{t}(\delta_{r_{s}}+\frac{1}{2}\sigma^{2}_{r_{s}})ds=\sum_{i\in\mathcal{S}}\pi_{i}d_{i}.$ ∎

Remark 2.1**.**

These estimations presented in Theorem 2.3 are dependent on the stationary distribution $\pi$ of the Markov chain $(r_{t})_{t\geqslant 0}.$ Especially, these estimations are also suitable for the case of $\delta\leqslant\frac{1}{2}\sigma^{2}$, which fills in the gap of the corresponding results in Wang et.al [3].

Theorem 2.4**.**

Assume (A) holds, then the solution $X(t)$ of model (1.2) with initial value (1.3) satisfies

[TABLE]

Remark 2.2**.**

Inequality (2.4) means that if the daily egg production rate goes to 0, then the limit inferior of the population of the sheep-blowfly will go to 0. In fact, one can see from (2.4) that $p_{i}\rightarrow 0$, $\forall\,i\in\mathcal{S}$ lead to $N^{*}\rightarrow 0$, which yields $\inf X(t)\rightarrow 0$ as $t\rightarrow\infty$ a.s.

Proof.

If $\liminf_{t\rightarrow\infty}X(t)>N^{*}$ a.s., then there exists a positive constant $\varepsilon$ such that $\liminf_{t\rightarrow\infty}X(t)=N^{*}+2\varepsilon$ a.s., and for this $\varepsilon$ there exists a constant $T>0$ so that $X(t)\geqslant N^{*}+\varepsilon$ for all $t>T$ a.s. It follows from the It$\hat{o}$ formula and equation (1.2) with (1.3) that

[TABLE]

Thus

[TABLE]

which yields $\limsup_{t\rightarrow\infty}X(t)=0$, a contradiction. ∎

Theorem 2.5**.**

Assume (A) holds, then the solution $X(t)$ of equation (1.2) with initial value (1.3) has the properties: if $\forall\,i\in\mathcal{S}$, $p_{i}=0$, then

[TABLE]

if $\tau_{i}=0$, $\forall\,i\in\mathcal{S}$, then

[TABLE]

further if $\lambda>0$, then $\lim_{t\rightarrow\infty}X(t)=0$; and if $\lambda=0$, then $\limsup_{t\rightarrow\infty}X(t)\leqslant 1$ a.s.

Theorem 2.6**.**

Assume $\forall\,i\in S$, $\tau_{i}=\tau$ and $a_{i}=a$ are independent of the Markov chain $(r_{t})_{t\geqslant 0}$, and $\hat{\mu}>0$, then the solution $X(t)$ of equation (1.2) with initial value (1.3) has the properties:

[TABLE]

where $\kappa\in(0,\hat{\mu})$ is the unique root of the equation $\kappa\vartheta\tau e^{\kappa\tau}+\kappa=\hat{\mu},$ and $\vartheta>\check{p}$ is a constant.

Proof.

Let $e^{\kappa t}V(x,t)=e^{\kappa t}[x^{2}+\vartheta\int^{t}_{t-\tau}g(x(s))ds]$, where $g(x)=x^{2}e^{-2ax}$, by applying It$\hat{o}$ formula we have

[TABLE]

Integrating two sides of this inequality from [math] to $T>0$ and taking expectation we obtain

[TABLE]

On the other hand, we have

[TABLE]

which together with (2.9) gives

[TABLE]

Hence,

[TABLE]

which implies $\limsup_{t\rightarrow\infty}\frac{1}{t}\log(\mathbb{E}[X^{2}(t)])\leqslant-\kappa.$ The rest of the proof uses the same arguments as those in the proof of Mao [7] (p.373-374). For completeness, we include the details here. For arbitrary $\epsilon\in(0,\kappa/2)$, there must exist a constant $K>0$ so that

[TABLE]

For $n=1,2,...$, by the Doob martingale inequality and Holder inequality, we get

[TABLE]

where $C$ is positive constant independent of $n$, and allowed to be different in deferent lines. Therefore, we obtain $\mathbb{E}[\sup_{n\leqslant t\leqslant n+1}X^{2}(t)]\leqslant Ce^{-(\kappa-\epsilon)n},$ which together with Chebyshev’s theorem yields

[TABLE]

It follows from the Borel-Cantelli lemma that for almost all $\omega\in\Omega$ there is a random integer $n_{0}(\omega)$ such that for all $n\geqslant n_{0},$

[TABLE]

which yields $\limsup_{t\rightarrow\infty}\frac{1}{t}\log X(t)\leqslant-\kappa/2+\epsilon$ a.s. The arbitarity of $\epsilon$ implies the desired assertion. ∎

3 Examples

In this section an example is given to check the results obtained in section 2 by numerical simulation. Let $\mathcal{S}=\{1,2,3\}$, $Q=[-10,4,6;2,-3,1;3,5,-8]$, $\delta=[2,1,4]$, $p=[4,2,8]$, $\tau=[1,1,1]$, $a=[0.4,0.2,0.3]$, $\sigma=[1.5,2,3]$, and $\phi(t)=1$, $t\in[-1,0]$, then the Markov chain $r_{t}$ is irreducible and has a unique stationary distribution $\pi=(0.1845,0.6019,0.2136)$ and one can see that $\delta_{1}>\sigma_{1}^{2}/2,$ $\delta_{2}=\sigma_{2}^{2}/2,$ and $\delta_{3}<\sigma_{3}^{2}/2$. It follows from Theorem 2.1 that the solution of equation (1.2) is global existence and positive almost surely on $\mathbb{R}^{+}_{0}$, Figure 2(left) shows this. Meanwhile, we get $\limsup_{t\rightarrow\infty}\mathbb{E}[X(t)]\leqslant 1.1221$ and $\limsup_{t\rightarrow\infty}\mathbb{E}[X^{7/5}(t)]\leqslant 0.3713$. Furthermore, from Theorem 2.3 we have $C\thickapprox[2.1022,3.6788,10.3229]$, $d=[3.125,3,8.5]$, $N^{*}=1.1883$ and the estimation of the sample Lyapunov exponent: $-4.1978\leqslant\liminf_{t\rightarrow\infty}\frac{1}{t}{\log X(t)}\leqslant\limsup_{t\rightarrow\infty}\frac{1}{t}{\log X(t)}\leqslant 2.4035$. One can see from Figure 2(right) that the value of $\frac{1}{t}\log X(t)$ oscillates up and down at zero as $t\rightarrow\infty$, which yields that the sample path of $X(t)$ oscillates up and down at 1 as $t\rightarrow\infty$. Thus, $X(t)\nrightarrow 0$ as $t\rightarrow\infty$, that is, the population will neither go to extinction finally nor explode in finite time. In order to verify the effectiveness of Theorem 2.6, we reset $p=[0.2,0.2,0.4]$, $a=[0.4,0.4,0.4]$, $\sigma=[1.5,1,2.5]$ then $\hat{\mu}>0$ holds. One can see from Figure 1(right) that $X(t)\rightarrow 0$ and $\frac{1}{t}\log X(t)<0$ as $t\rightarrow\infty.$

Figure 1.   Left: The solution $X(t)$ of equation (1.1) with initial value $\phi(t)=1$ for $t\in[-1,0]$ and the step $\Delta=1e-4$. Right(up): Simulation of the processes $(X(t),r_{t})$ with $r_{0}=3$ and $\hat{\mu}>0$. Right(down): $[\log X(t)]/t$ of the sample path $X(t)$.

Figure 2.   Here the step $\Delta=1e-4$. Left: Simulation of the processes $(X(t),r_{t})$ with $r_{0}=3,$ and $\phi(t)=1$ for $t\in[-1,0]$ (In order to see the line clearly, we only present the first $5\times 10^{3}$ states of $r_{t}$, and just show the values of $X(t)$ on interval $[800,1000]$). Right: $[\log X(t)]/t$ of the sample path $X(t)$ (For clarity, we just print the values of $[\log X(t)]/t$ on interval $[800,1000]$).

4 Conclusion

In this paper, we consider a stochastic Nicholson’s blowflies delay differential equation with regime switching. To the best of our knowledge it is the first time that regime switching is considered in Nicholson’s blowflies model. By choosing several suitable $V$ functions we prove that the solution $X(t)$ of the stochastic Nicholson’s blowflies model is globally existent, which is also positive on $\mathbb{R}^{+}_{0}$ almost surely, and present some estimation of the solution $X(t)$ to model (1.2)-(1.3). The results obtained in this paper extend some corresponding results in literature, especially we removed the restricted condition $\delta>\sigma^{2}/2$, which is needed in [3]. It is worth to point out that model (1.2) can be extended to a general case:

[TABLE]

where $G:\mathbb{R}\times\mathcal{S}\rightarrow\mathbb{R}^{+}_{0}$ is uniformly bounded, and the methods applied in Theorems 2.1$-$2.5 are also suitable for it. For model (1.2) with initial value (1.3), there are still some interesting problems left:

1. Find conditions for the almost surely strong (weak) persistence of the population, that is, under what conditions $\liminf_{t\rightarrow\infty}$ $X(t)>0$ a.s. ($\limsup_{t\rightarrow\infty}X(t)>0$ a.s.); or more accurately, estimate the maximum value of $K\in\mathbb{R}^{+}$ such that $\liminf_{t\rightarrow\infty}X(t)\geqslant K$ a.s.
2. In the case of $p>0$, $\tau>0$, the delay $\tau$ and parameter $a$ vary along the Markov chain $(r_{t})_{t\geqslant 0}$, that is, they depend on the Markov chain $(r_{t})_{t\geqslant 0}$. Find conditions for the upper boundedness of the population, that is, under what conditions there exists a $M\in\mathbb{R}^{+}$ such that $\limsup_{t\rightarrow\infty}X(t)\leqslant M$ a.s., and give the estimation of this upper bound $M$. Find conditions for the negativeness of the Lyapunov exponent of $X(t)$, that is $\limsup_{t\rightarrow\infty}\frac{1}{t}\log X(t)<0$ a.s.