Global behavior of a multi-group SEIR epidemic model with age structure and spatial diffusion

Abstract: Different epidemic models with one or two characteristics of multi-group, age structure and spatial diffusion have been proposed, but few models take all three into consideration. In this paper, a novel multi-group SEIR epidemic model with both age structure and spatial diffusion is constructed for the first time ever to study the transmission dynamics of infectious diseases. We first analytically study the positivity, boundedness, existence and uniqueness of solution and the existence of compact global attractor of the associated solution semiflow. Based on some assumptions for parameters, we then show that the disease-free steady state is globally asymptotically stable by utilizing appropriate Lyapunov functionals and the LaSalle’s invariance principle. By means of Perron-Frobenius theorem and graph-theoretical results, the existence and global stability of endemic steady state are ensured under appropriate conditions. Finally, feasibility of main theoretical results is showed with the aid of numerical examples for model with two groups which is important from the viewpoint of applications.


Introduction
Since the pioneering work of Kermack and McKendrick [1], many mathematical models have been proposed attempting to gain a better understanding of disease transmission, especially for the control strategy and dynamical behavior of infectious diseases [2][3][4][5][6][7][8]. Simple models with assumption that individuals are well mixed, which implies each individual has the same probability to be infected, are beneficial in that one can obtain analytical results easily but may be lack of realism. Epidemic models with population structures, like age, sex and patch (such as communities, cities, or counties), may be a more realistic way to describe complex disease dynamics. As a matter of fact, the total population should be classified into different groups and the vital epidemic parameters should vary among different population groups. In addition, at different age stages, the effects of infectious transmission are various, which is another important and key factor that needs necessarily to be included in model-ing this infectious transmission process. Thus, considering multi-group and age structure in epidemic models is very necessary and reasonable. Some recent developments on the transmission dynamics of multi-group and age structured epidemic models have been discussed in [8][9][10][11][12][13].
Since that the population distribute heterogeneously in different spatial location in the real life and they will move or diffuse for many reasons, in epidemiology, there is increasing evidence that environmental heterogeneity and individual motility have significant impact on the spread of infectious diseases [14,15]. In recent years, global behavior of spatial diffusion systems, which are suitable for diseases such as the rabies and the Black Death, has been attracted extensive attention of researchers and has been one of the hot topics [16][17][18][19][20][21][22][23][24][25]. Among these works, few take age structure or multi-group into consideration. Yang et al. [24] proposed a novel model incorporated with both age-since-infection and spacial diffusion of brucellosis infection, and the basic reproduction number and global behaviors of this system were completely investigated. Fitzgibbon et al. [19] considered a diffusive epidemic model with age structure where the disease spreads between vector and host populations. Then, the existence of solutions of the model was studied based on semigroup theory and the asymptotic behavior of the solution was analyzed. Luo et al. [21] incorporated spatial heterogeneity in n-group reactiondiffusion SIR model with nonlinear incidence rate to investigate the global dynamics of the disease-free and endemic steady states for this model. Zhao et al. [25] modeled host heterogeneity by introducing multi-group structure in a time delay SIR epidemic model and showed that basic reproduction number determines the existence of traveling waves of this system. To determine how age structure, multigroup population and diffusion of individuals affect the consequences of epidemiological processes, Ducrot et al. [18] formulated a multi-group age-structured epidemic model with the classical Fickian diffusion and studied the existence of travelling wave solutions for this model.
To the best of our knowledge, epidemic models established by researchers except for [18] only include one or two characteristics of multi-group, age structure and spatial diffusion. All the three characteristics are incorporated into epidemic model in [18], however, this model does not include the class of latent individuals. For some epidemic diseases like malaria, HIV/AIDS and West Nile virus, latent individuals may take days, months, or even years to become infectious. Moreover, the travel of latent individuals showing no symptoms can spread the disease geographically which makes disease harder to control. Motivated by the above discussion, in this paper, we investigate a diffusive version of multi-group epidemic system with age structure which is generalization of the model studied in [26] for the first time to allow for individuals moving around on the spatial habitat x ∈ Ω ⊂ R n with smooth boundary ∂Ω. The organization of this paper is as follows. Firstly, we present our model in the next section. In section 3, some preliminaries including the positivity, boundedness, existence and uniqueness of solution, and the existence of compact global attractor of the associated solution semiflow, are established. In section 4, the sufficient conditions on the existence and global stability of disease-free and endemic steady states are stated and proved. In section 5, we conduct numerical simulations to illustrate the validity of our theoretical results. In section 6, a brief conclusion is given.

The model
In 2015, Liu et al. [26] introduced age-of-latent and age-of-relapse into epidemic model which is appropriate for diseases such as tuberculosis and herpes virus infection. For these diseases, latent individuals may take days, months, or even years to become infectious and the treatment efficacy may decline with time for recovered individuals then cause recurrence of disease. In order to study the global dynamics o for these diseases, they formulated the following SEIR epidemic system with continuous age-dependent latency and relapse for t ≥ 0 and with initial conditions is the space of functions on (0, +∞) that are nonnegative and Lebesgue integrable. At time t, the densities of susceptible individuals, latent individuals with latent age a, infectious individuals and removed individuals with relapse age b are denoted by S (t), e(t, a), I(t), r(t, b), respectively. σ(a) and γ(b) denote the conversional rate from the latent class and the relapse rate in the removed class, which depend on age a and age b, respectively. Furthermore, β is the transmission rate of the disease between susceptible and infectious individuals, Λ is the density of the recruitment into the susceptible class (including the births and immigration), µ is the natural death rate of all individuals, δ 1 and δ 2 are the additional death rate induced by the infectious diseases, and c is the recovery rate from the infectious class. All parameters are assumed to be positive.
It is clear that the variations of different epidemic parameters between or within different groups can be well realized according to the description of multi-group epidemic models. Hence, Liu and Feng [27] extended model (2.1) to the situation in which the population is divided into n groups according to different contact patterns and derived the following multi-group SEIR epidemic model e k (t, 0) = n j=1 β k j S k (t)I j (t), r k (t, 0) = c k I k (t), for t ≥ 0 and with initial conditions for a, b > 0. Λ k , µ k and c k denote the recruitment rate of the susceptible class, the per-capita natural death rate and the recovery rate from the infectious class in group k, respectively. β k j denotes the transmission rate of the disease between susceptible individuals in group k and infectious individuals in group j. δ 1k and δ 2k denote the additional death rates of exposed and infectious individuals induced by the infectious diseases in group k, respectively. σ k (a) denotes the conversional rate from the latent class in group k, which depends on age a and γ k (b) denotes the relapse rate from the removed class into the infectious class in group k, which depends on age b. Spatial diffusion is an intrinsic characteristic for investigating the roles of spatial heterogeneity on diseases mechanisms and transmission routes and can lead to rich dynamics. Based on this fact, we generalize (2.2) by taking account of the case that individuals move or diffuse around on the spatial habitat x ∈ Ω ⊂ R n with smooth boundary ∂Ω. Let S k (t, x) and I k (t, x) be the densities of susceptible individuals and infectious individuals at time t and location x ∈ Ω in group k, respectively, where the habitat Ω is bounded and connected. And let e k (t, a, x) and r k (t, b, x) denote the densities of individuals in the latent class with age a and the removed class with age b at time t and location x in group k, respectively. Hence, the n-group diffusive SEIR epidemic model with age-dependent latent and relapse has the following form for x ∈ Ω, a, b ∈ R + = (0, +∞), with the homogeneous Neumann boundary conditions and initial functions d 1k , d 2k , d 3k , d 4k denote the diffusion coefficients of susceptible individuals, exposed individuals, infectious individuals and removed individuals in group k, respectively. And the other parameters have the same biological meanings as in (2.2). The homogeneous Neumann boundary conditions imply that there is no population flux across the boundary ∂Ω.
We define the functional spaces X = C(Ω, R) and Y = L 1 (R + , X) for model (2.3) equipped, respectively, with the norms The positive cones are denoted by X + and Y + . In addition, we define a vector space Throughout this paper, for convenience, we always denote S = (S 1 , S 2 , ..., S n ), e = (e 1 , e 2 , ..., e n ), for all t > 0 and φ ∈ C(Ω, R), where Γ ik (t, x, y) is the Green function. We have that T ik , i = 1, 2, 3, 4, k = 1, 2, ..., n are compact and strongly positive for each t > 0 by the Corollary 7.2.3 in [29]. Integrating the second equation in model (2.3) along the characteristic line t − a = c, where c is a constant, we obtain To study the asymptotic behaviors of the dynamics of model (2.3), we require the following assumptions on the model parameters.

Preliminaries
Proof. From the first equation of (2.3), we have The positivity of C which means the positivity for A k and B k , k = 1, 2, ..., n is established by constructing Picard sequences as follows.
Solving equation I k for system (2.3), we have 3) and the definitions of A k and B k , we have and Then it is obvious that From the positivity of β k j , σ k and γ k , together with the positivity of Γ 2k , Γ 3k and Γ 4k , it follows that Hence applying mathematical induction, we show that the sequence {C (m) } ∞ 0 is monotonically increasing.
Next, applying the contraction mapping principle, we show the sequence {C (m) } ∞ 0 converges to C(t, x) for any (t, x) ∈ [0, T ] × Ω as m approaches infinity. To this end, we define a variablê

By the definitions of A (m)
k and B (m) k , we havê |S k (t, ·)| X , and Hence, . Therefore, for any We choose λ sufficiently large such that n j=1 β k jŜ k (σ j +γ j ) λ 2 < 1 and c k (σ k +γ k ) k , r 0 k , A k and B k , we conclude that e k (t, a, x) and r k (t, a, x) are positive. For the positivity of I k , we prove this by contradiction. Suppose that there exist x 0 ∈ Ω and t 0 = inf{t ∈ R + |I k (t, x 0 ) = 0} such that By the third equation of system (2.3), we can easily obtain This leads to a contradiction. Hence, for any t ∈ [0, T ], we have (S (t, ·), e(t, ·, ·), I(t, ·), r(t, ·, ·)) ∈ X n + × Y n + × X n + × Y n + .
Let N k (t) = Ω S k (t, x) + ∞ 0 e k (t, a, x)da + I k (t, x) + ∞ 0 r k (t, b, x)db dx denotes the total population at time t in group k and region Ω.
for all t > 0, x ∈ Ω, the region Π defined by is positively invariant for system (2.3).
Proof. Following condition (3.4) and the equations of system (2.3), we have Noting the Neumann boundary conditions of system (2.3) and using the Gauss formula, we derive It follows that Thus if N k (t) > Λ k µ k |Ω|, then dN k (t) dt < 0. Moreover, we observe the ordinary differential equation where N k (0) means the initial value of total population in group k and region Ω. By applying the standard comparison theorem, we have for all t ≥ 0, Hence, Π is positive invariant for system (2.3).
The existence and uniqueness of the solution of model (2.3) follow from Banach-Picard fixed point theorem.
Lemma 3.1. [31]. The following problem admits a unique positive steady state ω * = Λ µ , which is globally attractive in X. Theorem 3.4. Let initial functions satisfy (S 0 , e 0 , I 0 , r 0 ) ∈ X n + × Y n + × X n + × Y n + . Then the system (2.3) has a unique solution (S (t, ·), e(t, ·, ·), I(t, ·), r(t, ·, ·)) ∈ X n + × Y n + × X n + × Y n + for t ∈ R + . Proof. To extend the domain of existence from t ∈ [0, T ] to t ∈ R + , it suffices to show that the solution does not blow up in finite time. In fact, by Theorem 3.1, we have that for all t > 0 and x ∈ Ω. From Lemma 3.1 and the comparison principle, we get that S k (t, x) is bounded above by the upper solution Λ k µ k . We now claim that e k (t, a, x) < +∞ for all t > 0, a > 0, x ∈ Ω and k = 1, 2, ..., n. From (2.4), it is sufficient to show that e k (t, 0, x) < +∞ for all t > 0 and x ∈ Ω. Suppose on the contrary that there exist t e > 0 and x e ∈ Ω such that lim t→t e −0 e k (t, 0, x e ) = +∞.
Proof. We construct a Lyapunov function Taking the derivation of V 1k (t, x) along the trajectory of (2.3) with respect to t, we have Recalling (2.4), we can rewrite V 2k (t, x) as follows Thus, we calculate ∂V 2k ∂t along the solution of system (2.3) and get ∂V 2k × e k (t, a, x)da. (4.8) Similarly, Thus, Hence, combining (4.7)-(4.10), we calculate the derivative of V k (t) along the solution trajectory of (2.3) as Using (4.6), we yield that Thus, from (4.5), we have χ k (0) < 1, ψ k (0) < α k , which implies the global asymptotic stability of disease-free steady stateĒ. = 0 for all t > 0, x ∈ Ω, then the space-independent steady state E * is globally asymptotically stable.
Proof. We define g(x) = x − 1 − ln x, clearly, g(x) is always positive for x > 0 and g (x) = 1 − 1 x . Consider a matrix D = (β k j ) n×n with entryβ k j = β k j L k S * k I * j and a digraph G = (U, H) which contains a set U = {1, 2, ..., n} of vertices and a set H of arcs (k, j) leading from initial vertex k to terminal vertex j, then, we denote a weighted digraph as (G, D) for which each arc ( j, k) is assigned a positive weight β k j . Furthermore, we denoteD as the Laplacian matrix of matrix (G, D). Then, the irreducibility of matrix (β k j ) n×n implies thatD is also irreducible. Let q k denote the cofactor of k-th diagonal element ofD. And we construct a Lyapunov function as the following form For convenience, we denote J k (t, x) = n j=1 β k j I j (t, x) and J * k = n j=1 β k j I * j . Taking the derivation of W 1k along the trajectory of (2.3) with respect to t, we have Calculating the derivative of W 2k along the solution of system (2.3) yields Similarly, Then, we take the initial functions as then we have ∞ 0 α 1 σ 1 (θ)π 11 (θ)dθ ≈ 0.00051 < 1, ∞ 0 γ 1 (θ)π 21 (θ)dθ ≈ 0.0817 < 1, ∞ 0 α 2 σ 2 (θ)π 12 (θ)dθ ≈ 0.00027 < 1, ∞ 0 γ 2 (θ)π 22 (θ)dθ ≈ 0.13319 < 1. In this case, from Theorem 4.2, we expect the disease-free steady stateĒ to be globally asymptotically stable. In fact, in Figure 1, the infective population I 1 (t, x) and I 2 (t, x) converge to zero over time. then we have ρ(M 0 ) ≈ 3.89217 > 1. In this case, from Theorem 4.3, we expect the space-independent steady state E * to be globally asymptotically stable. In fact, in Figure 2, the infective population I 1 (t, x) and I 2 (t, x) converge to the positive distribution over time.  Figure 2. Time evolution of infective population I 1 (t, x) and I 2 (t, x) for system (2.3) with parameters (5.2).

Conclusions
In this paper, as an additional structure of the system, we focus on the spatial diffusion of the population. Models with spatial diffusion allow individuals to move to adjacent positions through a random walk process, this is a key factor in considering the geographical spread of infectious diseases. Firstly, we propose the n-group diffusive SEIR epidemic model with age-dependent latent and relapse, it is a generalization of the model in [27] to a spatially diffusive system. Then, we investigate the positivity, boundedness, existence and uniqueness of solution and the existence of compact global attractor of the associated solution semiflow for this system. For these results, we use the method of constructing Picard sequences, Banach-Picard fixed point theorem and theories of partial functional differential equations. Thereafter, we establish the existence of disease-free and endemic steady states based on Perron-Frobenius theorem. we utilize appropriate Lyapunov functionals, graph-theoretical results and the LaSalle's invariance principle to prove the global stability of disease-free and endemic steady states. Thus, we presented the results of numerical simulations to verify the validity of our main theorems. This is important from the viewpoint of applications.
In this epidemic model, we are concerned with two kinds of spatial heterogeneity: the patch structure and spatial diffusion. Furthermore, age-of-latent and age-of-relapse are included into the epidemic model which is appropriate for diseases such as tuberculosis and herpes virus infection. Dynamical results obtained in this paper provide theoretical foundation for seeking effective measures to prevent, control and study disease transmission.
The expressions of basic reproduction number and endemic steady state depends on space are not analyzed in this paper owing to the complexity of model. In addition, how to improve the sufficient conditions that ensure the stabilities of steady states and make them be depended on basic reproduction number is also need to investigate. We leave these issues for future research.