Abstract:
This paper investigates the global dynamics for a class of multigroup SIR epidemic model with time fractional-order derivatives and reaction–diffusion. The fractional order considered in this paper is in (0, 1], which the propagation speed of this process is slower than Brownian motion leading to anomalous subdiffusion. Furthermore, the generalized incidence function is considered so that the data itself can flexibly determine the functional form of incidence rates in practice. Firstly, the existence, nonnegativity, and ultimate boundedness of the solution for the proposed system are studied. Moreover, the basic reproduction number R0 is calculated and shown as a threshold: the disease-free equilibrium point of the proposed system is globally asymptotically stable when 
while when 
the proposed system is uniformly persistent, and the endemic equilibrium point is globally asymptotically stable. Finally, the theoretical results are verified by numerical simulation.
Keywords: SIR epidemic model, multigroup, reaction–diffusion, fractional order, asymptotic stability.
Global dynamics for a class of reaction–diffusion multigroup SIR epidemic models with time fractional-order derivatives*
Vilniaus Universitetas
Esta obra está bajo una Licencia Creative Commons Atribución 4.0 Internacional.
Recepción: 03 Noviembre 2020
Publicación: 01 Enero 2022
As we all know, mathematical models play an important role in researching the dynamical behavior of infectious diseases. In the classical epidemic model, it is generally considered that individuals are completely mixed, and everyone has the same possibility of infection. However, due to the differences in age, geographical distribution, and other factors, it is more realistic to divide the total population into several different populations, that is, to establish a multigroup epidemic model. Lajmanovich et al. first proposed the SIS multigroup systems and researched the stability of the endemic equilibrium point [11]. Subsequently, there are many research efforts devoted to investigating the importance of multigroup epidemic models [6, 8, 14]. Guo et al. were the first to successfully establish the complete global dynamics of the multigroup epidemic model based on the basic reproduction number [7]. Boosted by the work of Guo et al., many researchers discussed the stability of various multigroup systems [3, 15, 20, 21, 25].
Meanwhile, individual diffusion behavior is widespread in the actual propagation of infectious diseases. With the development of global transportation, individuals in incu- bation period can easily travel from one place to another, which is thought to be one of the main reasons of the global pandemic of infectious diseases. For instance, SARS first appeared in China’s Guangdong Province in November 2002 and then quickly spread to other parts of China and even the world [26]. Also, COVID-19 was first detected at the end of December 2019 with successive cases occurring worldwide. Therefore, in order to better understand the impact of population mobility on the spread of infectious diseases, it is necessary to incorporate human movement into epidemic model to provide more theoretical guidance for epidemic control. Li et al. analyzed the stability and the uniform persistence of a SIRS epidemic model with diffusion [12]. Xu et al. studied the stability and the existence of traveling wave solutions of a SIS epidemic model with diffusion [28]. Recently, many diffusive epidemic models have been used to model within-group and inter-group interactions in spatially environments, for example, Wu et al. investigated a multigroup epidemic model with nonlocal diffusion and obtained the asymptotic behav- ior of traveling wave fronts [27].
It is worth noting in real life that the spread of infectious diseases not only depends on its current state, but also on its past state. Actually, it can be achieved that current state of fractional-order epidemic models depends on the past information since any fractional derivative contains a kernel function [30]. Furthermore, Smethurst et al. found that the patient waits for the doctor’s time to follow a power law model
[24]. More importantly, Angstmann et al. proposed a infectivity SIR model with fractional- order derivative, and they showed how fractional-order derivative arise naturally by con- tinuous time random walk [2]. As generalized of classical integers ones, Hethcote firstly proposed a fractional-order SIR model with a constant population [8]. Then Almeida et al. considered the local stability of two equilibrium points of a fractional SEIR epidemic model [1].
Typically, the reaction term describes a birth-death reaction occurring in a habitat or reactor. The diffusion term simulates the movement of the individual in the environment in real-world applications. The diffusion is often described by a power law
where 
is the diffusion coefficient, and 
is the elapsed time. In normal diffusion, the order 
But if 
particle undergoes superdiffusion, which mainly describes the process of active cell transport; if 
, this phenomenon is called subdiffusion, which can be the diffusion of proteins within cells or the diffusion of viruses between individuals [29]. And it results in a Caputo time-fractional reaction– diffusion system with fractional order 
. Meanwhile, it is pointed out in [19] that long waiting times model particle sticking, and the density of this process spreads slower than normal diffusion. Also, as shown in [19], Caputo time-fractional reaction– diffusion curve has a sharper peak and heavier tails, which can be used to describe the ability to control the transmission of the disease when only a small number of people are infected, such as COVID-19. The study of subdiffusion system has attracted widespread interest in recent years. Mahmoud et al. studied the Cauchy problem of the fractional- order evolution equation and obtained the expression of the solution of the time fractional- order reaction diffusion system [18]. The subdiffusive predator–prey system is discussed, and the analytical solution of the system is studied in [29]. However, few works have been devoted to studying the subdiffusion epidemic model. Motivated by this, in this work, we focus on time-fractional reaction–diffusion epidemic system, which means the spread of infectious diseases is slower than a Brown motion.
Based on the above discussion, the dynamics of the multigroup SIR epidemic model with generally incidence rates is investigated in this paper. Particularly, the susceptible individuals, infective individuals, and recovered individuals are assumed to follow Fickian diffusion.
The organization of this paper is as follows. A class of diffusive SIR epidemic model with time fractional-order derivatives is formulated and some preliminaries are introduced in Section 2. In Section 3, global dynamics of the proposed model are studied, and numerical simulations are presented to illustrate theoretical results in Section 4. Finally, a brief discussion is given in Section 5.
Before presenting a class of multigroup reaction–diffusion SIR epidemic model with time fractional-order derivatives, some necessary preliminaries are presented.
This section begins with some notations, definitions, and results.
Notation. Let 
be a continuous function; 
be the positive cone of 
with the norm 
where 
and 
be a open set of 
such that 
where 
is the boundary of 
be the positive cone of 
Definition 1. (See [22].) Caputo fractional derivative of order 
for a function 
is defined by
where 
and 
Lemma 1. Let 
be nonnegative and 
be the solution to the following system, respectively:
p2
Proof. Let 
then 
satisfies the following system:
Based on 
and 
we have 
Therefore, it can be deduced that 
Lemma 2. (See [29].) Consider the following system:
(1)
Suppose 
is mixed quasimonotonous and satisfies the local Lipschitz condition
where 
is constant, and 
where 
is a given constant. If the upper solution 
and the lower solutions 
satisfy 
system (1) has a unique solution in 
Lemma 3. The system with time fractional-order derivatives
(2)
has a unique global asymptotic stability of constant equilibrium 
Proof. Define the Lyapunov function
Calculating the fractional derivative of 
along the trajectories of system (2), one has
and .. = 0 if and only if . = ... Then according to [4], there exists a unique global asymptotic stability of constant equilibrium .. = b/µ for system (2).
Lemma 4. Consider the following system:
(3)
where 
satisfies the local Lipschitz condition, and 
Then for 
one has 
with 
Proof. According to the definition of Caputo fractional-order derivative, one has
Let 
then
By calculation,
Hence, 
is a solution of system (9). Further, 
By the uniqueness of the solution it is deduced that 
In [10], Korobrinikov et al. studied a multigroup SIR model as follows:
(4)But individual movement is not be considered in system (4) that is unrealistic, then Wu et al. considered the following SIR epidemic model with diffusion [27]:
(5)Based on the previous analysis, since fractional order has the long-term memory, which can describe the spread of infectious diseases more accurate. In addition, it is traditionally assumed that the incidence of disease transmission is bilinear with respect to the number of susceptible individuals and the number of infected individuals. But in reality, it is often difficult to obtain detailed information on the spread of infectious diseases because they may change with the surrounding environment. Therefore, the general incidence rates will be chose in this paper. Motivated by the above work, as an extension of system (5), a class of multigroup SIR epidemic model are investigated as follows:
(6) Where 
implies Caputo fractional-order operator 
denotes the Laplace operator;
denotes the outward normal derivative on the smooth bound- ary 
and 
represent the number of the susceptible, infective, and recovered individuals in 
group at time 
and spatial location 
respectively; 
imply the nature death rates of 
and 
in 
, respectively; 
denotes the disease-related death rates of 
in 
; 
represents the recruitment rate of the total population; 
implies the recovery rate of the infected individuals in kth group; 
denotes the diffusion rate of 
and 
in 
group;
represents the infection rate of 
infected by 
Furthermore,
and 
are positive constants for 
and 
are nonnegative constants for 
Before giving the main results, hypothesis in terms of generalized incidence rates
and
is made as follows:
and 
satisfy the local Lipschitz condition and 
is strictly monotone increasing on
and
is strictly monotone increasing on
for all 
for all 
where 
is nonnegative and irreducible. Furthermore,
if and only if 
for
and
if and only if 
for all 
Remark 1. Note that under hypothesis (H), many existing models can be regarded as a special form of system (6), such as 
and other nonlinear incidence rate in [16].
Some dynamical behavior of system (6) are investigated in this section. Here it can be found that the susceptible class 
and the infected class 
are not effected by the recovered class
of system (6). Hence, we will focus our attention on the following reduced system:
(7)Then some basic properties of system (7) are discussed in following parts.
It is significant to demonstrate the existence, uniqueness, and boundedness of a nonnega- tive solutions for system (7) before implementing its stable process. Thus, this subsection moves to the discussion of proprieties mentioned above.
Theorem 1. Under hypothesis (H), there exists a unique nonnegative solution
of system (7), and it is also ultimately bounded for any given initial function 
where 
and 
Proof. Consider these two function 
and
According to condition (i) in hypothesis (H), it is obvious that 
and 
are mixed quasimonotonous. Consider the following auxiliary system:
It is obvious that
is a pair of the lower solution to system (7). Then, according to Lemma 1, one has 
Furthermore, the following auxiliary system is introduced:
(8)then the above system (8) has a solution as follows:
Therefore,
then there exists a constant 
satisfied 
Further, consider the following auxiliary system:
(9)then the solution for the above system (9) is
Similarly, 
Since
then
(10)However, it is easy to see that
(11)It can be deduced from Eqs. (10) and (11) that
Similar to the above analysis, it can be obtained the following equation:
Based on the above analysis and Lemma 2, system (7) has a unique nonnegative global solution. Furthermore, the expression for the solution of system (7) is
where
with
represents a probability density;
represent generated strong continuous operator semigroups by 
denoting generated strong continuous operator semigroups by
can be rewritten by [23]
where
is the Green function yielded
with 
be the eigenvalue of
with the eigenfunction 
satisfying
Hence, by the boundedness of the eigenfunction
one has
According the upper solution
of system (7), one has 
which implies 
is ultimate bounded. Further, the ultimate bounded of 
will be analyzed. Let 
and 
Adding the first two equations of system (7) and integrating it on 
one has
Therefore, by [13] one has
then there exist two constants
and
satisfying
According to [17], the operator families 
is uniformly bounded. Hence, there exist two constants
and 
satisfying
for 
Finally, the uniformly boundedness of the infected group 
can be studied as follows:
where 
thus 
is ultimate bounded. Therefore, there exists a unique positive global solution
of system (7), and it is also ultimately bounded.
In this subsection, the global stability analysis of system (7) will be discussed. It is easy to find that the disease-free equilibrium point 
of system (7) always exists where 
Define the following fuction:
where 
and 
is the spectral radiuses of the matrix 
Lemma 5.The basic reproduction number 
Proof. Linearizing system (7) at the disease-free equilibrium point 
, one has
Let 
and 
Then it is easy to find 
Obviously, we have
By the definition of the basic reproduction number [5] one has
Thus, by the properties of matrix eigenvalues it can be deduced that 
Therefore,
is considered as a threshold parameter in place of 
In the following, the uniqueness and the global stability of 
are studied.
Theorem 2. Under hypothesis (H) and 
there exists the unique equilibrium point 
of system (7), and it is globally asympototically stable in domain 
where
Proof. Let 
and 
where 
Define
It is clear to find from 
that 
Then one has 
Since 
is irreducible, it can be obtained that 
and 
are irreducible. 
is also irreducible. If 
, the inequality 
holds. Further, it can be deduced that
if 
and 
Thus,
has a only trivial solution 
This shows that 
is the unique equilibrium of system (7) when 
Further, 
is positive, then 
is an eigenvalue of the matrix 
, and 
has a nonnegative eigenvector corresponding to 
Let 
be the positive left eigenvector of 
corresponding to the spectral radius 
that is,
Define Lyapunov function
Calculating the time fractional derivative of 
along the trajectories of system (7), one has
Let ..(.)
Let 
which is a positive definition function in 
Then it is concluded from [4] that 
is globally asymptotically stable in domain 
Theorem 3. Under hypothesis (H) and 
system (7) is uniform persistence, that is, for any initial value 
the solution 
satisfies
where
is a constant.
Proof. Define a set
and
Let
be the solution of system (7) under the initial value 
For any 
it can be known that all nonnegative solutions
generate a solution semiflow 
with
Thus, we have 
and it is obvious that
where 
is the identity matrix. It can be deduced from Lemma 4 that
Then 
Based on the above analysis, one has 
-semigroup on 
Obviously,
is compact for 
and point dissipative in 
The following system is considered:
It can be found from Lemma 3 that 
is globally asymptotically stable. Thus, system (7) is globally asymptotically stable at the disease-free equilibrium point
It can be deduced that the disease-free equilibrium 
in is 
a global attractor of 
which implies 
Let 
where 
Considering 
there exists a sufficiently small constant 
such that 
where
If 
there exists a solution 
of system (7) with the initial value such that 
then there exists a constant 
such that 
and 
Since 
is irreducible,
is irreducible. Let
be the positive left eigenvector of
corresponding to the spectral radius 
that is,
Define the following arbitrary function:
Calculating the time fractional derivative of 
along the trajectories of system (7), one has
which leads to a contradiction with 
Therefore,
Thus, it can be deduced from [25] that
is uniformly persistent. It is concluded that system (7) is uniformly persistent.
The ultimate boundedness and the uniform persistence imply the existence of a pos- itive equilibrium point of system (7). Therefore, the existence and global stability of the positive endemic equilibrium point of system (7) can be further discussed.
Theorem 4. Under hypothesis (H) and 
system (7) has at least one endemic equilibrium 
satisfying
Furthermore, if
system (7) is globally asymptotically stable at the endemic equilibrium point 
Proof. According to Theorem 1, for any given initial condition
the corresponding solution
is ultimately bounded, and system (7) is uniformly persistent when
Therefore, there exists a positive equilibrium point 
of system (7) that satisfies Eqs. (12), (13).
Next, the global stability of 
will be analyzed. Define the Lyapunov function
Where
and the coefficients
will be determined in Eq. (20). Calculating the time fractional- order derivative of 
along the trajectories of system (7), it can be conclude that
(15)For each 
it can be deduced from the divergence theorem that
Thus, Eq. (15) can be deduced that
(16)Let 
then
(17)Substituting Eq. (17) into Eq. (16), the following inequality holds:
Calculating the fractional-order derivative of 
along any solution of system (7), one has
(18)Since 
is the endemic equilibrium point of system (7), one has
(19)where 
denotes the cofactor of the 
diagonal entry of e 
, where
(20)with
It can be deduced from [7] that
exists a unique positive solution
Therefore,
(21)Thus, substituting Eq. (19) into Eq. (21), it can be obtained that
(22)Further, it is concluded from Eqs. (18) and (22) that
Based on [4], the endemic equilibrium point 
of system (7) is globally asymptotically stable.
Corollary 1.When 
the endemic equilibrium point 
is globally asymptotically stable if hypothesis (H) and 
satisfied.
Remark 2. It can be seen that Corollary 1 is similar with Theorem 6 of [17] when 
Remark 3. Not considering infection between populations, that is, when 
the reproduction number 
of group k is 
Furthermore,the disease-free equilibrium point 
is globally asymptotically stable when 
, and the endemic equilibrium point 
is globally asymptotically stable when 
In order to verify theoretical results numerically, numerical simulations are presented in this section. We consider system (7) with two-group case, which is suitable for infectious diseases transmitted between two cities or communities. Furthermore, system (7) with two groups (n = 2) can be calculated by the central difference method in L1-type space and Alikhanov-type discretization in time [9]. Furthermore, we consider the following incidence rate as an example:
which τ is a positive parameter measuring the psychological or inhibitory effect. Obviously,
satisfy hypothesis (H). The corresponding system can be expressed as
(23)Let assign the following values to the parameters of system (23):
It is easy to calculate that 
Based on by Theorem 2, the disease-free equilibrium point 
of system (23) is global stable which is verified by Figs. 1 and 2. Further, the following parameters is chose:
System (23) has a unique equilibrium point 
It can be calculated that 
and Eq. (14) is satisfied. Based on the above analysis, the endemic equilibrium point 
of system (23) is global stable, which is verified by Figs. 3 and 4.
Further, with regard to the disease-free equilibrium point of the first group, the influ- ence of different fractional order . on the stability of the infected are discussed. The error
Figure 1
The first group stability of the disease-free equilibrium E0.
Figure 2
The second group stability of the disease-free equilibrium E0.
Figure 3
The first group stability of the endemic equilibrium E∗.
images of the infected of 
are described in Fig. 5, respectively. It is easy to seen from Fig. 5 that although the infected will disappear, different order 
will have a sensitive effect on the change of solution. Further, when 
tends to 1, the numerical solutions of system (7) are also convergent to the solutions of the classical ones [17]. But the relationship between the change of the solution for system (7), and fractional order 
is not discussed.
Figure 4
The second group stability of the endemic equilibrium E∗.
Figure 5
The error about different fractional order of the disease-free equilibrium point.
In this article, incorporating the population diffusion and time fractional-order derivatives, theory analysis of a class of multigroup SIR epidemic model are investigated. Firstly, the existence and uniqueness of the nonnegative solution for system (7) are established. By using Lyapunov functions the global stability of the disease-free equilibrium point E0 is obtained when the basic reproduction number 
Besides, when 
the uniform persistence and the global stability of the endemic equilibrium point E∗ are discussed. The proposed model, a more accurate epidemic model, can help us to understand some dynamical behaviors of infectious diseases. Moreover, theoretical results may provide some useful guidance for making effective countermeasures on infectious diseases. However, the relationship between system (7) and fractional order α is still an open question, which will be our future work.
Figure 1
The first group stability of the disease-free equilibrium E0.
Figure 2
The second group stability of the disease-free equilibrium E0.
Figure 3
The first group stability of the endemic equilibrium E∗.
Figure 4
The second group stability of the endemic equilibrium E∗.
Figure 5
The error about different fractional order of the disease-free equilibrium point.
