The Global Dynamics of an Age-Structured Hand-Foot-Mouth Disease Model with Saturation Incidence and Time Delay

,


Introduction
Hand-foot-mouth disease (HFMD) is a viral infection primarily caused by enteroviruses.Coxsackie Asckievirus (Cox A16) and enterovirus (EV71) are the most prevalent HFMD virus strains [1].HFMD mostly afects children younger than the age of ten.But older children and adults may also be infected with HFMD [2].One of the most substantial components of death in children is HFMD, which has become a main public health worry that may have unquantifable societal and economic consequences.Tus, researches concerns about the spread of HFMD have risen in recent years.
Because there is no anti-HFMD medical therapy or vaccination available, a huge number of scholars have built mathematical models to study the transmission dynamics and illness control of HFMD.Mathematical models of HFMD infection have advanced signifcantly over the last few decades.A simple SIR (susceptible-infected-recovered) model was developed to forecast the number of infectious individuals and the length of the outbreak in [3].After that, more real factors have been included into this simple SIR model, and more complicated mathematical models were proposed in studying the dynamics of HFMD.For instance, Roy and Halder [4] established an SEIR model in 2010 by dividing the infectious group into asymptomatic (E (t)) and symptomatic (I (t)) infectious individuals.Furthermore, Liu [5] formulated an SEIQR model to investigate the efect of quarantine Q (t) on the disease development.Considering the asymptomatic individuals A and contaminated environments W, Wang et al. [6] discussed a SEIARW system and found that both factors are critical to delay and avoid epidemic outbreaks.Tan and Cao [7] developed a SEIVT epidemic model with vaccination to investigate HFMD's transmission in 2018.Recently, Sun et al. in [8] proposed a difusive HFMD model with a fxed latent period and nonlocal infections.Besides, these mathematical models are also used in the research of other infectious diseases, such as COVID-19 [9][10][11]; cofee berry disease [12]; and banana black Sigatoka disease [13].
In the transmission of HFMD, individuals who are latently infected have no symptoms and cannot convey the virus to others until several days after exposure (usually 3-7 days [23]).As a result, time delay (latent period) is an important factor in the research of the dynamic behaviors of HFMD systems.Many scholars have investigated delayed HFMD models, for instance [14,24,25].Te active HFMD individuals will disseminate the bacterium in the population once the latently infectious persons become active.Tus, infected people's infectivity is determined by their infection age, which is a signifcant factor for infectious diseases, particularly the ones with a prolonged latency and fuctuating infection rate like HFMD.
Note that as time progresses, the disease develops within the individuals with the diferent infectivity or many times one may easily observe vary infectiousness of an infected individual at diferent stages of infection.Particularly in the transmission of HFMD, latent individuals have no symptoms and cannot convey the virus to others.When asymptomatic individuals can infect others (that is to say they have infectiousness), depending on how long they have been infected, which is usually called infection age [26].Moreover, the infectivity of the host might continuously change with time and infection age.Tis means that infection age may be one of the informative factors to model for the diseases like HFMD.Hence, age-structured epidemic models should be extensively examined in order to provide better knowledge and further insights into transmission mechanisms.Te authors of [26][27][28][29][30][31][32][33][34] have widely explored the dynamic behavior of age-structured epidemic systems.
Inspired by the preceding discussion, in this work, we will investigate the HFMD model with age structure and time delay.Let E(a) be the infected individuals' transmission rate at age a.To establish the fundamental form of E(a), we utilize the data (see Table 1) from the China Center for Disease Control (CDC) [35] in 2018 and the least square approach to match the transmission rate E(a) at varied infection age a. Te numerical simulation shows that the transmission rate E(a) of HFMD follows the exponential function 6.275e − 0.3873a with a ≥ τ (see Figure 1), after the shortest period τ necessary during the initial infection to become infectious class.Terefore, when we discuss HFMD, the transmission rate E(a) is taken as the following exponential function: where τ represents the shortest time required from initial infection to become infectious class, E(a) ∈ L ∞ + ((0, +∞), R), κ and α * are constant coefcients of (1), κ > 0 and α * > 0. Te equation (1) means that if the infection age a of the infected individuals is less that the time delay τ, the infected individuals do not have the infectivity; otherwise, the transmission rate is α * e − κa .
Te form of the age-structured HFMD model with saturation incidence and time delay is investigated as follows: where δ(a): � μ + E(a), and with the initial conditions where S(t), I(t), Q(t), and R(t) represent the number of susceptible persons, infectious persons, quarantined persons, and recovered persons at time t, respectively.Te number of latent persons with age a is expressed by (t, a).As a result, system (2) is an SEIQR-type HFMD system, and the infection process fowchart is shown in Figure 2.  2) are considered biologically relevant, and all of them are positive.Te biological interpretations of the parameters of the system are described in Table 2.
It is worth pointing out here that, in comparison to the model in [14], system (2), which considers the infection age and infectious delay, is obviously more general and realistic.Te asymptotic behavior of solutions to system (2) is the focus of this work.Firstly, we will determine the basic reproduction number R 0 , and then study the infection-free equilibrium E 0 's asymptotic stability for system (2) when R 0 < 1 and R 0 > 1, respectively, through discussing the position of the corresponding characteristic equations' eigenvalues.Meanwhile, when R 0 > 1, the global stability of a positive equilibrium E * is also discussed by constructing a suitable Lyapunov function.Finally, we are particularly interested in how age-structure, time delay, and saturation incidence afect the dynamics of the system under consideration.
Te following is the remainder of this work.In Section 2, we give some preliminary conclusions and prove the wellposedness of the system (2).In Section 3, the existence of the equilibria, particularly the existence and uniqueness of the endemic equilibrium, as well as the linearized system of (2) surrounding the equilibria are discussed.Te stability of the equilibrium is demonstrated, in Section 4. Ultimately, in Section 5, we summarize this paper and provide some numerical simulations to demonstrate our theoretical results.

Preliminaries and Well-Posedness
For this second, model (2) will be converted as an abstract Cauchy problem (ACP); then we will prove the well-posedness of the system and discuss the nonnegativity and boundedness of solutions.For these intensives, we will frst gather some background information about linear operators and C 0 -semigroup theory, and some notations will be utilized throughout the paper.
We start by transforming the system (2) into an abstract evolution equation.Defne that and give the linear operator's formulation { }, which implies that D(A) is not dense in X.A nonlinear operator L: D(A) ⟶ X given by the following: and let T(t) � (S(t), E(t, •), I(t), Q(t), R(t), 0) T .Hence, based on the above notations, system (2) can be rewritten into an abstract Cauchy problem as follows: with   Generally speaking, solving an abstract diferential problem with a strong solution such as (7) is challenging.As a result, we solve (7) in integrated form as follows: Set In Section 3, the operator (A, D(A)) will be shown to be a Hille-Yosida operator, and thus, it produces a C 0 -semigroup on the closure of its domain, according to Teorem 3.8. in [36].Next, the well-posedness result of the model ( 7) can be obtained.

From biological interpretation, S(t), I(t), Q(t), and R(t)
represent the number of susceptible persons, infectious persons, quarantined persons, and recovered persons at time t, respectively.E(t, a) represents the number of latent persons with age a. Terefore, only nonnegative solutions (S(t), E(t, a), I(t), Q(t), R(t)) of the system (2) are meaningful.Te following result demonstrates that the nonnegative initial value solutions (S(t), E(t, a), I(t), Q(t), R(t)) of ( 2) remain nonnegative and bounded eventually.
Theorem .For any t ≥ 0, all the solutions of system (2) with nonnegative initial values stay nonnegative and are eventually bounded.
Proof.To begin, by integrating the second equation of ( 2) along the characteristic line, we can get the following result: Similar to the proof of Lemma 2.2 in [27], it is easy to check that E(t, a) remains nonnegative for nonnegative initial data.
Next, before proving Tus, the third equation of model ( 2) can be replaced as follows: Ten we have Trough the fourth equation of system (2), we have thus, Similarly, from the ffth equation, we can obtain Ten, we demonstrate for ∀t ≥ 0, S(t) is nonnegative.If we assume there ∃t 1 > 0 such that S(t 1 ) � 0, and S(t) > 0 for ∀t ∈ (0, t 1 ).In fact, by the frst equation of system (2), there is S ′ (t 1 ) � Λ + rR(t 1 ) > 0, which obviously is inconsistent.Tus, S(t) ≥ 0 for any t ≥ 0.
Finally, the boundness of the solutions of model (2) will be proved.Te total number of latent persons is represented by Tere is a biologically defned maximum age; hence, lim a ⟶ +∞ E(t, a) � 0 is a plausible assumption.From system (2), we deduce Terefore, As a result, the omega limit set of model ( 2) is constrained in the below feasible region: which is obviously a positive invariant for model (2), implying that system (2) is ultimately bounded.

Equilibria and Properties of Linearization at Equilibria
Te nonlinear system (2) will be linearized around the equilibrium solutions in this section.In order to do this, we frst consider whether or not the equilibrium exist.System (2) always maintains a disease-free equilibrium E 0 � (S (0) , 0, 0, 0, 0) with S (0) � Λ/μ.For the purpose of fnding the nontrivial equilibrium 2), the following equations are proposed: Discrete Dynamics in Nature and Society We solve the above equation in (17), assume p � μ + ι 1 + −  a 0 δ(θ)dθ da, and then obtain On the other hand, by the frst equation ( 17) and conducting some computations, we obtain that Tus, We propose the defnition of the basic reproduction number as follows: R 0 is defned to be the expected number of secondary cases produced in a completely susceptible population by a typical infected individual during its entire period of infection [23].
R 0 can also be calculated by the recipe in van den Driessche and Watmough [37].Moreover, Diekmann et al. [38] 2) has a single positive solution.As a result of the aforesaid study, we arrive at the following conclusion.
Theorem 3. Te model ( 2) always has a disease-free equilibrium E 0 � (S (0) , 0, 0, 0, 0).Furthermore, when Trough computations, the linearized system of (8) can be obtained as follows: in which 6 Discrete Dynamics in Nature and Society Tus, on X the compactness of the bounded linear DL(u) is obviously derived.Notice that Ω: � λ ∈ C: Re(λ) > − μ  , then the statement below can be proven.

Theorem 4. Te operator (A, D(A)) is a Hille
Tus, we know Integrating the second equation ( 26) with regard to the age variable a and summating all the equations, we get Discrete Dynamics in Nature and Society Tus, we have which implies that (A, D(A)) is a Hille ̶ Yosida operator.
By the proof of Teorem 4, the Hille ̶ Yosida is estimated; thus, we know Moreover, DL(u)S(t): X 0 ⟶ X is compact for any t > 0. Since which implies (T(t)) t≥0 is a quasi-compact C 0 -semigroup.Tus, Teorem B.1 in [39] deduces that, for some η > 0, e ηt ‖T(t)‖ ⟶ 0 as t ⟶ + ∞ when any eigenvalue of (A + DL(u)) has a negative real part.Now, we may derive the result below, based on the previous arguments.
(i) Te steady-state u is locally asymptotically stable if any eigenvalue of (A + DL(u)) has a strictly negative real part.(ii) If at least one of the eigenvalues of (A + DL(u)) has a strictly positive part, then the steady state u is not stable.

Stability of Equilibria
Depending on the preceding analysis, the local stability of disease-free equilibrium E 0 will be frstly analyzed in this section.
Proof.Based on the result of Teorem 6, here we only need to demonstrate that for every nonnegative solution (S(t), E(t, a), Indeed, (15) implies that As a result, for ∀ε > 0, there is t 0 such that S(t) ≤ S (0) + ε, for ∀t ≥ t 0 , and thus, Discrete Dynamics in Nature and Society Ten, consider the following linear system: We may conclude that system (39) permits a solution of the type where � E 0 (a), � I 0 , � Q 0 and � R 0 are positive and λ 0 is a root of the corresponding characteristic equation of this system by utilizing the similar way of proving Teorem 6. Te formulation (10) of solution E(t, a) for the second equation in the system (2) yields that E(t, a) ≤ � E(t, a), for t ≥ t 0 ; thus, Since R 0 < 1, all the eigenvalues of the system (39) have a negative real part, due to the continuous reliance of λ 0 on ε, there exists ε > 0 s.t.λ 0 < 0, indicating that lim Tus, the frst equation in ( 2) is asymptotic to equation Obviously, there is lim 0) .Using the asymptotic autonomous semifow theory (see Corollary 4.3 in [40]), we can fnd out that Te demonstration is now complete.
In the following, we primarily concern about the stability of endemic equilibrium E * � (S * , E * (a), I * , Q * , R * ).We linearize system (2) at endemic equilibrium E * to investigate its local stability.We do this by introducing the perturbation variables Discrete Dynamics in Nature and Society In order to analyze the stability of E * , we seek for solutions of the type x(t) � x 0 e λt , y(t, a) � y 0 (a)e λt , z(t) � z 0 e λt , m(t) � m 0 e λt , n(t) � n 0 e λt .Ten, the following eigenvalue issue is obtained: We obtain by solving the second and the third equations in (47): Discrete Dynamics in Nature and Society By the fourth and the ffth equations in (47), there are Te frst equation in (47) yields that By the last equation in (47), and then combining with (48)-(50), we obtain Tus, the characteristic equation ( 47) is derived as follows: (λ+δ(θ))dθ da < M ≤ 1, and r(c Tus, it is easy to demonstrate that the left side of the characteristic (52) holds for λ with Reλ ≥ 0: For the right side of (52), when Reλ ≥ 0, we obtain that where the notations of p, M, S * , R 0 defned at the beginning of Section 3. Terefore, for λ with Reλ ≥ 0 the left-hand side of ( 52) is strictly greater than 1, but the right-hand side of (52) is not greater than 1.Tat is an inconsistency, i.e., the solutions of the characteristic (52) must have a negative real part.As a result, the endemic equilibrium E * is locally asymptotically stable, which allows us to derive the following conclusion.
As the system ( 2) is a high-dimensional system, it is hard to estimate the global stability of the endemic equilibrium E * .However, we can verify the global stability of E * for the situation r � 0. In this case, the function g(x, y) and the Lyapunov functional are defned as follows: where Using ( 17), we take the derivative of each portion of the Lyapunov functional F described in (57) along with the solutions of system (2) separately.Firstly, we obtain Discrete Dynamics in Nature and Society By using (10), we can rewrite Secondly, using ( 18) and the fact xg x (x, y) + yg y (x, y) � g(x, y) yield (60) Denote that B(0) � S * /S (0) R 0 � 1, and Ten, we have Diferentiating F 3 ′ along the solutions of system (2) yields 14 Discrete Dynamics in Nature and Society Adding ( 58), (62), and (63) together yields Terefore, we have where For any x > 0, l(x) ≤ 0 with equality holding if and only if x � 1. Tus, H(t, a) ≤ 0. Hence, F ′ ≤ 0. Ten, it can be confrmed that the singleton E *   is the biggest invariant set of F ′ � 0. According to [41], we know that E * is globally asymptotically stable, as shown by the compact global attractor A 0 � E *  .Te preceding considerations lead to the following conclusion.
Theorem 9. Suppose R 0 > 1, and r � 0, the endemic equilibrium E * of system ( 2) is globally asymptotically stable.Remark 1.For R 0 > 1, and r > 0, it is difcult to construct a suitable Lyapunov function to verify the global stability of the endemic equilibrium E * .Here, we take r � 0.025, 0.05 and use numerical simulation to display the endemic equilibrium E * .Compared with the case r � 0, we fnd that the endemic equilibrium E * is still globally asymptotically stable (see Figure 3).Moreover, it can be seen that with the increasing of, the number of latent and infectious persons increase.Hence, the rate of loss of immunity r has a negative impact on the control of the disease.

Conclusion and Discussions
In this research, an age-structured SEIQR model of HFMD has been developed and studied, with the infection age and time lag.Te well-posedness of the system (2) is explored frst, and following that the solutions' positivity and boundedness for system (2) is discussed.Next, the basic reproduction number R 0 is derived, which has been shown to be the threshold for determining disease extinction or survival.Tat is to say, when R 0 < 1, the disease-free equilibrium E 0 is globally asymptotically stable; otherwise, the system is unstable.Moreover, we explore the local stability of the infected equilibrium E * by analyzing the distribution of roots to the related characteristic equation and the global stability by constructing suitable Lyapunov functions for the epidemic model (2) when r � 0, and demonstrate the E * is still globally stable when r ≠ 0 by numerical simulation.Tus, by managing the basic reproduction number R 0 , we may control the spread of HFMD in the population.
In the following, we utilize numerical simulations to demonstrate the theoretical results of the system (2) by Matlab software and reveal the infuences of age-structure, saturation incidence σ and time delay τ on the disease dynamics.Firstly, according to Chinese real data about HFMD, we choose the parameters in the system (2).Te entire population of China in 2020 is 1412120000 and the birth rate is 8.52 per thousand [42]; thus, Λ � 1002605.2.And the natural mortality rate of China is 7.07 per thousand, thus μ � 0.00707.From [35], the number of infected cases in 2018 is about 2353310.We assume that q � 0.45, thus I(0) � 1294321, Q(0) � 1058989.Te report shows that the exposure period for HFMD is about 4-7 days; thus, the initial latent population is E(0) � 392218.We just supply an a priori estimate for the initial recovered population R(0) � 4532432.Te ratio of the entire youthful population aged from 1 to 14 is 17.9; therefore, S(0) � 252770000.We may derive the death rate of infection ι 1 � 0.000025.Te decubation of infection was about 2 weeks [43,44], that is, the recovery rate of infectious individuals c 1 � 0.0714.Besides, in all the following fgures, we set E(t): �  +∞ 0 E(t, a)da represents the total number of latent individuals.
(1) Te maximum infection age we chose is 100, β � 0.0005, ι 2 � 0.000015, c 2 � 0.075, r � 0.025, τ � 0.5, θ � 0, α * � 0.35 reference from [23].Accordingly, we obtain R 0 < 1. Te numbers of the susceptible persons, the latent persons, the infectious persons, the quarantined persons, and the recovered persons with three diferent initial values are displayed in Figure 4, implies that the disease-free equilibrium E 0 of system (2) is globally asymptotically stable when R 0 < 1.It is compatible with Teorem 6, which indicates that if we control the fundamental regeneration number R 0 to be kept below 1, the illness will be eradicated.
(2) When R 0 < 1, the impact of delay "τ" on the solutions of system (2) are displayed in Figure 5.We take τ � 0, 0.5, 0.8, respectively, and keep the other parameters the same with Figure 4. Figure 5 states that as the delay τ becomes longer, the time for achieving disease-free equilibrium E 0 increases.(3) Under R 0 > 1, the stability of positive equilibrium E * will be shown next.When we take the parameters β � 0.22, and others keep the same with Figure 4, then we get R 0 > 1. Figure 6 illustrates that the solutions of system (2) with three distinct initial values will trend to E * as t approaches infnity, revealing that the positive solution E * probably be globally asymptotically stable when R 0 > 1. (4) Fourthly, Figure 7 shows the infuence of delay "τ" on the solutions of system (2) when R 0 > 1. Te delay τ does not infuence the stability of E * , but from Figure 7 we can see that as the delay τ (in the incidence of latent infections) increases, the number of HFMD infections, quarantines, and recoveries decreases.6) Finally, we examine that the efect of the saturation incidence rate parameters of the system (2) on the dynamics of HFMD when R 0 > 1. Te impact of the saturation incidence rate on the solutions of the system (2) are displayed in Figure 9.When σ � 0.195, 0.2, 0.205, respectively, and the other parameters are the same as in Figure 6, then R 0 > 1.
From Figure 9, we know that as the saturation incidence rate's parameters σ increases, the number of the susceptible persons S(t) raises, as well as the number of the latent persons E(t), infectious persons I(t), quarantined persons Q(t), and recovered persons R(t) all decrease.In summary, the higher the saturation incidence rates's parameter σ is, the faster the HFMD can be controlled.
Comprehensive to is it can know the growth of time delay τ adverse to the formation of disease-free equilibrium E 0 when R 0 < 1, and causes the increasing of the number of HFMD infections, quarantines and recoveries decrease when R 0 > 1. Terefore, the increase in the saturation incidence rate σ is beneft for the control of HFMD.

Figure 2 :
Figure 2: Flow chart of the SEIQR model for HFMD.

Figure 1 :
Figure 1: Te ftting curve of the transmission rate E(a) of latent HFMD individuals who are 0-15 years old.

Figure 3 :
Figure 3: Te solutions of system (2) with diferent r when R 0 > 1. S(t) represents the susceptible individuals, E(t) �  +∞ 0 E(t, a)da represents the total number of latent individuals, I(t) represents the infectious individuals, Q(t) represents the quarantined individuals, and R(t) represents the recovered individuals.

Figure 4 : 8 Figure 5 :Figure 6 :
Figure 4: Te graph of disease-free equilibrium E 0 of system (2), where S(t) represents the susceptible individuals, E(t) �  +∞ 0 E(t, a)da represents the latent individuals, I(t) represents the infectious individuals, Q(t) represents the quarantined individuals, and R(t) represents the recovered individuals.Te three diferent color curves correspond to three sets of diferent initial values.

( 5 )
Accordingly, the distribution of latent individuals with respect to infection age at the endemic equilibrium E * , E * (a) is shown in Figure 8(a), which corresponding to the second solution line in Figure 6, and the distributions with both infection age and time, E(t, a) is shown in Figure 8(b).(

Table 1 :
Te number of latent HFMD individuals in China.

Table 2 :
Te biological meanings of parameters in the system (2).
have proposed that R 0 is the spectral radius of the next-generation matrix.In fact, each phrase in R 0 has a distinct epidemiological + ι 1 + c 1 + q is the average infection period.β is the transmission rate of the infectious individual.S (0) is the overall number of susceptible individuals.σ is the saturation incidence.As a result, R 0 represents new cases generated by the average number of typical infected members during the course of the infection period.As determined by the preceding study, the basic reproduction number of system (2) is R 0 .