Dynamical behavior of a coupling SEIR epidemic model with transmission in body and vitro, incubation and environmental effects

1.   Introduction

Infectious diseases have serious impacts on human health, social stability, economic development, happiness of families, and even national security. In particular, humans beings are facing serious threats from various kinds of viruses, e.g., AIDS, influenza, hepatitis B, dengue, cholera, MERS, SARS, Ebola and COVID-19 ^{[1,2,3,4,5,6]}. Therefore, it is vitally important to study the dynamic properties of viral spread, which could provide evidently theoretical guidance for the prediction and control of disease transmission.

During the past few decades, many epidemic models of various viruses have been extensively investigated by medical scientists, and bio-mathematicians, and many very important results have been obtained, e.g., ^{[7,8,9,10,11,12,13,14,15,16,17,18,19,20,21]} and the references cited therein. In particular, the dynamics of infectious diseases driven by environmental viruses have gradually become central issues. Therein, the epidemiological processes among several compartmental classes of the population, and the virus infection processes in hosts are coupled by two kinds of models at quick and slow time scales, respectively, which are very different from most existing models that study these two processes separately and do not explicitly connect them. During the disease transmission process, the environmental contamination rate depends not only on the number of infected individuals but also on the density of the virus, which makes it more necessary to study the interaction between these two different kinds of transmitters of disease by using coupling models, i.e., the transmission model of a virus in a single host and the disease spread model between hosts (see, e.g., ^{[7,8,9,11,12,13,14,15,16,20,21]} and the references cited therein).

Recently, epidemic models with two processes of within-host virus infection and between-host disease transmission have been investigated in many articles. (See ^{[9,11,12,13,14,15,20,21]}). In Particular, Wen et al. ^[13] investigated a discrete time environmentally-driven coupling dynamic model of within-host virus infection and between-host disease transmission. Wen et al. ^[22] further improved the above discrete time model into a more practical form with saturation incidence. In ^[20,21], the authors also investigated age-structured and reaction-diffusion epidemic models for coupling within- and between-host dynamics in environmentally-driven infectious diseases.

In many infectious diseases, such as influenza, COVID-19, the incubation period cannot be neglected. These diseases have no obvious symptoms in the early stages of infection, some symptoms only appear after a period of incubation, and some of them are contagious not only during the infectious period but also during the incubation period^[23,24,25]. In particular, COVID-19 is contagious during its incubation period. With the idea of infectivity in the incubation period, Jiao et al. ^[26] proposed a SEIR epidemic model with homestead-isolation on the susceptible population. However, for some diseases, such as tuberculosis, measles and AIDS, that are transmissible up on adequate contact with an infected individual, a susceptible individual becomes lurk, that is, infected but not infective. This individual remains in the incubation period before becoming infective.

We already know that the immune response has a great role in controlling the spread of disease during the virus infection process ^{[27,28,29,30]}. Particularly, when antibodies in a host do not work or the effect of antibodies is weak, viruses will continue to spread. As a result, the spread of disease between hosts will be difficult to control. Otherwise, when a host has an extensive antibody response, the production rate of B cells in the host will increase, and environmental virus invasion into the host will be prevented; then, the transmission of disease is well controlled. Humoral immunity is much more efficient than cellular immunity in some infections^[31]. Murase et al. ^[32] introduced a model that described the interaction among target cells, pathogens and immune responses. Wang et al. ^[27] applied this model to investigate the dynamic behavior of in host virus infection models with humoral immunity and intracellular delays.

Furthermore, many works have realized that the interaction driven by the environment between two kinds of spread routes, i.e., virus transmission (or replication) within the body and between hosts, has a significant effect on the spread of diseases ^{[7,8,9,10,11,12,13,14,15,16,20,21]}. In Particular, asymptomatic infected individuals will make the control of disease spread more difficult^[23,25]. Based on the above consideration, in this paper, we propose an SEIR epidemic model coupling virus spread within the body and between hosts, incubation, antibody response and environmental effects. The paper is organized as follows: In Section 2, the virus infection model with humoral immunity in the host, and the SEIR epidemic model of in vitro transmission between hosts are proposed with incubation at fast and slow time scales respectively, and both models are coupled by the environmental concentration of the virus. In Section 3, some properties and important results of the virus infection model with humoral immunity are given. In Section 4, the positivity and boundedness of solutions for the SEIR model without antibody responses are discussed, criteria on the existence and local stability of equilibria are established, and the uniform persistence of the model is obtained. In Section 5, we further study the SIER model with an antibody response. The basic reproduction number is calculated, and the existence and global stability of equilibria are studied. In Section 6, the main conclusions are illustrated by numerical examples and discussion. In Section 7, we give a brief conclusion.

2.   Model description

In this section, we provide a detailed description of the SEIR epidemic model coupling virus transmission in the body and between hosts, incubation and environmental effects. We first proposed the following assumptions.

$(A_1)$ The transmission of diseases is mainly caused by susceptible individuals being exposed to a virus in the environment. This kind of infection is commonly called indirect infection. The environmental contamination rate is related to the quantity of infected individuals and the virus load $V$ of hosts, with the form $\theta{V}I$.

$(A_2)$ Susceptible individuals have an incubation period after being infected by the virus. Usually, different lurks have different incubation periods. To facilitate discussion, we take an average incubation period provided that the infected individuals will not transmit disease in the incubation period.

$(A_3)$ The infected individuals usually will either die or be cured as a remover after being infected. To facilitate discussion, we assume that the removed individuals will not become susceptible again or that the disease will not recur.

$(A_4)$ The Viruses in the environment primarily comes from unique releasers, i.e., infected individuals. The amount of virus released by different infected individuals is related to the total load of virus within the bodies of hosts. To facilitate discussion, we take the average number of viruses within the infected individuals.

Based on the above assumptions, we can first establish the following SEIR epidemic model of viral infection:

where $S(t)$, $E(t)$, $I(t)$ and $R(t)$ indicate the numbers of susceptible, latent, infected and removed individuals at time $t$, $U(t)$ is the contamination rate of the virus in the environment at time $t$, $\beta$ indicates the probability of susceptible individuals being exposed to environmental viruses, $\Lambda$ indicates the recruitment rate of susceptible individuals, $\mu$ indicates the natural mortality rate of the entire population (including susceptible, latent, infected, and removed individuals), $\alpha$ indicates the conversion rate of lurks into infected individuals, $\sigma$ indicates the cure rate of infected individuals, $\zeta$ indicates the disease-related mortality rate of infected individuals, $\theta$ indicates the emission rate of virus to the environment released by each infected individual, $\theta VI(t)(1-U(t))$ indicates the increase in virus concentration in the environment per unit time as a whole, $\xi$ indicates the decay rate of virus in the environment, and $\gamma$ is the per capita virus clearance rate in the environment.

In the above model, $V$ represents the amount of virus carried in infected individuals. To obtain the specific form of $V$, we need to further study the dynamic process of virus infection in host. Dynamic models of virus infection in different types of hosts have been investigated ^{[9,11,12,13,20,21]}. In this article, we further propose a virus infection model with humoral immunity in the host as follows:

where $s$ represents the dynamic evolution time of cells, viruses, and $B$ cells in hosts, with $T = T(s)$, $T^* = T^*(s)$, $V = V(s)$ and $B = B(s)$ being the density of healthy cells, infected cells, virus load and $B$ cells at a fast time scale $s$, respectively. The parameters $\Lambda_c$, $k$, $m$, $d$, $p$ and $c$ are all positive constants. $\Lambda_c$ indicates the recruitment rate of healthy cells, $k$ indicates the infection rate of cells, $m$ indicates the natural mortality rate of cells, $d$ indicates the infection mortality rate of cells, $p$ indicates the rate of virus release by infected cells, $c$ indicates the virus clearance rate in hosts, $q$ indicates the neutralization rate of $B$ cells, $h$ indicates the production rate of $B$ cells, and $\omega$ is the mortality rate of $B$ cells.

In this article, the influence of environmental viruses is also introduced. The number of viruses in the body of a host will increase due to contact with environmental viruses, breathing and diet. Therefore, we improve model (2.2) and change it into the following form:

where function $g(U(t))$ represents the increase in the number of viruses in the host caused by environmental viruses invading into the host per unit time. Virus invasion usually indicates that the virus invades the host through breathing, eating, or contact with air, water, food or other objects contaminated by environmental viruses. Obviously, as the number of environmental viruses increases, the number of viruses that invade the host will also increase.

Thus, we finally establish the SEIR model coupling virus transmission both in the body and between hosts, incubation and environmental effects as follows:

There are two time scales in model (2.4). One is the dynamic evolution time $s$ of susceptible cells, infected cells, viruses and B cells within the host. The other is the dynamic evolution time $t$ of susceptible, latent, infectious, and recovered individuals and the level of environmental in vitro contamination of hosts. Time $s$ is a fast time scale within hosts, and time $t$ is a slow time scale between hosts. Usually, $s$ is faster than $t$; then, we can assume that $t = \omega{s}$, where $\omega$ is a very small positive number. However, there is a variable $g(U(t))$ in subsystem (2.3) within hosts and a variable $V(s)$ in subsystem (2.1) between hosts. In this case, model (2.4) can be seen as a coupling models (2.1) and (2.3) in the body and between hosts.

3.   Virus infection model with humoral immunity in a host

For the virus infection model (2.3), since the virus infection change process in the host is much faster than the disease transmission process between hosts, we can assume that the environmental contamination rate $U(t)$ in model (2.3) remains constant $U\; (0\leq U \leq 1)$. Thus, model (2.3) becomes an isolated virus infection dynamics model.

We further suppose that environmental contamination is mainly related to the density of viruses living in polluted environment, and the host is infected through the ingestion of contaminated food. Therefore, function $g(U)$ is assumed to satisfy the following condition:

$({\bf{H}}_{\bf{1}})$ $g'(U) > 0$, $g''(U)\leq 0$ for all $U\geq 0$, and $g(0) = 0$.

There are some special forms of $g(U)$ satisfy condition $({\bf{H}}_{\bf{1}})$. For example, the linear function $g(U) = aU$ used in ^[9,11,12], means that the amount of environmental virus entering the host linearly increases with $U$. The nonlinear function $g(U) = \frac{a(U)}{1+bU}$ is used in ^[12], where $a > 0$, $b > 0$, which means that the amount of environmental virus entering the host will reach saturation.

From the biological significance of model (2.3), any solution $(T(s), T^*(s), V(s), B(s))$ is assumed to satisfy the initial conditions:

Based on the positivity and boundedness of solutions, the conclusion is as follows.

Theorem 1. Any solution of model (2.3) with initial conditions (3.1) is positive and ultimately bounded for all $s\geq 0$.

The proof of Theorem 1 is not difficult, so we omit it here.

When $U = 0$, model (2.3) degrades into model (2.2). The basic reproduction number of virus infection and antibody response for model (2.2) are defined, respectively, by

From ^{[9,11,12,13,20,21]}, we can summarize the following conclusions.

Theorem 2. Assume that $U = 0$ in model (2.3).

(1) Model (2.3) always has infection-free equilibrium $A_0 = (T_0, 0, 0, 0)$, where $T_0 = \frac{\Lambda_c}{m}$;

(2) If $R_{w0} > 1$, then model (2.3) has no antibody response infection equilibrium $A_1 = (T_1, T^*_1, V_1, 0)$, where, $T_1 = \frac{T_0}{R_{w0}}, \; T^*_1 = \frac{mT_0}{m+d}(1-\frac{1}{R_{w0}})$ and $V_1 = \frac{m(R_{w0}-1)}{k};$

(3) If $R_{w0} > 1$ and $R_{w1} > 1$, then model (2.3) has antibody response infection equilibrium $A_2 = (T_2, T^*_2, V_2, B_2)$, where, $T_2 = \frac{h\Lambda_c}{k\omega+mh}, \; T^*_2 = \frac{k\omega\Lambda_c} {(m+d)(k\omega+mh)}, \; V_2 = \frac{\omega}{h}$ and $B_2 = \frac{pk\Lambda_ch -c(m+d)(k\omega+mh)}{q(m+d)(k\omega+mh)}.$

Theorem 3. Assume that $U = 0$ in model (2.3). Then, the following conclusions hold:

(1) If $R_{w0}\leq 1$, then infection-free equilibrium $A_0 = (T_0, 0, 0, 0)$ is globally asymptotically stable;

(2) If $R_{w0} > 1$ and $R_{w1}\leq 1$, then the antibody-free infection equilibrium $A_1 = (T_1, T^*_1, V_1, 0)$ is globally asymptotically stable;

(3) If $R_{w0} > 1$ and $R_{w1} > 1$, then the antibody response infection equilibrium $A_2 = (T_2, T^*_2, V_2, B_2)$ is globally asymptotically stable.

When $U > 0$, the basic reproduction number of the antibody response to the infection for model (2.3) is defined by

Obviously, when $U = 0$, $R_w = R_{w1}$. The existence of equilibria for model (2.3) is given in the following conclusion.

Theorem 4. Assume that $U > 0$ in model (2.3).

(1) Model (2.3) always has antibody-free infection equilibrium $A_3(U) = (T_3(U), T^*_3(U), V_3(U), 0)$, where $T_3(U)$, $T^*_3(U)$ and $V_3(U)$ are given below. Furthermore,

(2) If $R_w > 1$, then model (2.3) has antibody response infection equilibrium $A_4(U) = (T_4(U), T^*_4(U), V_4(U), B_4(U))$, where $T_4(U)$, $T^*_4(U)$, $V_4(U)$ and $B_4(U)$ are given below. Furthermore,

Proof. Clearly, equilibrium $A_3(U) = (T_3(U), T^*_3(U), V_3(U), 0)$ satisfies the following equations:

We directly obtain

and $T_{3}(U) = \frac{m}{m+kV_3(U)}T_0 < T_0$. Substituting (3.4) into the first equation of (3.3), we further obtain the quadratic equation with $T_3(U)$ as a root:

where $a_1 = \frac{(m+d)g(U)}{pm}+T_0(1+\frac{1}{R_{w0}})$ and $a_2 = \frac{T^2_0}{R_{w0}}.$ By direct calculation, Eq (3.5) has a unique root $T_{3}(U) = \frac{1}{2}(a_1-\sqrt{a^2_1-4a_2})$, which satisfies $T_{3}(U) < T_0$. Thus, model (2.3) has a unique no antibody response infection equilibrium $A_3(U) = (T_3(U), T^*_3(U), V_3(U), 0)$ with

When $R_w > 1$, with calculation, model (2.3) has a unique antibody response infection equilibrium $A_4(U) = (T_4(U), T^*_4(U), V_4(U), B_4(U))$ with

This completes the proof.

Regarding the stability of equilibria $A_3(U)$ and $A_4(U)$, the following conclusion is established.

Theorem 5. Let $U > 0$ in model (2.3).

(1) If $R_w\leq 1$, then theantibody-free infection equilibrium $A_3(U)$ is globally asymptotically stable.

(2) If $R_w > 1$, then the antibody response infection equilibrium $A_4(U)$ is globally asymptotically stable.

Proof. For equilibrium $A_3(U)$, we define the Lyapunov function $L_3$ as follows:

For equilibrium $A_4(U)$, we define the Lyapunov function $L_4$ as follows:

Computing the derivatives of $L_3$ and $L_4$ along the solution of model (2.3), we can obtain

When $R_w\leq 1$, $V_3\leq\frac{\omega}{h}$. Hence, for any $(T, T^*, V, B)\in R_+^4$, $\frac{dL_3}{dt}\leq 0$. Moreover, $\frac{dH_3}{dt} = 0$ implies that $T = T_3$, $T^* = T^*_3$ and $V = V_3$. Furthermore, $B = 0$ when $V = V_3$. Thus, from the LaSalle invariance principle ^[33], equilibrium $A_3(U)$ is globally asymptotically stable.

When $R_w > 1$, $\frac{dL_4}{dt}\leq 0$ for any $(T, T^*, V, B)\in R_+^4$. Moreover, $\frac{dH_4}{dt} = 0$ implies that $T = T_4$, $T^* = T^*_4$ and $V = V_4$. In addition, $B = B_4$ when $V = V_4$. Thus, from the LaSalle invariance principle ^[33], equilibrium $A_4(U)$ is globally asymptotically stable. This completes the proof.

4.   SEIR model without antibody response

4.1.   Positivity and boundedness

Here, we discuss the nonnegativity and boundedness of solutions for model (2.4) in the general case, so we discuss neither the antibody response, nor the antibody dose response. Since the viral infection change process in the host is much faster than the disease transmission process between the hosts, we can assume that the state of virus infection in hosts has reached its equilibrium while the state of disease transmission between the host has not changed. Therefore, we can assume that the fast time variable $V(s) = \bar{V}(U)$ in model (2.4); furthermore we have the following form of the coupled model:

According to Theorem 4, we have $\bar{V}(U)$ with the following expression,

For the convenience of statements, we denote $R_+^{5}: = \{(x_1, x_2, x_3, x_4, x_5): x_i\geq 0, \; i = 1, 2, 3, 4, 0\leq x_5\leq 1\}$. Due to the biological significance of model (4.1), any solution $(S(t), E(t), I(t), U(t))$ of model (4.1) satisfies the following initial condition:

where $(S_0, E_0, I_0, R_0, U_0)\in R_+^{5}$ and $S > 0$.

Regarding the nonnegativity and boundedness of solutions for model (4.1), we have the following conclusions.

Theorem 6. For any initial point $x_0 = (S_0, E_0, I_0, R_0, U_0)\in R^{5}_{+}$, model (4.1) has a unique nonnegative solution $u(t) = (S(t), E(t), I(t), R(t), U(t))\in R_+^5$ defined on $[0, \infty)$, and the solution is also ultimately bounded. Furthermore, $0\leq U(t)\leq 1$ for all $t\geq 0$.

Proof. Let solution $(S(t), E(t), I(t), R(t), U(t))$ be defined for $t\in[0, \tau_{\infty})$ with $\tau_{\infty}\leq\infty$. We first prove the nonnegativity of the solution. Assume that $(S_0, E_0, I_0, R_0, U_0) > 0$. If we define $m(t) = \min\{S(t), E(t), I(t), R(t), U(t)\}$, then $m(0) = \min\{S(0), E(0), I(0), R(0), U(0)\} > 0$. Now, we need to verify $m(t) > 0$ for any $t > 0$. If we assume that there is a $t_0\in(0, \tau_{\infty})$ satisfying $m(t_0) = 0$ and $m(t) > 0$ for any $t\in [0, t_0)$. If $m(t_0) = S(t_0)$, then from the first equation of model (4.5), we easily have

It follows that $S(t_0)\geq S(0)e^{-\int_0^{t_0}(\beta U(u)+\mu)du} > 0,$ which leads to a contradiction. Similarly, if $m(t_0) = E(t_0)$, $m(t_0) = I(t_0)$, $m(t_0) = R(t_0)$ or $m(t_0) = U(t_0)$, we also obtain a contradiction. Therefore, $m(t) > 0$ for any $t\in[0, \tau_{\infty})$. If we assume that $(S_0, E_0, I_0, R_0, U_0)\geq0$, then, from the continuity of solutions with respect to initial values, we obtain $(S(t), E(t), I(t), R(t), U(t))\geq 0$ for any $t\in[0, \tau_{\infty})$.

Next, we prove the boundedness of the solutions. Let $N(t) = S(t)+E(t)+I(t)+R(t)$; then, we have

This shows that $N(t)$ is bounded for any $t\in[0, \tau_{\infty})$. Now we consider $U(t)$. Assume that $U(0) < 1$. If there is a $t_0\in(0, \tau_{\infty})$ satisfying $U(t_0) = 1$ and $U(t) < 1$ for any $t\in [0, t_0)$, then $\frac{dU(t_0)}{dt}\geq 0$. However, from the fifth equation of model (4.1), we have $\frac{dU(t_0)}{dt} = -(\xi+\gamma)$, which causes a contradiction. Therefore, $U(t) < 1$ for any $t\in[0, \tau_{\infty})$. Furthermore, assume that $U(0)\leq 1$; then, by the continuity of solutions with respect to initial values, we can obtain $U(t)\leq 1$ for any $t\in[0, \tau_{\infty})$. Thus, we finally determine that solution $(S(t), E(t), I(t), R(t), U(t))$ is defined for all $t\geq 0$. That is, $\tau_{\infty} = \infty$.

Again, from (4.4), we have $\limsup_{t\rightarrow \infty}N(t)\leq \frac{A}{\mu}$. Specifically, we also obtain $N(t)\leq \frac{A}{\mu}$ for all $t\geq0$ while $N(0)\leq \frac{A}{\mu}$. This shows that solution $u(t)$ is also ultimately bounded. This completes the proof.

4.2.   Basic reproduction number and equilibria

Now, we consider the SEIR model without an antibody response. Because the removed $R$ does not emerge in the other four equations of model (4.1), we can only investigate the following sub-model

For the convenience of discussion, we always assume that $R_{w0} > 1$ and $R_{w}\leq1$ in this section. Thus, from (4.2), we can determine that the fast time variable $\bar{V}(U)$ in model (4.5) has the form $\bar{V}(U) = \frac{1}{c}(g(U)+pT^*_3(U))$. Moreover, from (4.2), we also have $\bar{V}(0) = \frac{m(R_{w0}-1)}{k}$. Now, we define the basic reproduction number for model (4.5) as follows:

Obviously, if $R_{w0} < 1$, then $R_{b0} < 0$; if $R_{w0} = 1$, then $R_{b0} = 0$; if $R_{w0} > 1$, then $R_{b0} > 0$; and if $R_{b0} > 1$, then $R_{w0} > 1$.

We rewrite $R_{b0} = \frac{1}{(\mu+\sigma+\zeta)}\cdot\theta\cdot\bar{V}(0)\cdot\frac{1}{(\xi+\gamma)}\cdot\beta\cdot S_0\cdot\frac{1}{(\mu+\alpha)}\cdot\alpha$. Where, $\frac{1}{(\mu+\sigma+\zeta)}$ indicates the average infected period of an infected individual. $\theta$ indicates the emission rate of the virus to the environment released by each infected individual. $\bar{V}(0)$ indicates the host's virus load in the initial stage of environmental pollution. $\frac{1}{(\xi+\gamma)}$ indicates the average survival time of the virus in the environment. $\beta$ indicates the probability of infection of susceptible individuals after contact with environmental viruses. $S_0$ indicates the total number of susceptible individuals in the initial stage of infection. $\frac{1}{\mu+\alpha}$ indicates the average latent period of latent individuals. $\alpha$ indicates the conversion rate of latent individuals into infected individuals.

Therefore, the basic regeneration number $R_{b0}$ denotes the number of new cases in which in the initial stage of infection, an infected individual releases viruses into the environment during the infection period. These viruses infect susceptible individuals during the survival period and make susceptible individuals become latent individuals, and last these latent individuals become infected individuals in the incubation period.

Define the functions $F(U) = \frac{1-U}{c}\big[g(U)+\frac{pm}{m+d}\big(T_0-T_3(U)\big)\big]$, $G(U) = \frac{1}{R_{b1}}+\frac{\beta U}{\mu R_{b1}}$ and $H(U) = F(U)-G(U)$, where $R_{b1} = \frac{\beta{S_0}\alpha\theta}{(\mu+\alpha)(\mu+\sigma+\zeta)(\xi+\gamma)}$. Clearly, $R_{b0} = R_{b1}\bar{V}(0)$.

Model (4.5) always has disease-free equilibrium $\tilde{W}_0 = (S_0, 0, 0, 0)$ with $S_0 = \frac{A}{\mu}$. Let $\tilde{W} = (\tilde{S}, \tilde{E}, \tilde{I}, \tilde{U})$ be the positive equilibrium of model (4.5) that satisfies the following equations:

Easily, we can obtain $\tilde{S} = \frac{A}{\mu+\beta\tilde{U}}$, $\tilde{E} = \frac{A\beta\tilde{U}}{(\mu+\alpha)(\mu+\beta\tilde{U})}$, $\tilde{I} = \frac{A\beta\alpha\tilde{U}}{(\mu+\alpha)(\mu+\beta\tilde{U})(\mu+\sigma+\zeta)}$ and $\tilde{U}$ is a zero point of equation $H(U) = 0$ in $(0, 1)$. Let $H_M = \max_{0\leq U \leq 1}{H(U)}$. Combining the basic reproduction number $R_{b0}$, assuming the existence of equilibria, we can obtain the following conclusion.

Lemma 1. (1) Model (4.5) always has disease-free equilibrium $\tilde{W}_0 = (S_0, 0, 0, 0)$ with $S_0 = \frac{A}{\mu}$.

(2) Model (4.5) has a unique endemic equilibrium $\tilde{W} = (\tilde{S}, \tilde{E}, \tilde{I}, \tilde{U})$ if and only if one of the conditions given below holds:

$(a) \; R_{b0} > 1; \; \; (b)\; R_{b0} = 1\; and\; H_M > 0; \; \; (c)\; R_{b0} < 1\; and\; H_M = 0$.

(3) Model (4.5) has two positive equilibria, $\tilde{W_1} = (\tilde{S_1}, \tilde{E_1}, \tilde{I_1}, \tilde{U_1})$ and $\tilde{W_2} = (\tilde{S_2}, \tilde{E_2}, \tilde{I_2}, \tilde{U_2})$ if and only if the condition given below holds:

$(d) \; R_{b0} < 1\; and\; H_M > 0$.

(4) Model (4.5) has only disease-free equilibrium $\tilde{W}_0 = (S_0, 0, 0, 0)$ if and only if one of the conditions given below holds:

$(e) \; H_M < 0;\; \; (f)\; R_{b0} = 1\; and\; H_M = 0$.

Proof. With calculation we obtain $H(0) = \bar{V}(0)(1-\frac{1}{R_{b0}})$ and $H(1) = -\frac{1}{R_{b1}}-\frac{\beta}{\mu R_{b1}} < 0.$ The second order derivative of $H(U)$ is given as follows:

From the assumption $({\bf H}_{\bf{1}})$, $a'_{1}(U) = g'(U)\frac{m+d}{pm} > 0$ and $a''_{1}(U) = \frac{(m+d)g''(U)}{pm}\leq 0$. In addition, $a_2 > 0$, where $a_1$ and $a_2$ are defined in (3.5). This implies that $T'_{3}(U) = \frac{1}{2}a'_{1}(U)(1-\frac{a_{1}(U)}{\sqrt{a^2_{1}(U)-4a_{2}}}) < 0$ and $T''_{3}(U) = \frac{1}{2}a''_{1}(U)[1-\frac{a_{1}(U)}{\sqrt{a^2_{1}(U)-4a_{2}}}]+\frac{1}{2}(a'_{1}(U))^2\frac{4a_{2}}{[a^2_{1}(U)-4a_{2}]^{\frac{3}{2}}} > 0$. Therefore, we further obtain $H''(U) < 0$ for any $U\geq0$. $H(U)$ is an upper convex function that has at most two zeros.

If condition $(a)$ is true, then due to $R_{b0} > 1$, $H(U) = 0$ has a unique positive root $\tilde{U}\in(0, 1)$. Therefore, endemic equilibrium $\tilde{W} = (\tilde{S}, \tilde{E}, \tilde{I}, \tilde{U})$ exists and is unique.

If condition $(b)$ holds, $H(0) = 0$ when $R_{b0} = 1$. From $H_M > 0$, $H(U) = 0$ has unique positive root $\tilde{U}\in(0, 1)$. Thus, endemic equilibrium $\tilde{W} = (\tilde{S}, \tilde{E}, \tilde{I}, \tilde{U})$ exists and is unique.

Under condition $(c)$, $H(0) < 0$ when $R_{b0} < 1$. $H(U) = 0$ has unique positive roots $\tilde{U}\in(0, 1)$ with $H(\tilde{U})$ when $H_M = 0$. Therefore, a unique endemic equilibrium $\tilde{W} = (\tilde{S}, \tilde{E}, \tilde{I}, \tilde{U})$ exists.

Let condition $(d)$ hold; then, due to $H(0) < 0$ and $H_M > 0$, $H(U) = 0$ has only two positive roots. Thus, model (4.5) has only two positive equilibria: $\tilde{W_1} = (\tilde{S_1}, \tilde{E_1}, \tilde{I_1}, \tilde{U_1})$ and $\tilde{W_2} = (\tilde{S_2}, \tilde{E_2}, \tilde{I_2}, \tilde{U_2})$.

If condition $(e)$ is true, owing to $H_M < 0$, $H(U) = 0$ has no roots in $(0, 1)$. Hence, only disease-free equilibrium $\tilde{W}_0 = (\frac{A}{\mu}, 0, 0, 0)$ exists.

Last, for condition $(f)$, then, due to $H(0) = 0$, $H(U) < 0$ for all $U\in(0, 1]$. Obviously, $H(U) = 0$ has only one root $U = 0$. Hence, only disease-free equilibrium $\tilde{W}_0 = (\frac{A}{\mu}, 0, 0, 0)$ exists. This completes the proof.

Remark 1. Conclusion (3) in Lemma 1 indicates that model (4.5) generates backward bifurcation at disease-free equilibrium $W_0$ when $R_{b0} < 1$.

4.3.   Stability of equilibria

We first consider the stability of the disease-free equilibrium $\tilde{W}_0 = (S_0, 0, 0, 0)$ of the model (4.5). The following conclusion is established.

Theorem 7. (a) If $R_{b0} < 1$, then equilibrium $\tilde{W}_0$ is locally asymptotically stable;

(b) If $R_{b0} > 1$, then equilibrium $\tilde{W}_0$ is unstable.

Proof. Calculating the Jacobian matrix of model (4.5) at any equilibrium $\tilde{W} = (\tilde{S}, \tilde{E}, \tilde{I}, \tilde{U})$, we obtain

where $J_{44}(\tilde{W}) = \theta\tilde{I}\big(\bar{V}'(\tilde{U})(1-\tilde{U})-\bar{V}(\tilde{U})\big)-(\xi+\gamma)$.

From (4.8), the Jacobian matrix at equilibrium $\tilde{W}_0$ is

Furthermore, the characteristic equation of $J(\tilde{W}_0)$ is

where

When $R_{b0} < 1$, we obtain $a_i > 0$ for $i = 0, 1, 2$. By calculation we further obtain

Thus, from the Routh-Hurwitz criterion ^[35], all roots of equation $f(\lambda) = 0$ have negative real parts. Hence, $\tilde W_0$ is locally asymptotically stable. When $R_{b0} > 1$, $a_0(\tilde{W}_0) < 0$; then, equation $f(\lambda) = 0$ has at least one positive root. Therefore, $\tilde W_0$ is unstable. This completes the proof.

Next, the stability of the endemic equilibrium $\tilde{W} = (\tilde{S}, \tilde{E}, \tilde{I}, \tilde{U})$ of model (4.5) is considered. From (4.8), we can see that the Jacobian matrix of model (4.5) at equilibrium $\tilde{W}$ is

where $J_{44}(\tilde{W}) = \theta\tilde{I}\big(\bar{V}'(\tilde{U})(1-\tilde{U})-\bar{V}(\tilde{U})\big)-(\xi+\gamma)$. From the fourth equation of equations (4.7), we can obtain $\theta{\tilde{I}} = \frac{(\xi+\gamma)\tilde{U}}{\bar{V}(\tilde{U})(1-\tilde{U})} = \frac{(\xi+\gamma)}{F(\tilde{U})}\tilde{U}$. Then, we have $J_{44}(\tilde{W}) = -(\xi+\gamma)\frac{F(\tilde{U})-\tilde{U}F'(\tilde{U})}{F(\tilde{U})}$. Denote $K(\tilde{U}) = F(\tilde{U})-\tilde{U}F'(\tilde{U})$; then, $K'(\tilde{U}) > 0$ for $\tilde{U} > 0$. Since $K(0) = \bar{V}(0)\geq 0$, $K(\tilde{U}) > 0$ for $\tilde{U} > 0$.

with calculation we can obtain the characteristic equation $J(\tilde{W})$ as follows:

where

Obviously, $b_i(\tilde{W}) > 0$ for $i = 2, 3$. Moreover, we can propose the following conclusion.

Theorem 8. Assume that $F'(\tilde{U})\leq 0$ and

and either conditions $(a)$ or $(b)$ in Lemma 1 holds. Then, unique endemic equilibrium $\tilde{W}$ is locally asymptotically stable.

Proof. Equations (4.7) imply $\tilde{S} = \frac{(\mu+\alpha)\tilde{E}}{\beta\tilde{U}} = \frac{(\mu+\alpha)(\mu+\sigma+\zeta)\tilde{I}}{\beta\tilde{U}\alpha} = \frac{(\mu+\alpha)(\mu+\sigma+\zeta)(\xi+\gamma)}{\beta\alpha\theta F(\tilde{U})}$. Then, we obtain

From $F'(\tilde{U})\leq 0$, $K(\tilde{U}) = F(\tilde{U})-\tilde{U}F'(\tilde{U})\geq F(\tilde{U})$ and $\frac{K(\tilde{U})}{F(\tilde{U})}-1\geq 0$. This implies that $b_i(\tilde{W}) > 0$ for $i = 0, 1$.

Next, we verify $b_3(\tilde{W})b_2(\tilde{W})-b_1(\tilde{W}) > 0$. which calculation, we obtain

where $A = \beta\tilde{S}\alpha\theta\bar{V}(\tilde{U})(1-\tilde{U})$. From the expression of $\tilde{S}$ and $F(\tilde{U}) = \bar{V}(\tilde{U})(1-\tilde{U})$, $A = (\mu+\alpha)(\mu+\sigma+\zeta)(\xi+\gamma) > 0$. Hence, we can obtain $b_3(\tilde{W})b_2(\tilde{W})-b_1(\tilde{W}) > 0$. Thus, from the Routh-Hurwitz criterion ^[35], all roots of characteristic equation $f_1(\lambda) = 0$ have negative real parts. Therefore, equilibrium $\tilde{W}$ is locally asymptotically stable. This completes the proof.

Remark 2. In the above conclusions, we show that the unique endemic equilibrium is locally asymptotically stable based on conditions $(a)$ and $(b)$ in Lemma 1. However, the same conclusion is not obtained for condition $(c)$ in Lemma 1. Therefore, we give an interesting open problem of whether the unique endemic equilibrium in condition $(c)$ is also locally asymptotically stable.

Remark 3. In Theorem 8, in addition to conditions $R_{b0} > 1$ or $R_{b0} = 1, H_M > 0$, additional hypotheses $F'(\tilde{U})\leq 0$ and condition (4.9) are also needed. Therefore, an interesting open question is whether (4.9) can be directly verified. The other open question is whether condition $F'(\tilde{U})\leq 0$ can be removed.

Remark 4. In Theorem 8, only the local asymptotic stability of equilibrium $\tilde{W}$ is established. Whether we can further construct the appropriate Lyapunov function and use the LaSalle invariant principle ^[33] or the geometric method ^[34] to obtain the global stability of equilibrium $\tilde{W}$ is still an open problem.

Now, the stability of positive equilibrium $\tilde{W}_1 = (\tilde{S_1}, \tilde{E_1}, \tilde{I_1}, \tilde{U_1})$ of model (4.5) is considered in the following theorem.

Theorem 9. If $R_{b0} < 1$ and $H_M > 0$, then positive equilibrium $\tilde{W_1}$ is unstable.

Proof. Using the same calculation used in Theorem 8, the characteristic equation of model (4.5) at equilibrium $\tilde{W_1}$ is given in the following form:

where

From conditions $R_{b0} < 1$ and $H_M > 0$, $H'(\tilde{U_1}) > 0$ and $H(\tilde{U_1}) = 0$. Furthermore, $F'(\tilde{U_1}) > G'(\tilde{U_1}) > 0$ and $F(\tilde{U_1}) = G(\tilde{U_1})$. Hence, since

we can obtain

Hence, the characteristic equation $f_2(\lambda) = 0$ has at least one positive root. Thus, positive equilibrium $\tilde{W_1}$ is unstable. This completes the proof.

Finally, the stability of the positive equilibrium $\tilde{W}_2 = (\tilde{S_2}, \tilde{E_2}, \tilde{I_2}, \tilde{U_2})$ of model (4.5) is discussed. Using the same calculation used in Theorem 8, we obtained the characteristic equation of model (4.5) at equilibrium $\tilde{W_2}$, which is given in the following form:

where

Obviously, $b_i(\tilde{W_2}) > 0$ for $i = 2, 3$. Furthermore, we can give the following theorem.

Theorem 10. Assume that $F'(\tilde{U}_2)\leq 0$ and

and condition $(d)$ in Lemma 1 holds. Then, positive equilibrium $\tilde{W_2}$ is locally asymptotically stable.

Proof. From the equations of which equilibrium $\tilde{W_2}$ is satisfied, $\tilde{S_2} = \frac{(\mu+\alpha)\tilde{E_2}}{\beta\tilde{U_2}} = \frac{(\mu+\alpha)(\mu+\sigma+\zeta)\tilde{I_2}}{\beta\tilde{U_2}\alpha}$. The fourth equation of model (4.5) indicates that $\theta{\tilde{I_2}} = \frac{(\xi+\gamma)\tilde{U_2}}{\bar{V_2}(\tilde{U_2})(1-\tilde{U_2})} = \frac{(\xi+\gamma)}{F(\tilde{U_2})}\tilde{U_2}$. In addition, $\tilde{S_2} = \frac{(\mu+\alpha)(\mu+\sigma+\zeta)(\xi+\gamma)}{\beta\alpha\theta F(\tilde{U_2})}$. Hence, we can obtain

From the assumption $F'(\tilde{U_2})\leq 0$, $K(\tilde{U_2}) = F(\tilde{U_2})-\tilde{U_2}F'(\tilde{U_2})\geq F(\tilde{U_2})$ and $\frac{K(\tilde{U_2})}{F(\tilde{U_2})}-1\geq 0.$ Finally, $b_i(\tilde{W_2}) > 0$ for $i = 0, 1$.

Next, we verify $b_3(\tilde{W_2})b_2(\tilde{W_2})-b_1(\tilde{W_2}) > 0$. With calculation, we obtain

where $B = \beta\tilde{S_2}\alpha\theta\bar{V}(\tilde{U_2})(1-\tilde{U_2})$. From the expression of $\tilde{S_2}$ and $F(\tilde{U}_2) = \bar{V}(\tilde{U_2})(1-\tilde{U_2})$, $B = (\mu+\alpha)(\mu+\sigma+\zeta)(\xi+\gamma) > 0$. Hence, we can obtain $b_3(\tilde{W_2})b_2(\tilde{W_2})-b_1(\tilde{W_2}) > 0$. Thus, from the Routh-Hurwitz criterion ^[35], all roots of characteristic equation $f_3(\lambda) = 0$ have negative real parts. Therefore, equilibrium $\tilde{W}_2$ is locally asymptotically stable. This completes the proof.

Remark 5. In Theorem 10, we see that in addition to conditions $R_{b0} < 1$ or $H_M > 0$, the additional hypotheses $F'(\tilde{U}_2)\leq 0$ and condition (4.10) are also required. Therefore, an interesting open question is whether we can get that (4.10) directly holds when $R_{b0} < 1$, $H_M > 0$ and $R_w\leq 1$. The other open question is whether condition $F'(\tilde{U}_2)\leq 0$ can be removed.

4.4.   Uniform persistence

Now, we discuss the uniform persistence of positive solutions for model (4.5). We can establish the following conclusion.

Theorem 11. Assume that $R_{b0} > 1$; then, there is a constant $\varepsilon > 0$ such that for any initial value $x_0 = (S_0, E_0, I_0, U_0)\in R^4_{+}$ with $E_0\ne0$, $I_0\ne0$ and $U_0\ne0$ solutions, the $u(t) = (S(t), E(t), I(t), U(t))$ of model (4.5) satisfies

Proof. For any initial point $x_0 = (S_0, E_0, I_0, U_0)\in R^4_{+}$ with $E_0\ne0, I_0\ne0$ and $U_0\ne0$, let $u(t) = (S(t), E(t), I(t), U(t))$ be the solution to model (4.5) satisfying initial condition $u(0) = x_0$. From Theorem 6,

By comparison principle, $\liminf_{t\rightarrow \infty}S(t)\geq\frac{A}{\mu+\beta}.$ This shows that $S(t)$ is uniformly persistent.

Define the set $X = \{ x = (S, E, I, U)\in R^4_{+}: E > 0, I > 0, U > 0\}.$ The boundary of $X$ is $\partial{X} = \{(S, E, I, U)\in R^4_{+}: E = 0\; or\; I = 0\; or\; U = 0\}.$ Denote $M_\partial = \{x_0\in R^4_{+}: u(t)\in\partial{X}, \forall\; t\geq0\}.$

Let $M_0 = \{W_0\}$. Clearly, $M_0\subset\cup_{x_0\in{M_\partial}}\omega(x_0)$, where $\omega(x_0)$ is the $\omega$-limit set of solution $u(t)$ with initial value $u(0) = x_0$. Since $x_0\in M_{\partial}$, and $u(t)\in\partial{X}$ for all $t\geq0$, we have $E(t)\equiv0$, $I(t)\equiv0$ or $U(t)\equiv0$. If $E(t)\equiv0$, then from the third and fourth equations of model (4.5), we know that $I(t)\equiv0$ and $U(t)\equiv0$. Thus, model (4.5) is reduced to the following equation:

From this, we can obtain $\lim_{t\rightarrow \infty}S(t) = S_0$, which implies that $\omega(x_0) = \{W_0\}$. If $I(t)\equiv0$, then from the second and fourth equations of model (4.5), we know that $I(t)\equiv0$ and $U(t)\equiv0$. Thus, model (4.5) is also reduced to Eq (4.11), which implies that $\omega(x_0) = \{W_0\}$. If $U(t)\equiv0$, then similar to the above discussions, we obtain $\omega(x_0) = \{W_0\}$. Therefore, we finally obtain $M_0 = \cup_{x_0\in{M_\partial}}\omega(x_0)$. Moreover, $M_0$ is isolated and noncyclic in $\partial{X}$.

Now we prove that $K^s(W_0)\cap{X} = \emptyset$, where $K^s(W_0)$ is the stable set of $W_0$. By contradiction, we assume that there is a $x_0\in X$ such that $\lim_{t\rightarrow \infty}u(t) = W_0$. Since $R_{b0} > 1$, we can choose a sufficiently small constant $\varepsilon > 0$ such that

Thus, there is a $t^* > 0$ such that $S(t)\geq{S_0-\varepsilon}$, $U(t) < \varepsilon$, $E(t) < \varepsilon$ and $I(t) < \varepsilon$ for all $t\geq t^*$. Furthermore, from $\lim_{U\to 0}\bar{V}(U) = \bar{V}(0)$, we can also obtain $\bar{V}(U(t)) > \bar{V}(0)-\varepsilon$ for all $t\geq t^*$. The following function is defined:

Then $\lim_{t\rightarrow \infty}L(t) = 0$. When $t\geq{t^*}$,

Clearly, $\frac{dL(t)}{dt} > 0$ for all $t\geq{t^*}$. Therefore, $L(t)$ is an increasing function of $t\geq{t^*}$. $\lim_{t\rightarrow \infty}L(t)\ne0$, which leads to a contradiction. Then, $K^s(W_0)\cap{X} = \emptyset$. According to theory of persistence for dynamic systems, there is a constant $\varepsilon$ such that for any $x_0\in{X}$,

Therefore, the uniform persistence of model (4.5) is obtained. This completes the proof.

5.   SEIR model with antibody response

In this section, the SEIR model with an antibody response is considered. We always assume that $R_{w0} > 1$ and $R_{w} > 1$ in this section. Hence, from (4.2), the fast time variable $\bar{V}(U) = \frac{\omega}{h}$ in model (4.5) is a constant that does not depend on $U$. Thus, model (4.5) takes the following form:

The positivity and boundedness of model (5.1) has been established in Theorem 6 in Section 4.1. The basic reproduction number for model (5.1) is defined as

With the existence of equilibrium of model (5.1), the conclusion is given below.

Lemma 2. (i) If $R_{0} \leq 1$, then model (5.1) has only disease-free equilibrium $W_0 = (S_0, 0, 0, 0)$ with $S_0 = \frac{A}{\mu}$.

(ii) If $R_{0} > 1$, then model (5.1) has a unique endemic equilibrium $\hat{W} = (\hat{S}, \hat{E}, \hat{I}, \hat{U})$ with

The proof of Lemma (2) is simple. Hence, we omit it here. Based on the global stability of equilibria $W_{0}$ and $\hat{W}$ for model (5.1), the conclusion is given below.

Theorem 12. (i) If $R_{0}\leq1$, then disease-free equilibrium $W_{0}$ is globally asymptotically stable.

(ii) If $R_{0} > 1$, then endemic equilibrium $\hat{W}$ is globally asymptotically stable.

Proof. For conclusion $(i)$, we define the Lyapunov function $L_0$ as follows:

The derivative of $L_0(t)$ along with the solution $W_{0}$ of model (5.1) is given by

Clearly, if $R_0\leq1$, then $\frac{dL_{0}}{dt} < 0$ and $\frac{dL_{0}}{dt} = 0$ if and only if ($S(t), E(t), I(t), U(t)$) = $(S_0, 0, 0, 0)$. According to the Lyapunov theorem ^[39] and LaSalle invariable principle ^[33], $W_0$ is globally asymptotically stable.

For conclusion $(ii)$, we define the Lyapunov function $L_1$ as follows:

The derivative of $L_1(t)$ along with any solution $\hat{W}$ of model (5.1) is given by

From the equations in which equilibrium $\hat{W}$ is satisfied, we can easily obtain $\mu+\sigma+\zeta = \frac{\alpha{\hat{E}}}{\hat{I}}$ and $\frac{(\xi+\gamma)\hat{U}}{\theta{V}(1-\hat{U})} = \hat{I}$. Furthermore, we have $\frac{(\mu+\alpha)(\mu+\sigma+\zeta)\hat{I}}{\alpha}+\frac{(\xi+\gamma)\hat{U}(\mu+\alpha)(\mu+\sigma+\zeta)}{\alpha\theta\frac{\omega}{h}(1-\hat{U})} = 2(\mu+\alpha)\hat{E}$. Moreover, from $\xi+\gamma = \frac{\theta{I_1}V(1-\hat{U})}{\hat{U}}$, we also have $-\frac{(\mu+\alpha)(\mu+\sigma+\zeta)U}{\alpha\theta\frac{\omega}{h}(1-\hat{U})}(\xi+\gamma) = -(\mu+\alpha)\hat{E}\frac{U}{\hat{U}}$ and then $-\frac{(\mu+\alpha)(\mu+\sigma+\zeta)I}{\alpha}[-1+\frac{1-U}{1-\hat{U}}-\frac{(1-U)\hat{U}}{(1-\hat{U})U}]$ = $-\frac{(\mu+\alpha)(\mu+\sigma+\zeta)}{\alpha}I$ $\frac{(U-\hat{U})^2}{U(1-\hat{U})}-(\mu+\alpha)\hat{E}\frac{I\hat{U}}{\hat{I}U}$. Finally, we can obtain

Clearly, if $R_0 > 1$, then $\frac{\mathrm{d}L_1}{\mathrm{d} t}\leq0$ with $\frac{\mathrm{d}L_1}{\mathrm{d} t} = 0$ only at $\hat{W}$. According to the Lyapunov theorem ^[39] and LaSalle invariable principle ^[33], $\hat{W}$ is globally asymptotically stable. This completes the proof.

Let $\psi = \frac{h}{\omega}$ indicate the total quantity of B cells in the host. At this point, the basic reproduction number becomes the following:

We rewrite $R_0 = \frac{1}{(\mu+\sigma+\zeta)}\cdot\theta\cdot\frac{1}{\psi}\cdot\frac{1}{(\xi+\gamma)}\cdot\beta\cdot S_0\cdot\frac{1}{(\mu+\alpha)}\cdot\alpha$. Combining the explanation of basic reproduction number $R_{b0}$ in Subsection 4.2, $R_0$ decreases with the increase in the quantity of B cells in hosts, the decrease in the incubation period of latent individuals and the survival period of environmental viruses, the decrease in the infection rate of susceptible individuals infected with environmental viruses, and the decrease in the emission rate of infected individuals releasing the virus into the environment.

When $R_0 = 1$, we obtain the critical value of $\psi$ as follows:

Thus, according to Theorem 12, we obtain the following conclusions.

Corollary 1. If the number $\psi$ of B cells in the host satisfies $\psi > \psi_0$, then the disease between the hosts will be extinct. In contrast, if $\psi < \psi_0$, then the disease between the hosts will be persistent.

Remark 6. Corollary 1 shows that the increasing number of B cells in infected and susceptible individuals can effectively control the spread of diseases between the hosts. However, the main way to increase the number of B cells in susceptible and infected individuals is to carry out active and effective treatment and vaccination. Therefore, when a certain infectious disease appears in an area, actively carrying out large-scale effective vaccination, timely treatment of patients who have been infected with the disease, and effective treatment are important means of preventing the disease in a timely and effective manner.

6.   Numerical examples

In this section, the theoretical results obtained in Theorems 7, 8 and 10 are illustrated by the following numerical examples.

Example 1. In model (2.4), we choose the parameters $A = 3$, $\beta = 0.0351$, $\mu = 0.044$, $\gamma = 0.0151$, $\theta = 1.4\times10^{-6}$, $\Lambda_c = 5995$, $K = 1.51\times10^{-6}$, $m = 0.29$, $d = 0.16$, $p = 954$, $c = 61$, $\sigma = 0.05$, $\alpha = 0.004$, $\zeta = 0.032$, $\xi = 0.035$ and $\omega = 0.005$.

(a) and (b) satisfy the initial values $(S(k), E(k), I(k), U(k)) = (75+14*k, 50+15*k, 50+16*k, 0.08+17*k)$, where $k = 1, ..., 20$. (c) satisfies the initial values $(S(k), E(k), I(k)) = (175+8*k, 80+2*k, 60+2*k, 0.19+20*k)$, where $k = 1, ..., 20$.

By calculation, $R_{w0} = 1.0849 > 1$ and $R_{b0} = 0.7208 < 1$. Then, from Figure 1, the disease-free equilibrium $W_0 = (68.1818, 0, 0, 0)$ of model (4.5) is locally asymptotically stable, which means that Theorem 7 is true.

Example 2. In model (2.4), we choose the parameters $A = 30$, $\beta = 0.0351$, $\mu = 0.044$, $\gamma = 0.0151$, $\theta = 1.4\times10^{-6}$, $\Lambda_c = 5995$, $K = 1.51\times10^{-6}$, $m = 0.29$, $d = 0.16$, $p = 954$, $c = 61$, $\sigma = 0.05$, $\alpha = 0.004$, $\zeta = 0.032$, $\xi = 0.035$, $b4 = 1$ and $w = 4\times10^4$.

(a) and (b) satisfy the initial values $(S(k), E(k), I(k), U(k)) = (482+0.02*k, 280+0.07*k, 0.17+0.4*k, 0.29+0.5*k)$, where $k = 1, ..., 20$. (c) satisfies the initial values $(S(k), E(k), I(k)) = (412+0.9*k, 247+0.9*k, 70+0.9*k, 0.29+0.9*k)$, where $k = 1, ..., 20$.

By calculation, $R_{w0} = 1.0849 > 1$ and $R_{b0} = 7.2085 > 1$. Model (4.5) has a unique endemic equilibrium $\tilde{W} = (\tilde{S}, \tilde{E}, \tilde{I}, \tilde{U})$ = $(412.1898,247.1593, 7.8463, 0.82)$. Furthermore, $F'(\tilde{U}) = -20658 < 0$, $H_M = 14036.9165 > 0$ and $b_3b_2b_1-b_4b_1^2-b_3^2b_0 = 1.7112e-04 > 0$. Then, Figure 2 shows that the unique endemic equilibrium $\tilde{W} = (\tilde{S}, \tilde{E}, \tilde{I}, \tilde{U})$ is locally asymptotically stable, which means that the Theorem 8 is true.

Example 3. In model (2.4), we only need to choose parameters $\beta = 0.0357$, $\mu = 0.042$, $\theta = 1.5\times10^{-8}$, $d = 0.2$, $c = 30$ and $w = 4\times10^6$; the other parameters are the same as those given in Example 2.

(a) and (b) satisfy initial values $(S(k), E(k), I(k), U(k)) = (500+3*k, 100+3*k, 0.19+0.5*k, 0.23+0.5*k)$, where $k = 1, ..., 20$. (c) satisfies the initial values $(S(k), E(k), I(k)) = (680+5*k, 360+5*k, 150+1*k, 0.19+0.9*k)$, where $k = 1, ..., 20$.

By calculation, $R_{w0} = 2.0258 > 1$ and $R_{b0} = 1.0548 > 1$. Model (4.5) has a unique endemic equilibrium $\tilde{W} = (\tilde{S}, \tilde{E}, \tilde{I}, \tilde{U})$ = $(668.8068, 41.5242, 1.3395, 0.08)$. Furthermore, $F'(\tilde{U}) = 4434.3 > 0$, $H_M = 10232.2459 > 0$ and $b_3b_2b_1-b_4b_1^2-b_3^2b_0 = 3.5990e-06 > 0$. Then, Figure 3 shows that the unique endemic equilibrium $\tilde{W}$ is locally asymptotically stable, which means that Theorem 8 is true even if $F'(\tilde{U}) > 0$.

Example 4. In model (2.4), we choose the parameters $A = 4$, $\beta = 0.0358$, $\mu = 0.042$, $\gamma = 0.0151$, $\theta = 1.45\times10^{-6}$, $\Lambda_c = 5995$, $K = 1.51\times10^{-6}$, $m = 0.29$, $d = 0.16$, $p = 954$, $c = 61$, $\sigma = 0.05$, $\alpha = 0.004$, $\zeta = 0.062$, $\xi = 0.035$, $b4 = 1$ and $w = 4\times10^5$.

(b) and (c) satisfy the initial values $(S(k), E(k), I(k), U(k)) = (80+0.06*k, 10+0.07*k, 0.21+0.08*k, 0.21+0.09*k)$, where $k = 1, ..., 20$. (d) satisfies the initial values $(S(k), E(k), I(k)) = (80+1*k, 20+1*k, 10+2*k, 0.02+0.1*k)$ and $(S(k), E(k), I(k)) = (92+0.5*k, 2+0.5*k, 5+2*k, 0.04+0.1*k)$, $k = 1, ..., 20$, where the second set of initial values was taken near the unstable equilibrium $\tilde{W_1} = (\tilde{S_1}, \tilde{E_1}, \tilde{I_1}, \tilde{U_1})$.

By calculation, $R_{w0} = 1.0849 > 1$ and $R_{b0} = 0.9081 < 1$. Model (4.5) has two positive equilibria which are $\tilde{W_1} = (\tilde{S_1}, \tilde{E_1}, \tilde{I_1}, \tilde{U_1})$ = $(92.0980, 2.867116, 0.0745, 0.04)$ and $\tilde{W_2} = (\tilde{S_2}, \tilde{E_2}, \tilde{I_2}, \tilde{U_2})$ = $(79.0639, 14.7678, 0.3836, 0.24)$. Furthermore, $F'(\tilde{U_2}) = -370.0812 < 0$, $H_M = 926.4939 > 0$ and $b_3b_2b_1-b_4b_1^2-b_3^2b_0 = 6.7554e-06 > 0$. Then, Figure 4 shows that the endemic equilibrium $\tilde{W_1}$ is unstable and that $\tilde{W_2}$ is locally asymptotically stable, which means that Theorem 10 is true.

Example 5. In model (2.4), we only need to choose parameters $\beta = 0.0356$ and $\theta = 1.4\times10^{-6}$, and the other parameters are the same as the parameters given in Example 4.

(b) and (c) satisfy initial values $(S(k), E(k), I(k), U(k)) = (80+0.06*k, 7+0.07*k, 0.15+0.08*k, 0.12+0.09*k)$, where $k = 1, ..., 20$. (d) satisfies initial values $(S(k), E(k), I(k)) = (180+1*k, 100+1*k, 10+2*k, 0.2+0.1*k)$ and $(S(k), E(k), I(k)) = (89+0.5*k, 4+0.5*k, 5+2*k, 0.07+0.1*k)$, $k = 1, ..., 20$, where the second set of initial values was taken near the unstable equilibrium $\tilde{W_1} = (\tilde{S_1}, \tilde{E_1}, \tilde{I_1}, \tilde{U_1})$.

By calculation, $R_{w0} = 1.0849 > 1$ and $R_{b0} = 0.8817 < 1$. Model (4.5) has two positive equilibria which are $\tilde{W_1} = (\tilde{S_1}, \tilde{E_1}, \tilde{I_1}, \tilde{U_1})$ = $(89.8473, 4.9221, 0.1278, 0.07)$ and $\tilde{W_2} = (\tilde{S_2}, \tilde{E_2}, \tilde{I_2}, \tilde{U_2})$ = $(81.9001, 12.1782, 0.3163, 0.19)$. Furthermore, $F'(\tilde{U_2}) = 6093.1 > 0$, $H_M = 316.3866 > 0$ and $b_3b_2b_1-b_4b_1^2-b_3^2b_0 = 5.9487e-06 > 0$. Then, from Figure 5 shows that endemic equilibrium $\tilde{W_1}$ is unstable and that $\tilde{W_2}$ is locally asymptotically stable, which means that Theorem 10 is true even though $F'(\tilde{U_2}) > 0$.

7.   Conclusions

In this article, we investigate an SEIR epidemic model (2.4) coupling virus transmission in the body and vitro of hosts, incubation, humoral immunity and environmental effect, which are characterized and linked by two subsystems, i.e., (2.1) and (2.3).

With respect to the virus transmission process with humoral immunity in the body of hosts, we assume that the environmental contamination rate $U(t)$ in model (2.3) remains constant at $U\; (0\leq U \leq 1)$ since the spread of virus infection within hosts is much faster than that in vitro of hosts. The basic reproduction number $R_w$ with an antibody response is defined, by which antibody-free infection equilibrium $A_3(U)$ is globally asymptotically stable if $R_w\leq 1$, while infection equilibrium $A_4(U)$ is globally asymptotically stable if $R_w > 1$.

As is the SEIR model without an antibody response, we assume that the virus load in the body of hosts will tend to reach equilibrium, i.e., the fast time variable $V(s)$ in model (2.4) satisfies $V(s) = \bar{V}(U) = \frac{1}{c}(g(U)+pT^*_3(U))$ in Theorems 4. The basic reproduction number $R_{b0}$ is defined, from which we find that disease-free equilibrium $\tilde{W}_0$ is locally asymptotically stable if $R_{b0} < 1$, while the unique positive equilibrium $\tilde{W}$ is locally asymptotically stable if $R_{b0} > 1$ (or $R_{b0} = 1$ and $H_{M} > 0$) and $F'(\tilde{U})\leq 0$. When $R_{b0} < 1$ and $H_{M} > 0$, the system has two different positive equilibria $\tilde{W_1}$ and $\tilde{W_2}$ with $U_1 < U_2$, which means that system (4.5) experiences backward bifurcation at $R_{b0} = 1$. Meanwhile, $\tilde{W_2}$ is locally asymptotically stable if $F'(\tilde{U_2})\leq 0$ and $\tilde{W_1}$ are unstable.

Furthermore, for the SEIR model with an antibody response, we assume that the fast time variable $V(s)$ in model (2.4) satisfies $V(s) = \bar{V}(U) = \frac{\omega}{h}$ in Theorems 4. The basic reproduction number $R_{0}$ is defined, by which we find the disease-free equilibrium $W_0$ is globally asymptotically stable when $R_0\leq1$, while the unique endemic equilibrium $\hat{W}$ is globally asymptotically stable if $R_0 > 1$.

From the numerical examples we know that $F'(\tilde{U})\leq 0$ and $F'(\tilde{U_2})\leq0$ are pure mathematical conditions, and can only be used to prove the local stability of endemic equilibria $\tilde{W}$ and $\tilde{W_2}$. Generally, we hope that the local stability of model (4.5) can only be determined by the basic reproduction number $R_{b0}$, but we used some additional conditions to do that. So we have summarized several open problems. The first is whether condition $F'(\tilde{U})\leq 0$ and $F'(\tilde{U_2})\leq0$ can be taken off in the proof of Theorems 8 and 10. The second is whether we can get that (4.9) directly holds. The third one is whether we can further establish an appropriate Lyapunov function to obtain the global stability of $\tilde{W}$. We will continue to investigate these questions in the future.

The results obtained in this paper show that the strength of antibodies in hosts has a great effect on the spread of diseases between hosts. When the antibodies in hosts do not work or are weak, the results obtained in Section 4 show that backward bifurcation could occur. Even if the basic reproduction number $R_{b0}$ is less than 1, the disease will continue to spread, which will results in the control of the spread of disease between hosts being very difficult, making it difficult to effectively treat the disease.

When hosts have an extensive antibody response, the results obtained in Section 5 show that the spread of disease between hosts will be easy to control. Under antibody action, the basic reproduction number $R_0$ could be less than 1, and the diseases between hosts could be effectively controlled if we could decrease the incubation period of latent individuals, increase the production rate of B cells in hosts, and reduce the virus concentration (or load) both in the environment or in the body of hosts. Therefore, based on the above discussion with the explanation of the basic reproduction number $R_0$, we can take the following prevention and control measures.

1. Vaccination and antibody immunotherapy. After vaccination, susceptible individuals will widely produce antibodies, which could effectively prevent the virus from invading. Even if a susceptible individual is infected by the virus, due to the effect of the antibody response, the amount of virus in hosts can only reach $\bar{V}(U) = \frac{\omega}{h}$, so most of the infected individuals will only have some minor symptoms, and will not experience long-term infection. Furthermore, if antibody immunotherapy is performed on susceptible individuals, it will increase the number of antibody cells $\psi$ in infected hosts, and eventually make the basic reproduction number $R_0$ less than 1, then the spread of disease will be effectively controlled until the disease eliminated.

2. Timely strict isolation and treatment of infected individuals, isolation control of close human contacts. These measures could effectively prevent close contact between infected and susceptible individuals, and reduce the release of virus into the environment by infected individuals, further reducing the infection rate of susceptible individuals and the virus emission rate of infected individuals, and finally making the basic reproduction number $R_0$ less than 1, this could effectively control the disease transmission between hosts until its disappearance.

3. Legitimately expand the range of management, supervision and publicity, e.g., wearing masks (surgical masks, goggles, etc.), practicing safe social distancing in public or crowded places, frequent hand-washing, exercising of the body, maintain a healthy diet, and obtaining sufficient rest. These measures could effectively reduce the infection rate of susceptible individuals, so that the basic reproduction number $R_0$ will be less than 1, effectively controlling the spread of disease.

4. Manually eliminating the virus in the environment. These measures could not only effectively improve the clearance rate of environmental viruses but also reduce the average survival period of environmental viruses. Therefore, the basic reproduction number $R_0$ would be less than 1, the spread of disease would be effectively controlled, and the disease would eventually become extinct.

In this paper we proposed the SEIR model (2.4) coupling virus transmission both in body and vitro of hosts. Since the virus infection changes process in host is much faster than the disease transmission process between hosts, to simplify discussions, we assumed that the environmental contamination rate $U(t)$ in the fast time model (2.3) remains constant at $U$, and then, the fast time variable $V(s)$ in the slow time model (2.1) reaches its equilibrium $\bar{V}(U)$ while the state of disease transmission between the host is unchanged. Thus, we established the coupled model (4.1). Obviously, model (4.1) is a special limit state coupled model. In this paper, we mainly investigated the dynamical behavior of model (4.1). However, it is more realistic to directly investigate the coupling model (2.4) with both fast time $s$ and slow time $t$. Particularly, in model (2.4) we can assume that fast and slow times satisfy $t = \omega s$ is sufficiently, where $\omega$ is an enough small positive number. Therefore, an interesting open problem is to investigate coupling model (2.4) with assumption $t = \omega s$. We will discuss this problem in the future.

Acknowledgement

This authors were supported by the National Natural Science Foundation of China (11771373, 11861065, 12261087), the Open Project of Key Laboratory of Applied Mathematics of Xinjiang Province (2021D04014) and the Scientific Research Programmes of Colleges in Xinjiang (XJEDU20211002).

Conflict of interest

The authors declare that they have no conflicts of interest.