# Pattern Dynamics of Nonlocal Delay SI Epidemic Model with the Growth of the Susceptible following Logistic Mode

Authors: Zun-Guang Guo; Jing Li; Can Li; Juan Liang; Yiwei Yan
Journal: Complexity
 (2020)

---

## Abstract

In this paper, we investigate pattern dynamics of a nonlocal delay SI epidemic model with the growth of susceptible population following logistic mode. Applying the linear stability theory, the condition that the model generates Turing instability at the endemic steady state is analyzed; then, the exact Turing domain is found in the parameter space. Additionally, numerical results show that the time delay has key effect on the spatial distribution of the infected, that is, time delay induces the system to generate stripe patterns with different spatial structures and affects the average density of the infected. The numerical simulation is consistent with the theoretical results, which provides a reference for disease prevention and control.


---

## Body

## 1. Introduction

Research on infectious disease models can be traced back to the pioneering work of Kermack and Mckendrick [1, 2], and the traditional model is usually described by ordinary differential equations as a result of the spatially homogeneous assumption [3–10]. Considering the instantaneous change, researchers have established epidemic models via impulsive differential equations [11–15]. However, besides the change of time, the spatial distribution of populations has an significant effect on diseases transmission. As a matter of fact, spatial distribution of the population is not homogeneous, but depend on the spatial location and the surrounding environment factors, such as the amount of food, the number of natural enemies, and the number of the infected [16, 17]. Therefore, considering the space in epidemic models is more suitable for the spread process of infectious diseases. Recently, many reaction-diffusion models have been proposed to study their dynamical behaviors and study the corresponding control measures [18–29].Since they can reach the current position from any other positions in the whole space when humans or animals are moving in space, which is more reflective for the actual situation, therefore, the nonlocal delay of infectious disease has become a research focus in epidemiology. Wang and Wu [30] studied the Kermack–Mckendrick SIR model with nonlocal delay and discussed the dynamic properties of disease spatial diffusion. Pan [31] studied the existence of wave front solutions for a class of infectious disease models with nonlocal diffusion and time delay by constructing upper and lower solutions. Zhen et al. [32] also obtained the traveling wave solutions for a class of SIR epidemic models with spatiotemporal delays. Using the upper and lower solutions and its related monotone iterative techniques, Tang et al. [33] studied the sufficient conditions for the global asymptotic stability of the disease-free equilibrium of the bird system.Most of the abovementioned studies focus on the traveling wave solutions and global stability of the nonlocal epidemic model. However, few literatures study the pattern formation of infectious disease models based on reaction-diffusion equation with nonlocal delay. Pattern dynamics can effectively characterize the spatial distribution of the infected, so as to provide decision-making guidance for the government. In fact, in the natural environments with limited resources, the assumption of logistic growth is more consistent with the actual situation [34–37]; then, we assume the susceptible possess the logistic growth in the absence of the infected. In addition, we take the standard incidence ratio to describe the infection between the susceptible and the infected. To this end, we propose the following model:(1)∂S∂t=DSΔS+rS1−SK−βS∫−∞t∫ΩQx−y,t−sIy,sdydsS+I−dS−AS,∂I∂t=DIΔI+βS∫−∞t∫ΩQx−y,t−sIy,sdydsS+I−μI−dI,with initial boundary value condition(2)n⟶⋅∇SI=0,Sx,0,Ix,0are given,and the descriptions of the symbols in model (1) are shown in Table 1. The kernel function Qx,t=1/4πtexp−x2/4t1/τexp−t/τ represents the weight from the other possible positions to the position HTML translation failed before time t. The nonlocal item.depicts the cumulative number of the infected reaching x position at t time who starts from any y position in the entire space at s time.Table 1
Symbol descriptions of model (1).
SymbolInterpretationsSx,tThe density of the susceptible at positionx(x∈R2) at time tIx,tThe density of the infected at positionx(x∈R2) at time tAThe emigration rate of the susceptibleΩThe diffusion domainn⟶Unit outward normal to∂Ω∂ΩThe closed boundary of the diffusion domainΩΔLaplace operator∂2/∂X2+∂2/∂Y2DSThe diffusion rate of the susceptibleDIThe diffusion rate of the infectedβThe infection rateμThe mortality rate due to infectious disease causeddThe natural mortalityrThe intrinsic growth rateKThe carrying capacity(3)∫−∞t∫R2Qx−y,t−sIy,sdsdy.In order to study the spatiotemporal dynamics of model (1), in Section 2, we first obtain the steady states of model (1), linearize the model at the endemic steady state, and then analyze the conditions of system generating Turing instability and parameter space of Turing pattern generated. In Section 3, different stripe patterns are shown for different time delays, which indicates that time delay affects the structures of Turing patterns. In Section 4, the summary and discussion are given.


## 2. Linear Analysis and Turing Patterns

In this section, we define a new variableVx,t to replace the nonlocal term ∫−∞t∫R2Qx−y,t−sIy,sdyds; then, system (1) can be transformed into 3-variable reaction-diffusion system:(4)∂c∂t=Fc+DΔc,where(5)c=SIV,Fc=rS1−SK−βSVS+I−d+ASβSVS+I−μ+dI1τI−V,D=DSDI1.By studying system (4) without diffusion, so the ordinary differential equation is given:(6)dStdt=rSt1−StK−βStVtSt+It−d+ASt,dItdt=βStVtSt+It−μ+dIt,dVtdt=1τIt−Vt,system (6) has two constant steady states E0=S0,I0,V0=Kr−A−d/r,0,0 and E1∗=S1∗,I1∗,V1∗, where(7)S1∗=Kr+μ−A−βr,I1∗=Kr+μ−A−ββ−d−μrd+μ,V1∗=Kr+μ−A−ββ−d−μrd+μ.According to the actually biological situation, the population number should be nonnegative. So, ifr>A+d holds, then S0>0, that is, E0 is disease-free steady states. In addition, if the condition(8)H1r>A+β−μ,β>μ+d,is satisfied, then S1∗>0, I1∗>0, and V1∗>0, which shows that E1∗ is the endemic steady state.The main focus of infectious diseases is to study the outbreak and prevalence, so as to provide theoretical guidance for the prevention and control of infectious diseases. Then, the following contents only study endemic steady stateE1∗.Through setting(9)w=S−S1∗I−I1∗V−V1∗,we substitute (9) into (4), expand Fc at the endemic steady state E1∗ based on Talyor expansion, and then obtain the following linear equation:(10)∂w∂t=DΔw+Jw,J=a11a12a13a21a22a23a31a32a33,where(11)a11=Aβ+β2+βd−βr−d2−2dμ−μ2β,a12=β−d−μd+μβ,a13=−d−μ,a21=β−d−μ2β,a22=−2βd+2βμ−d2−2dμ−μ2β,a23=d+μ,a31=0,a32=1τ,a33=−1τ.In order to solve the solution of system (10) satisfying the boundary condition (2), we first define Wx as the time-independent spatial eigenvalue problem, which satisfies the following system:(12)ΔW+k2W=0,n⟶⋅∇W=0,forxon∂Ω,where k is the eigenvalue (i.e., wavenumber).We now look for solutionswx,t of (10) in the form(13)wx,t=C1C2C3eλt+ik⟶·x,where k⟶ is the wave vector with magnitude k=k⟶. By inserting (13) into (10), the characteristic equation matrix is obtained:(14)J−λE−D=0,which is equivalent as follows:(15)λ3+a2kλ2+a1kλ+a0k=0,with(16)a2k=DS+DI+1k2−A−β+d+2μ+r+1τ,a1k=DS+DI+DSDIk4−1βτAβDIτ+β2DIτ−2βdDSτ+βdDIτ−2βDSμτ−βDIrτ+d2DSτ−d2DIτ+2dDSμτ−2dDIμτ+DSμ2τ−DIμ2τ+Aβτ+β2τ−βdτ−2βμτ−βrτ−βDS−βDIk2−1βτ2Aβdτ+2Aβμτ−Ad2τ−2Adμτ−Aμ2τ+3β2dτ+3β2μτ−2βd2τ−6βdμτ−2βdrτ−4βμ2τ−2βμrτ+d2μτ+d2rτ+2dμ2τ+2dμrτ+μ3τ+μ2rτ+Aβ+β2−βμ−βr,a0k=DSDIk6−1βτAβDIτ+β2DIτ−2βdDSτ+βdDIτ−2βDSμτ−βDIrτ+d2DSτ−d2DIτ+2dDSμτ−2dDIμτ+DSμ2τ−DIμ2τ−βDSDIk4−1βτ2Aβdτ+2Aβμτ−Ad2τ−2Adμτ−Aμ2τ+3β2dτ+3β2μτ−2βd2τ−6βdμτ−2βdrτ−4βμ2τ−2βμrτ+d2μτ+d2rτ+2dμ2τ+2dμrτ+μ3τ+μ2rτ+AβDI+β2DI−βdDS+βdDI−βDSμ−βDIr+d2DS−d2DI+2dDSμ−2dDIμ+DSμ2−DIμ2k2−d+μ−d+β−μA+β−μ−rβτ,where a2k>0 for any k because of condition (H1) satisfied.Moreover, one can obtain the corresponding characteristic equation for system (6) as follows:(17)λ3+a20λ2+a10λ+a00=0,where(18)a20=−A−β+d+2μ+r+1τ,a10=−1βτ2Aβdτ+2dμrτ+μ3τ+μ2rτ+Aβ+β2−βμ−βr+2Aβμτ−Ad2τ−2Adμτ−Aμ2τ+3β2dτ+3β2μτ−2βd2τ−6βdμτ−2βdrτ−4βμ2τ−2βμrτ+d2μτ+d2rτ+2dμ2τ,a00=−d+μβ−d−μA+β−μ−rβτ.On the basis of condition (H1), one can derivea20>0 and a00>0. According to Hurwitz criterion, the condition for the local asymptotic stability of E1∗ of system (6) is(19)H2a20a10−a00>0.Theorem 1.
The nonspatial system (6) gives rise to a Hopf bifurcation at endemic steady state E1∗ if only and if the condition(20)H3a20a10−a00=0,is established.Proof.
It is obvious thata20>0 and a00>0 due to the existence of endemic equilibrium E1∗ (i.e., (H1)). According to the theorem of [38], when conditions a20>0, a00>0, and a20a10−a00=0 hold, a Hopf bifurcation occurs for nonspatial system (6).
Next, we find the conditions under which system (4) generates Turing instability near endemic steady state E1∗: E1∗ is locally asymptotic stability for system (6), but E1∗ loses stability for system (6) with diffusion (i.e., system (4)). Furthermore, we already know a2k>0 for any k, while the signs of a0k and a2ka1k−a0k are uncertain. Now, we need to look for the condition that E1∗ of system (4) becomes unstable. Obviously, Hurwitz criterion is not satisfied, that is, a0k<0 or a2ka1k−a0k<0; then, the corresponding instability conditions are given. Based on such two cases, we derive the following theorems with respect to Turing instability.Theorem 2.
For system (4), if condition (H1) holds, and one of the following conditions(C1)h1<0 and H2z1<0(C2)h1>0, h2<0, h22−3h3h1>0 and H2z1<0is satisfied, then endemic steady state E1∗ is unstable for some k, where(21)z1=−h2+h22−3h3h13h3.Proof.
Letz=k2>0 and H2k2=a0k, then H2z=h3z3+h2z2+h1z+h0, where(22)h3=DSDI>0,h2=−1βτAβDIτ+β2DIτ−2βdDSτ+βdDIτ−2βDSμτ−βDSDI−βDIrτ+d2DSτ−d2DIτ+2dDSμτ−2dDIμτ+DSμ2τ−DIμ2τ,h1=−1βτ2Aβdτ+2Aβμτ−Ad2τ−2Adμτ−Aμ2τ+3β2dτ+3β2μτ−2βd2τ−6βdμτ−2βdrτ−4βμ2τ−2βμrτ+d2μτ+d2rτ+2dμ2τ+2dμrτ+μ3τ+μ2rτ+AβDI+β2DI−βdDS+βdDI−βDSμ−βDIr+d2DS−d2DI+2dDSμ−2dDIμ+DSμ2−DIμ2,h0=−d+μ−d+β−μA+β−μ−rβτ>0.
Considering the properties of cubic polynomials, the coefficient of the third order ofH2z is greater than zero and H20=h0>0. By solving the first-order derivative H2z about z, namely, H2′z=3h3z2+2h2z+h1=0, we can deduce(23)z2=−h2−h22−3h3h13h3≤z1=−h2+h22−3h3h13h3.(i)
If conditionh1<0 holds, then one can obtain z2<0<z1, which further shows that z1 is the minimum point. Then, we can get a2k<0 for some k when H2z2,min<0 combining with H20=h0>0 (see Figure 1(a)).(ii)
If conditionsh1>0, h2<0, and h22−3h3h1>0 hold, then we can get 0<z2<z1, and z2 and z1 are maximum and minimum points, respectively. Furthermore, according to H2z1<0, one can derive that a2k<0 for some k (see Figure 1(b)).
Obviously, Hurwitz criterion does not hold for Case (i) or Case (ii), and then the endemic steady stateE1∗ is unstable for system (4).
Thus, the conditions of system (4) generating Turing instability is given:(24)H1,H2,C1or,C2.
In addition, we setH3z=a2ka1k−a0k and z=k2>0, that is,(25)H3z=hq3z3+hq2z2+hq1z+hq0,where(26)hq3=DI+1DS+1DS+DI>0,hq2=−1τβ2AβDSDIτ+AβDI2τ+2β2DSDIτ+β2DI2τ−2βdDS2τ−2βdDSDIτ+βdDI2τ−2βDS2μτ−4βDSDIμτ−2βDSDIrτ−βDI2rτ+d2DS2τ−d2DI2τ+2dDS2μτ−2dDI2μτ+DS2μ2τ−DI2μ2τ+2AβDSτ+2AβDIτ+2β2DSτ+2β2DIτ−2βdDSτ−2βdDIτ−4βDSμτ−2βDSrτ−4βDIμτ−2βDIrτ+Aβτ+β2τ−βdτ−βDS2−2βDSDI−βDI2−2βμτ−βrτ−2βDS−2βDI,hq1=1τ2βA2βDIτ2+2Aβ2DIτ2−4AβdDSτ2−2AβdDIτ2−4AβDSμτ2−4AβDIμτ2−2AβDIrτ2+2Ad2DSτ2+4AdDSμτ2+2ADSμ2τ2+β3DIτ2−5β2dDSτ2−3β2dDIτ2−5β2DSμτ2−5β2DIμτ2−2β2DIrτ2+5βd2DSτ2+14βdDSμτ2+4βdDSrτ2+2βdDIμτ2+2βdDIrτ2+9βDSμ2τ2+4βDSμrτ2+3βDIμ2τ2+4βDIμrτ2+βDIr2τ2−d3DSτ2+d3DIτ2−5d2DSμτ2−2d2DSrτ2+3d2DIμτ2−7dDSμ2τ2−4dDSμrτ2+3dDIμ2τ2−3DSμ3τ2−2DSμ2rτ2+DIμ3τ2+A2βτ2+2Aβ2τ2−2Aβdτ2−4Aβμτ2−2Aβrτ2+β3τ2−2β2dτ2−4β2μτ2−2β2rτ2+βd2τ2+4βdμτ2+2βdrτ2+4βμ2τ2+4βμrτ2+βr2τ2−2AβDSτ−2AβDIτ−2β2DSτ−2β2DIτ+2βdDSτ+βdDIτ+4βDSμτ+2βDSrτ+3βDIμτ+2βDIrτ−2Aβτ−2β2τ+βdτ+3βμτ+2βrτ+βDS+βDI,hq0=1τ2β2A2βdτ2+2A2βμτ2−A2d2τ2−2A2dμτ2−A2μ2τ2+5Aβ2dτ2+5Aβ2μτ2−5Aβd2τ2−14Aβdμτ2−4Aβdrτ2−9Aβμ2τ2−4Aβμrτ2+Ad3τ2+5Ad2μτ2+2Ad2rτ2+7Adμ2τ2+4Adμrτ2+3Aμ3τ2+2Aμ2rτ2+3β3dτ2+3β3μτ2−5β2d2τ2−15β2dμτ2−5β2drτ2−10β2μ2τ2−5β2μrτ2+2βd3τ2+11βd2μτ2+5βd2rτ2+18βdμ2τ2+14βdμrτ2+2βdr2τ2+9βμ3τ2+9βμ2rτ2+2βμr2τ2−d3μτ2−d3rτ2−4d2μ2τ2−5d2μrτ2−d2r2τ2−5dμ3τ2−7dμ2rτ2−2dμr2τ2−2μ4τ2−3μ3rτ2−μ2r2τ2+A2βτ+2Aβ2τ−2Aβdτ−4Aβμτ−2Aβrτ+β3τ−3β2dτ−5β2μτ−2β2rτ+βd2τ+4βdμτ+2βdrτ+4βμ2τ+4βμrτ+βr2τ−Aβ−β2+βμ+βr.
Through some analyses, we can give the following conclusion.Figure 1
The diagrams of the cubic functionH2z in Theorem 1: (a) case (i) and (b) case (ii).
(a)(b)Theorem 3.
For system (4), suppose that the following holds(C3)hq1<0, hq0>0, and H3z3<0(C4)hq1>0, hq0>0, hq2<0, hq22−3hq3hq1>0, and H3z3<0
Then, endemic steady stateE1∗ becomes unstable for some k, where(27)z3=−hq2+hq22−3hq3hq13hq3.Given that the cubic polynomial functionsH2z and H3z have similar structures and conditions, thus the proof process of Theorem 3 is similar to that of Theorem 2; then, one can easily obtain the above result.Therefore, we derive another condition for system (4) to give rise to Turing instability:(28)H1,H2,C3or,C4.


## 3. Numerical Results

In this section, we will carry out numerical simulations for system (4) with Neumann boundary conditions, which will be used to illustrate our theoretical analysis. For the convenience of numerical simulation, we select two-dimensional space area Ω=0,a×0,a and set the same space step h in X and Y directions of the space area. Simultaneously, time step is fixed as Δt. Then, the mesh generation is given:(29)xi=ih,i=0,1,2,…,J,yj=jh,j=0,1,2,…,J,Jh=a,tn=nΔt,n=0,1,2,….Fc in system (4) is relabeled as Fc=f1S,I,Vf2S,I,Vf3S,I,V, and the following finite difference scheme is adopted for system (4):(30)Si,jn+1−Si,jnΔt=DSSi+1,jn−2Si,jn+Si−1,jnh2+Si,j+1n−2Si,jn+Si,j−1nh2+f1Si,jn,Ii,jn,Vi,jn,Ii,jn+1−Ii,jnΔt=DIIi+1,jn−2Ii,jn+Ii−1,jnh2+Ii,j+1n−2Ii,jn+Ii,j−1nh2+f2Si,jn,Ii,jn,Vi,jn,Vi,jn+1−Vi,jnΔt=Vi+1,jn−2Vi,jn+Vi−1,jnh2+Vi,j+1n−2Vi,jn+Vi,j−1nh2+f3Si,jn,Ii,jn,Vi,jn,where Si,jn, Ii,jn, and Vi,jn represent the approximate value of Sxi,yj,tn, Ixi,yj,tn, and Vxi,yj,tn, respectively. In the simulation process, the initial value is a small random disturbance at endemic steady state E1∗, and the area range of individuals’ activity is Ω=0,100×0,100. The space region is evenly divided into 100×100 grids, and the corresponding space step is h=1. The time interval is 0,5000, and the time step is set to Δt=0.005. The first-order difference is used for the first-order derivative of the time, and the second-order derivative of the space is discretized by the center difference method in process of numerical simulation. For infectious diseases, the spatial distribution of infected individuals is a key issue for the control of infectious disease. Turing pattern can visually describe the spatial distribution of infected individuals. The focus of our simulation is how nonlocal delay affects spatial distribution of infected individuals. We thus take parameter value as K=1000, d=0.3, μ=0.2, r=0.8, DS=1, DI=10, A=0.4, and β=0.58; time delay τ is selected as the control parameter. Figures 2 and 3 are the time sequence diagram of spatial distribution of the density of the infected for time delays τ=0.1 and τ=0.6, respectively. It can be seen from these two figures that spatial distribution of the density of the infected evolves with time, which finally presents a regular nonuniform macroscopic structure in limited time, namely, stripe patterns, which indicates the number of infected individuals gradually increases over time and eventually gathers together. Under different time delay parameters, the infected finally forms different spatial distribution structures in the first row of Figure 4, which are called strip patterns. These patterns are further presented in three-dimensional space in the second row of Figure 4, which directly reflects the change of the density of the infected with space. Figure 5 shows the relationship between the average density of the infected and time delay τ, that is, the average density of the infected decreases with the increase of time delay τ (τ is less than the critical value), while the opposite case appears when τ is greater than the critical value.Figure 2
The spatial evolution of the density of the infected with time for time delayτ=0.1: (a) t=90000; (b) t=200000; (c) t=300000; (d) t=1000000.
(a)(b)(c)(d)Figure 3
The spatial evolution of the density of the infected with time for time delayτ=0.6: (a) t=80000; (b) t=180000; (c) t=400000; (d) t=1000000.
(a)(b)(c)(d)Figure 4
The first row shows stripe patterns with three different spatial structures in two-dimensional space for different time delay parameters: (a and d)τ=0.1; (b and e) τ=0.2; (c and f) τ=0.3. These patterns are further presented in a three-dimensional space in the second row.
(a)(b)(c)(d)(e)(f)Figure 5
The relationship between the average density of the infected and time delayτ.


## 4. Conclusions and Discussion

This paper investigates spatiotemporal dynamics of the SI epidemic model with nonlocal effect and the growth of the susceptible population following logistic mode, and we obtain two spatially homogeneous steady states including disease-free steady state and endemic steady state. We mainly analyze the dynamical behavior near endemic steady state-Turing instability. In general, the condition of system generating Turing instability is that the endemic steady state is stable in the system without diffusion (i.e., ordinary equation system), but becomes instability in the system with diffusion (i.e., reaction-diffusion system); then, the system finally gives rise to stable spatial structure, which is called Turing pattern. The conditions of system (4) generating Turing instability are (24) or (28) based on Hurwitz criterion, and the parameter space of the existence of Turing pattern is given by condition (24). While the corresponding parameter space is difficult to obtain because of the complexity condition (28), which will be carried out in the following work. Time delay widely exists in the process of disease transmission, such as incubation period and immune period. We thus study the influence of time delay on infectious diseases in numerical simulations, find that the average density of the infected firstly decreases, and then increase with the increase of time delay, which provides a theoretical support for disease control. At the same time, we simulate the evolution of spatial distribution of the density of the infected under different time delays, and the spatial distribution of the density of the infected finally forms stripe pattern with different spatial structure, which does not change with time. The high incidence area of infectious diseases is found, which also provides data information for prevention and control of infectious disease.In the research process of pattern dynamics of system (4), we consider a special form of kernel function so that the nonlocal delay system could be transformed into general reaction-diffusion system; then, pattern dynamics can be followed. In fact, there are some other forms of kernel functions in the biological system, but different forms shall induce nonlocal delay system to convert into different differential system; thus, this is a very meaningful work to study the effects of different kernel functions on the pattern dynamics in the future.



---

*Source: 9217387-2020-07-26.xml*