Deciphering role of inter and intracity human dispersal on epidemic spread via coupled reaction-diffusion models

Abstract

Human mobility has been significantly influencing public health since time immemorial. A susceptible-infected-deceased epidemic reaction diffusion network model using asymptotic transmission rate is proposed to portray the spatial spread of the epidemic among two cities due to population dispersion. Qualitative behaviour including global attractor and persistence property are obtained. We also study asymptotic behaviour of the whole network with the help of asymptotic behaviour at individual cities. The epidemic model shows up two equilibria, (i) the disease-free, and (ii) unique endemic equilibria. An expression that can be used to calculate the basic reproduction number for heterogeneous environment, for the entire network is obtained. We use graph theory to analyze the global stability of our diffusive two-city model. We also performed bifurcation analysis and discovered that endemic equilibrium changes stability via Hopf bifurcations. A significant reduction in the number of infectives were observed when proper migration rate is maintained between the cities. Numerical results are provided to illuminate and clarify theoretical findings. Simulation experiments for two-dimensional spatial models show that infectious populations will increase if contact heterogeneity is increased, but it will decline if infective populations perform more local random movement. We observe that infection risk may be understated if the parameters used to estimate the basic reproduction number remains unchanged through space or time.

Appendix

Proof of Proposition 5

We consider the following auxiliary system

$$\begin{aligned} -D_{S_1} \nabla ^2 S_1= & {} (\theta r_1(x)+(1-\theta ) r_{10}) S_1 \left( 1-\frac{S_1+I_1}{K_1}\right) \nonumber \\&-\frac{(\theta \beta _1(x)-(1-\theta )\beta _{10}) S_1 I_1}{S_1+I_1+c_1}, \end{aligned}$$

(74)

$$\begin{aligned} -D_{I_1} \nabla ^2 I_1= & {} \frac{(\theta \beta _1(x)-(1-\theta )\beta _{10}) S_1 I_1}{S_1+I_1+c_1}-(\theta a_1(x)-(1-\theta )a_{10}) I_1, \nonumber \\ \frac{\partial S_1}{\partial \nu }= & {} \frac{\partial I_1}{\partial \nu }=0. \end{aligned}$$

(75)

where $r_{10}, \beta _{10}, a_{10}$ are positive constant and the parameter $\theta \in [0,1]$. Problem (74–75) becomes problem (47) at $\theta =1$. We divide the proof into three parts for easy understanding.

Step 1. We find the upper bounds for any positive solution $(S_1,I_1)$ to (74–75). In view of (74), it holds

$$\begin{aligned} \int _{\varOmega }S_1 dx\le K_1 ~~~\text {and}~~~ \int _{\varOmega }I_1 dx \le K_1. \end{aligned}$$

(76)

Thus, we can find a positive constant C independent of $\theta \in [0,1]$ such that

$$\begin{aligned} (\theta r_1(x)+(1-\theta ) r_{10}) S_1 \left( 1-\frac{S_1+I_1}{K_1}\right) \le max(max ~~r_1(x),r_{10})\int _{\varOmega }S_1 dx \le C \end{aligned}$$

and

$$\begin{aligned} (\theta \beta _1(x)+(1-\theta ) \beta _{10}) \frac{S_1 I_1}{S_1+I_1+c_1} \le max(max~~\beta _1(x),\beta _{10})\int _{\varOmega }I_1 dx \le C \end{aligned}$$

The positive constant C does not depend on the parameter $\theta \in [0,1]$ and can be different depending on its position. Applying $L^1$ estimate theory for elliptic equations [5] to Equ. (74–75), we obtain $||S_1||_{W^{1,1({\varOmega })}} \le C$ and $||I_1||_{W^{1,1({\varOmega })}} \le C$. Application of Sobolev embedding theorem gives us,

$$\begin{aligned} W^{1,1({\varOmega })} \rightarrow L^p({\varOmega }),~~~ \forall 1\le p \le \frac{n}{n-1}~~\text {or}~~1\le p<\infty ~~\text {if}~~n=1. \end{aligned}$$

we have

$$\begin{aligned} ||S_1||_{L^p({\varOmega })},||I_1||_{L^p({\varOmega })}\le C, ~~~ \forall 1\le p \le \frac{n}{n-1}~~\text {or}~~1\le p<\infty ~~\text {if}~~n=1. \end{aligned}$$

(77)

Applying $L^p$ estimate for elliptic equations [12] to (74–75) leads to

$$\begin{aligned}&||S_1||_{W^{2,p}({\varOmega })},||I_1||_{W^{2,p}({\varOmega })} \le C, \nonumber \\&\quad \forall 1\le p \le \frac{n}{n-1}~~\text {or}~~1\le p<\infty ~~\text {if}~~n=1. \end{aligned}$$

(78)

Again we apply Sobolev embedding theorem, to get

$$\begin{aligned}&||S_1||_{L^{p^*}({\varOmega })},||I_1||_{L^{p^*}({\varOmega })}\le C, \nonumber \\&\quad \forall 1\le p^* \le \frac{n}{n-3}~~\text {or}~~1\le p^*<\infty ~~\text {if}~~n=1. \end{aligned}$$

(79)

Repeating the above process finitely many times, one can affirm that

$$\begin{aligned} ||S_1||_{L^\infty ({\varOmega })}\le C, ||I_1||_{L^\infty ({\varOmega })}\le C. \end{aligned}$$

(80)

Step 2. Now, we find lower bounds for any positive solution $(S_1,I_1)$ to (74–75). Integrating (75) over ${\varOmega }$ gives

$$\begin{aligned} \int _{\varOmega }[\theta \beta (x)+(1-\theta )\beta _{10}] \frac{S_1 I_1}{S_1+I_1+c_1} dx=\int _{\varOmega }[\theta a(x)+(1-\theta ) a_{10}]I_1 dx, \end{aligned}$$

(81)

Clearly, (81) indicates

$$\begin{aligned} c \int _{\varOmega }I_1 dx\le d \int _{\varOmega }S_1 dx, \end{aligned}$$

(82)

where $c=min(min~~a(x), a_{10})>0, d=max(max ~\beta (x), \beta _{10})>0$. One can then insert $\int _{\varOmega }I_1 dx \le K_1-\int _{\varOmega }S_1 dx$ into (82) to get

$$\begin{aligned} \int S_1 dx \ge \frac{cK_1}{c+d}. \end{aligned}$$

(83)

Notice that $S_1$ satisfies

$$\begin{aligned} {-\nabla ^2 S_1+\frac{1}{D_{S_1}}max ({max~~\beta _1(x),\beta _{10}})S_1>0, \forall x \in {\varOmega }}. \end{aligned}$$

(84)

Thus, together with (83) and Lemma 1 with $q=1$ concludes that

$$\begin{aligned} S_1(x) \ge C, \forall x \in {\varOmega }. \end{aligned}$$

(85)

We next take $min_{\varOmega }I_1(x)=I(x_0)$. According to [34], one can see

$$\begin{aligned} \frac{[\theta \beta _1(x_0)+(1-\theta )\beta _{10}]S_1(x_0)}{S_1(x_0)+I_1(x_0)+c_1} \le \theta a_1(x)+(1-\theta ) a_{10}. \end{aligned}$$

(86)

This leads to

$$\begin{aligned} min_{{\bar{{\varOmega }}}} I_1(x)=I_1(x_0)\ge \frac{min(min_{{\bar{{\varOmega }}}}~~\beta _1(x),\beta _{10}) S_1(x_0)}{max(max_{{\bar{{\varOmega }}}}~~a_1(x),a_{10}))}-c_1-S_1(x_0) \ge C. \end{aligned}$$

(87)

If $\beta _1(x_0) S_1(x_0)/a_1(x_0)-c_1-S_1(x_0)>0$, from the above analysis of step 1 and 2 we can always find a positive constant $C_*>1$, which is independent of $\theta \in [0,1]$, such that any positive solution $(S_1,I_1)$ of (74–75) satisfies

$$\begin{aligned} \frac{1}{C^*}<S_1(x), I_1(x)<C_*, \forall x \in {\bar{{\varOmega }}} \end{aligned}$$

(88)

Step 3. Finally, we find existence of positive solution to (74–75). Let us denote a set,

$$\begin{aligned} {\varTheta }=\{(S_1,I_1) \in C({\bar{{\varOmega }}}) \times C({\bar{{\varOmega }}}): \frac{1}{C^*}<S_1(x), I_1(x)<C_*, \} \end{aligned}$$

(89)

Thus, (74–75) has no positive solution $(S_1,I_1)\in \partial {\varTheta }$. For $\theta \in [0,1],$ we also define the operator

$$\begin{aligned} H(\theta ,(S_1,I_1))= & {} (-\nabla ^2 +I)^{-1}({\hat{h}} (\theta ,(S_1,I_1)),{\tilde{h}}(\theta ,(S_1,I_1))), \end{aligned}$$

(90)

$$\begin{aligned} {\hat{h}}= & {} S+D_{S_1}^{-1}\left( (\theta r_1(x)+(1-\theta r_{10})) S_1 \left( 1-\frac{S_1+I_1}{K_1}\right) \right. \nonumber \\&\qquad \qquad \left. -\frac{(\theta \beta _1(x)-(1-\theta )\beta _1(x)) S_1 I_1}{S_1+I_1+c_1} \right) ,\nonumber \\ {\tilde{h}}= & {} I+D_{I_1}^{-1} \nonumber \\&\left( \frac{(\theta \beta _1(x)-(1-\theta )) S_1 I_1}{S_1+I_1+c_1}-(\theta a_1(x)-(1-\theta ) a_{10}) I_1 \right) , \end{aligned}$$

(91)

Clearly, the existence of positive solutions of (47) is identical to the existence of fixed point of the operator H(1, .) in ${\varTheta }$. From standard elliptic regularity theory one can find that H is a compact operator from $[0,1] \times {\varTheta }$ to $C({{\bar{{\varOmega }}}}) \times C({\bar{{\varOmega }}})$. Furthermore, we have

$$\begin{aligned} (S_1,I_1) \ne H(\theta ,(S_1,I_1)), \forall \theta \in [0,1] ~~ \text {and}~~(S_1,I_1)\in \partial {\varTheta }. \end{aligned}$$

Therefore, the topological degree $deg(I-H(\theta ,.),{\varTheta })$ is well-defined and is independent of $\theta \in [0,1]$. Denote

$$\begin{aligned} S_{10}^*= & {} \frac{a_1(\sqrt{B}+(c_1 r_1+K_1(a_1-\beta _1+r_1)))}{2 \beta _1 r_1} \\ I_{10}^*= & {} \frac{-(a_1-\beta _1)^2 K_1-(\beta _1(c_1-K_1)+a_1(c_1+K_1))r_1+\sqrt{(a_1-\beta _1)^2 B}}{2 \beta _1 r_1} \end{aligned}$$

where $B=(a_1-\beta _1)K_1^2+2K_1(\beta _1(c_1-K_1)+a_1(c_1+K_1))r_1+(c_1+K_1)^2 r_1^2$. [54] have already proved that $(S_{10}^*,I_{10}^*)$ is linearly stable when $\beta _1(x), r_1(x), a_1(x)$ are constant. Using well-known Leray-Schauder degree index formula, we have

$$\begin{aligned} deg(I-H(0,.), {\varTheta })=index(I-H(0,.),(S_{10}^*,I_{10}^*))=1. \end{aligned}$$

Therefore, from homotopy invariance of the Leray-Schauder degree it follows that

$$\begin{aligned} deg(I-H(1,.),{\varTheta })=deg(I-H(0,.), {\varTheta })=1, \end{aligned}$$

which implies that H(1, .) has at least one fixed point in ${\varTheta }$. As a consequence, (47) has at least one positive solution. $\square$