# Dynamics of an HBV infection model with cell-to-cell transmission and CTL immune response

Hajar Besbassi<fnref n="1" />, Zineb Elrhoubari<fnref n="1" />, Khalid Hattaf<fnref n="1" />,<fnref n="2" />, Noura Yousfi<fnref n="1" />

Emails: hajar.besbassi@gmail.com, zineb.elrhoubari@gmail.com, k.hattaf@yahoo.fr, nourayousfi@gmail.com

**RÉSUMÉ.** Dans ce travail, nous proposons un modèle mathématique pour décrire la dynamique de l'infection par le virus de l'hépatite B (VHB) en tenant compte la guérison des cellules infectées, l'exportation des cellules lymphocytes T cytotoxiques (LTC) précurseurs du thymus et les deux modes de transmission qui sont l'infection de virus-à-cellule et la transmission de cellule-à-cellule. La stabilité locale de l'équilibre libre et de l'équilibre d'infection chronique est obtenue par des équations caractéristiques. En outre, la stabilité globale des deux équilibres est établie par l'utilisation de deux techniques, la méthode directe de Lyapunov pour l'équilibre libre et l'approche géométrique pour l'équilibre de l'infection chronique.

**ABSTRACT.** In this work, we propose a mathematical model to describe the dynamics of the hepatitis B virus (HBV) infection by taking into account the cure of infected cells, the export of precursor cytotoxic T lymphocytes (CTL) cells from the thymus and both modes of transmission that are the virus-to-cell infection and the cell-to-cell transmission. The local stability of the disease-free equilibrium and the chronic infection equilibrium is obtained via characteristic equations. Furthermore, the global stability of both equilibria is established by using two techniques, the direct Lyapunov method for the disease-free equilibrium and the geometrical approach for the chronic infection equilibrium.

**MOTS-CLÉS :** Infection par le VHB, immunité, transmission cellule-à-cellule, stabilité.

**KEYWORDS :** HBV infection, immunity, cell-to-cell transmission, stability.

<footnote n="1">¹Laboratory of Analysis, Modeling and Simulation (LAMS), Faculty of Sciences Ben M'sik, Hassan II University, P.O Box 7955 Sidi Othman, Casablanca, Morocco</footnote>

<footnote n="2">²Centre Régional des Métiers de l'Education et de la Formation (CRMEF), 20340 Derb Ghalef, Casablanca, Morocco</footnote>
# 1. Introduction

HBV infection is a major global health problem that can cause acute or chronic infection and puts people at high risk of death from cirrhosis and liver cancer. In 2015, hepatitis B resulted in 887000 deaths, mostly from complications (including cirrhosis and hepatocellular carcinoma) [1]. In addition, cytotoxic T lymphocytes (CTL) cells play an important role in antiviral defense by killing the infected cells. On the other hand, HBV can spread by two fundamental modes, one by virus-to-cell infection through the extracellular space and the other by cell-to-cell transfer involving direct cell-to-cell contact [2, 3]. For these reasons, we propose the following model :

$$ egin{cases} \dot{x} = \lambda - dx - f(x, y, v)v - g(x, y)y + \rho y, \\ \dot{y} = f(x, y, v)v + g(x, y)y - (a + \rho)y - pyz, \\ \dot{v} = ky - \mu v, \\ \dot{z} = s + \frac{cyz}{\omega + y} - bz, \end{cases} \quad (1) $$

where $x(t)$, $y(t)$, $v(t)$ and $z(t)$ represent the concentrations of uninfected cells, infected cells, free virus and CTL cells at time $t$, respectively. Susceptible host (healthy hepatocytes) cells are produced at rate $\lambda$, die at rate $dx$ and become infected either by free virus at rate $f(x, y, v)v$ or by direct contact with an infected cell at rate $g(x, y)y$. Hence, the term $f(x, y, v)v+g(x, y)y$ denotes the total infection rate of uninfected cells. Infected cells cured at rate $\rho y$, die at rate $ay$ and are killed by the CTL immune response at rate $pyz$. Free virus is produced by an infected cell at rate $ky$ and decays at rate $\mu v$. CTL cells expand in response to viral antigen derived from infected cells at rate $\frac{cyz}{\omega + y}$, where $c$ is HBV-specific CTL stimulation rate and $\omega$ represents virus load for half-maximal CTL cells stimulation [4] and decay in the absence of antigenic stimulation at the rate $bz$. The parameter $s$ denotes the export of precursor CTL cells from the thymus [4].

As in [5], the incidence functions $f(x, y, v)$ and $g(x, y)$ for the two modes are continuously differentiable and satisfy the following hypotheses :

(H$_0$) $g(0,y)=0$, for all $y \ge 0$; $\frac{\partial g}{\partial x}(x, y) \ge 0$ (or $g(x, y)$ is a strictly monotone increasing function with respect to $x$ when $f \equiv 0$) and $\frac{\partial g}{\partial y}(x, y) \le 0$, for all $x \ge 0$ and $y \ge 0$.

(H$_1$) $f(0, y, v) = 0$, for all $y \ge 0$ and $v \ge 0$,

(H$_2$) $f(x, y, v)$ is a strictly monotone increasing function with respect to $x$ (or $\frac{\partial f}{\partial x}(x, v, y) \ge 0$ when $g(x, y)$ is a strictly monotone increasing function with respect to $x$), for any fixed $y \ge 0$ and $v \ge 0$,

(H$_3$) $f(x, y, v)$ is a monotone decreasing function with respect to $y$ and $v$.

Biologically, the four hypotheses are reasonable and consistent with the reality. For more details on the biological significance of these four hypotheses, we refer the reader to the works [5–7]. Further, the general incidence functions $f(x, y, v)$ and $g(x, y)$ include various types of incidence rates existing in the literature.
The main objective of this work is to investigate the dynamical behavior of system (1). To do this end, we start with the existence, the positivity and boundedness of solutions, which implies that our model is well-posed. After, we determine the basic reproduction number and steady states of the model. Finally, the local and global stabilities of the disease-free equilibrium and the chronic infection equilibrium are established.

## 2. Positivity and boundedness of solutions

The first important step is to validate our model as model that represents the evolution of cells and virus. In the following result, we show that all cell and virus concentrations are non-negative and bounded.

**Theorem 2.1.** All solutions of system (1) starting from positive initial value $(x_0, y_0, v_0, z_0)$ remain bounded and positive for all $t > 0$. Moreover, we have

$$ (i) \quad T(t) \le T(0) + \frac{\lambda}{\delta}, $$

$$ (ii) \quad v(t) \le v_0 + \frac{k}{\mu} \|y\|_{\infty}, $$

$$ (iii) \quad z(t) \le z_0 + \frac{s}{b} + \frac{c}{p\omega} \left[ \frac{\lambda}{b} + x_0 + y_0 + \alpha_1 \|x\|_{\infty} + \alpha_2 \|y\|_{\infty} \right], $$

where $T = x + y$ that represents the total cells of liver, $\delta = \min\{a, d\}$, $\alpha_1 = \max\{0, 1 - \frac{d}{b}\}$ and $\alpha_2 = \max\{0, 1 - \frac{a}{b}\}$.

**Proof.** First, we prove that any solution starting in the first quadrant $\mathbb{R}_+^4 = \{(x, y, v, z) \in \mathbb{R}_+^4 : x \ge 0, y \ge 0, v \ge 0, z \ge 0\}$ stays in $\mathbb{R}_+^4$. In fact, for $(x(t), y(t), v(t), z(t) \in \mathbb{R}_+^4)$, we have

$$ \begin{aligned} \dot{x}|_{x=0} &= \lambda + \rho y \ge 0, & \dot{y}|_{y=0} &= f(x, 0, v)v \ge 0, \\ \dot{v}|_{v=0} &= ky \ge 0, & \dot{z}|_{z=0} &= s > 0, \end{aligned} \quad (2) $$

this immediately implies that all solutions of system (1) with initial condition $(x_0, y_0, v_0, z_0) \in \mathbb{R}_+^4$ stay in the first quadrant.

Next, we prove the boundedness of solutions. We have that $\dot{T} \le \lambda - \delta T$. Then

$$ T(t) \le T(0)e^{-\delta t} + \frac{\lambda}{\delta}(1 - e^{-\delta t}). \quad (3) $$

Since $0 \le e^{-\delta t} \le 1$ and $1 - e^{-\delta t} \le 1$, we deduce (i).

From the third equation of (1), we have

$$ \dot{v} = ky - \mu v, $$

which implies that

$$ v(t) = v(0)e^{-t\mu} + k \int_{0}^{t} y(\theta)e^{(\theta-t)\mu} d\theta. \quad (4) $$

Thus,

$$ v(t) \le v(0) + \frac{k}{\mu} \|y\|_{\infty} (1 - e^{-t\mu}). \quad (5) $$
Since $1 - e^{-t\mu} \le 1$, we have (ii).

From the fourth equation of (1), we get

$$ \dot{z} + bz \le s + \frac{c}{w}yz, \qquad (6) $$

Hence,

$$ \dot{z} + bz \le s + \frac{c}{pw}[\lambda - (\dot{x} + dx) - (\dot{y} + ay)], \qquad (7) $$

Then

$$ \begin{aligned} z(t)e^{bt} - z_0 & \le \left(\frac{s}{b} + \frac{c\lambda}{bp\omega}\right)(e^{bt} - 1) - \frac{c}{p\omega} \int_0^t e^{(b-d)\theta} \frac{d}{d\theta}(x(\theta)e^{d\theta})d\theta \\ & \quad - \frac{c}{p\omega} \int_0^t e^{(b-a)\theta} \frac{d}{d\theta}(y(\theta)e^{a\theta})d\theta. \end{aligned} \qquad (8) $$

By integration by parts, we have

$$ \begin{aligned} \int_0^t e^{(b-d)\theta} \frac{d}{d\theta}(x(\theta)e^{d\theta})d\theta &= [x(\theta)e^{b\theta}]_0^t - (b-d) \int_0^t x(\theta)e^{b\theta}d\theta, \\ \int_0^t e^{(b-a)\theta} \frac{d}{d\theta}(y(\theta)e^{a\theta})d\theta &= [y(\theta)e^{b\theta}]_0^t - (b-a) \int_0^t y(\theta)e^{b\theta}d\theta. \end{aligned} \qquad (9) $$

Thus,

$$ \begin{aligned} z(t) & \le \left[ \frac{c}{p\omega}(x_0 + y_0) + z_0 \right] e^{-bt} + \left( \frac{s}{b} + \frac{c\lambda}{bp\omega} \right) (1 - e^{-bt}) \\ & \quad + \frac{c}{p\omega} \left[ \int_0^t [(b-d)x(\theta) + (b-a)y(\theta)]e^{b(\theta-t)}d\theta - x(t) - y(t) \right]. \end{aligned} \qquad (10) $$

In order to find an upper bound to this integral, we study the following cases :
If $b - d \le 0$ and $b - a \le 0$, then

$$ z(t) \le z_0 + \frac{s}{b} + \frac{c}{p\omega} \left( \frac{\lambda}{b} + x_0 + y_0 \right). \qquad (11) $$

If $b - d \le 0$ and $b - a \ge 0$, then

$$ z(t) \le z_0 + \frac{s}{b} + \frac{c}{p\omega} \left[ \frac{\lambda}{b} + x_0 + y_0 + \left(1 - \frac{a}{b}\right) \|y\|_{\infty} \right]. \qquad (12) $$

If $b - d \ge 0$ and $b - a \le 0$, then

$$ z(t) \le z_0 + \frac{s}{b} + \frac{c}{p\omega} \left[ \frac{\lambda}{b} + x_0 + y_0 + \left(1 - \frac{d}{b}\right) \|x\|_{\infty} \right]. \qquad (13) $$

If $b - d \ge 0$ and $b - a \ge 0$, then

$$ z(t) \le z_0 + \frac{s}{b} + \frac{c}{p\omega} \left[ \frac{\lambda}{b} + x_0 + y_0 + \left(1 - \frac{d}{b}\right) \|x\|_{\infty} + \left(1 - \frac{a}{b}\right) \|y\|_{\infty} \right]. \qquad (14) $$

From (11)-(14), we can conclude (iii). ■
### 3. Basic reproduction number and Equilibria

Obviously, system (1) has always one infection free equilibrium $E_f(\frac{\lambda}{d}, 0, 0, \frac{s}{b})$. Then we define the basic reproduction number of (1) as follows:

$$R_0 = \frac{kf(\frac{\lambda}{d}, 0, 0) + \mu g(\frac{\lambda}{d}, 0)}{\mu(a + \rho + p\frac{s}{b})}, \quad (15)$$

which can be rewritten as $R_0 = R_{01} + R_{02}$, where

$$R_{01} = \frac{k}{a + \rho + p\frac{s}{b}} \times f(\frac{\lambda}{d}, 0, 0) \times \frac{1}{\mu},$$

and

$$R_{02} = g(\frac{\lambda}{d}, 0) \times \frac{1}{a + \rho + p\frac{s}{b}}.$$

In the formula (15), $\frac{1}{a + \rho + p\frac{s}{b}}$ denotes the average life expectancy of infected cells, which is less than $\frac{1}{a}$ because of the role of immune cells; $\frac{k}{a + \rho + p\frac{s}{b}}$ denotes the amount of virus generated from an infected during its survival period; $\frac{1}{\mu}$ is the average life expectancy of viruses; $\frac{\lambda}{d}$ denotes the number of susceptible cells at the beginning of the infectious process, which means that $f(\frac{\lambda}{d}, 0, 0)$ and $g(\frac{\lambda}{d}, 0)$ are the values of both incidence functions when all cells are uninfected. Hence, $R_{01}$ is the basic reproduction number corresponding to virus-to-cell infection mode, whereas $R_{02}$ is the basic reproduction number corresponding to cell-to-cell transmission mode. Therefore, $R_0$ describes the average number of newly infected cells generated from one infected cell at the beginning of the infectious process.

To find the other equilibrium of (1), we solve the following system

$$\lambda - dx - f(x, y, v)v - g(x, y)y + \rho y = 0, \quad (16)$$

$$f(x, y, v)v + g(x, y)y - (a + \rho)y - pyz = 0, \quad (17)$$

$$ky - \mu v = 0, \quad (18)$$

$$s + \frac{cyz}{\omega + y} - bz = 0. \quad (19)$$

By (16) to (19), we have $v = \frac{ky}{\mu}$, $z = \frac{s(\omega + y)}{b\omega + (b - c)y}$ and $x = \psi(y)$, where

$$\psi(y) = \frac{\lambda}{d} - y \left( \frac{a}{d} + \frac{ps(\omega + y)}{db\omega + dy(b - c)} \right). \quad (20)$$

Define a function $\Gamma$ on $[0, +\infty) - \{\frac{b\omega}{b-c}\}$ as follows

$$\Gamma(y) = kf\left(\psi(y), y, \frac{ky}{\mu}\right) + \mu g\left(\psi(y), y\right) - \mu \left(a + \rho + \frac{ps(\omega + y)}{b\omega + y(b - c)}\right). \quad (21)$$
When $R_0 > 1$, we have $\Gamma(0) = \mu(a + \rho + p\frac{s}{b})(R_0 - 1) > 0$ and

$$ \psi'(y) = -\frac{a}{d} - \frac{ps(\omega + y)}{b\omega + y(b-c)} - \frac{psd\omega cy}{[db\omega + dy(b-c)]^2} < 0, $$

$$ \Gamma'(y) = k \left( \psi'(y) \frac{\partial f}{\partial x} + \frac{\partial f}{\partial y} + \frac{k}{\mu} \frac{\partial f}{\partial v} \right) + \mu \left( \psi'(y) \frac{\partial g}{\partial x} + \frac{\partial g}{\partial y} \right) - \frac{\mu c \rho s \omega}{[b\omega + y(b-c)]^2} < 0. \quad (22) $$

Let $\alpha = \frac{b\omega}{c-b}$ be a pole of $\Gamma$. Hence, we discuss two cases :

(i) If $c > b$, then $\alpha > 0$. Since $z = \frac{s(\omega+y)}{b\omega+y(b-c)} \ge 0$, we have that $y < \alpha$, which means that there is no positive equilibrium point if $y \ge \alpha$. It is evident to show that

$$ \lim_{y \to \alpha^{-}} \Gamma(y) = -\infty. \quad (23) $$

Thus, there exists a unique $y^* \in (0, \alpha)$ such that $\Gamma(y^*) = 0$. Since $\psi(0) = \frac{\lambda}{d} > 0$ and $\lim_{y \to \alpha^{-}} \psi(y) = -\infty$. Then there exists a unique $\bar{y} \in (0, \alpha)$ such that $\psi(\bar{y}) = 0$. So, we have

$$ \Gamma(\bar{y}) = -\mu(a + \rho + \frac{ps(\omega + \bar{y})}{b\omega + \bar{y}(b-c)}) < 0. $$

Then we deduce that $0 < y^* < \bar{y}$ implying $\psi(\bar{y}) < \psi(y^*) < \psi(0)$. Thus, $0 < x^* < \frac{\lambda}{d}$. Also, it is clear that $v^*$ and $z^*$ are positive. Thus, model (1) has a unique chronic infection equilibrium $E^*(x^*, y^*, v^*, z^*)$, where $x^* \in (0, \frac{\lambda}{d})$, $y^* \in (0, \bar{y})$, $v^* > 0$ and $z^* > 0$.

(ii) If $c < b$, then $\alpha < 0$ and $\lim_{y \to +\infty} \psi(y) = -\infty$. As $\psi(0) = \lambda/d > 0$, then there exists a unique $\tilde{y} \in (0, +\infty)$ such that $\psi(\tilde{y}) = 0$. We have

$$ \Gamma(\tilde{y}) = -\mu(a + \rho + \frac{ps(\omega + \tilde{y})}{b\omega + \tilde{y}(b-c)}) < 0. \quad (24) $$

Similarly to above, we can show that $x^*, v^*$ and $z^*$ are positive. Therefore, model (1) has a unique chronic infection equilibrium $E^*(x^*, y^*, v^*, z^*)$, where $x^* > 0$, $y^* \in (0, \bar{y})$, $v^* > 0$ and $z^* > 0$.

The previous discussions can be summarized in the following result.

**Theorem 3.1.**

(i) When $R_0 \le 1$, the model (1) has always one infection-free equilibrium of the form $E_f(\frac{\lambda}{d}, 0, 0, \frac{s}{b})$.

(ii) When $R_0 > 1$, the model (1) has a unique chronic infection equilibrium of the form $E^*(x^*, y^*, v^*, z^*)$ with $x^* > 0$, $y^* > 0$, $v^* > 0$, and $z^* > 0$.
## 4. Local stability of equilibria

In this section, we discuss the local stability of both equilibria of model (1). Note that the Jacobian matrix of (1) is given by

$$ 
\begin{pmatrix}
-d - v \frac{\partial f}{\partial x} - y \frac{\partial g}{\partial x} & -v \frac{\partial f}{\partial y} - y \frac{\partial g}{\partial y} - g(x, y) + \rho & -v \frac{\partial f}{\partial v} - f(x, y, v) & 0 \\
v \frac{\partial f}{\partial x} + y \frac{\partial g}{\partial x} & v \frac{\partial f}{\partial y} + y \frac{\partial g}{\partial y} + g(x, y) - (a + \rho + p z) & v \frac{\partial f}{\partial v} + f(x, y, v) & -\rho y \\
0 & \frac{k}{cz\omega} & -\mu & 0 \\
0 & \frac{cs}{b\omega} & 0 & \frac{cy}{\omega + y} - b
\end{pmatrix}.
\quad (25) $$

Firstly, we get the following result.

**Theorem 4.1.** The infection-free equilibrium $E_f$ is locally asymptotically stable if $R_0 < 1$ and becomes unstable if $R_0 > 1$.

**Proof.** Evaluated (25) at $E_f$, we obtain

$$ J_{E_f} = \begin{pmatrix} -d & -g(\frac{\lambda}{d}, 0) + \rho & -f(\frac{\lambda}{d}, 0, 0) & 0 \\ 0 & g(\frac{\lambda}{d}, 0) - (a + \rho + p \frac{s}{b}) & f(\frac{\lambda}{d}, 0, 0) & 0 \\ 0 & k & -\mu & 0 \\ 0 & \frac{cs}{b\omega} & 0 & -b \end{pmatrix}. \quad (26) $$

Then the characteristic equation at $E_f$ is given by

$$ (\xi+d)(\xi+b) \left[ \xi^2 + \left( \mu + a + \rho + p \frac{s}{b} - g(\frac{\lambda}{d}, 0) \right) \xi + \mu \left( a + \rho + p \frac{s}{b} \right) (1-R_0) \right] = 0. \quad (27) $$

Hence, the roots of (27) are :

$$ \xi_1 = -d, $$

$$ \xi_2 = -b, $$

$$ \xi_3 = \frac{g(\frac{\lambda}{d}, 0) - \mu - a - \rho - p \frac{s}{b} - \sqrt{\Delta}}{2}, \quad (28) $$

$$ \xi_4 = \frac{g(\frac{\lambda}{d}, 0) - \mu - a - \rho - p \frac{s}{b} + \sqrt{\Delta}}{2}, $$

where $\Delta = [\mu + a + \rho + p \frac{s}{b} - g(\frac{\lambda}{d}, 0)]^2 - 4\mu(a + \rho + p \frac{s}{b})(1 - R_0)$. Clearly, $\xi_1$, $\xi_2$ and $\xi_3$ are negative. However, $\xi_4$ is negative if $R_0 < 1$ and is positive if $R_0 > 1$. Therefore, $E_f$ is locally asymptotically stable if $R_0 < 1$ and unstable if $R_0 > 1$. ■

Next, we focus on the local stability of the chronic infection equilibrium $E^*$.

**Theorem 4.2.** The chronic infection equilibrium $E^*$ is locally asymptotically stable if $R_0 > 1$.

**Proof.** Suppose that $R_0 > 1$. Evaluating (25) at $E^*$ gives

$$ \xi^4 + a_1 \xi^3 + a_2 \xi^2 + a_3 \xi + a_4 = 0, \quad (29) $$
where

$$a_1 = d + \mu + a + \rho + p z^* + \frac{s}{z^*} + C_1 - C_2,$$

$$a_2 = d(a + \rho + p z^* - C_2) + \left(\mu + \frac{s}{z^*}\right) \left(d - a - \rho - p z^* + C_1 - C_2\right) - (a + 2\rho + p z^*) C_1 - k C_3 \\ - \mu \left(\frac{c y^*}{\omega + y^*} - b\right) + \frac{c \rho \omega y^* z^*}{(\omega + y^*)^2},$$

$$a_3 = \left(\mu + b - \frac{cy^*}{\omega + y^*}\right) \left(d(a + \rho + p z^*) + C_2 + (a + p z^*) C_1\right) - k \left(d + b - \frac{cy^*}{\omega + y^*}\right) C_3 \\ + \mu \frac{s}{z^*} \left(d + a + \rho + p z^* + C_1 - C_2\right) + \frac{c \rho \omega y^* z^*}{(\omega + y^*)^2} (d + \mu + C_1),$$

$$a_4 = d \left( \frac{cy^*}{\omega + y^*} - b \right) \left( kC_3 + \mu C_2 - \mu(a + \rho + p z^*) \right) + \mu C_1 \left( (a + p z^*) \frac{s}{z^*} + \frac{c \rho \omega y^* z^*}{(\omega + y^*)^2} \right),$$

with

$$C_1 = v^* \frac{\partial f}{\partial x}(x^*, y^*, v^*) + y^* \frac{\partial g}{\partial x}(x^*, y^*),$$

$$C_2 = v^* \frac{\partial f}{\partial y}(x^*, y^*, v^*) + y^* \frac{\partial g}{\partial y}(x^*, y^*) + g(x^*, y^*),$$

$$C_3 = v^* \frac{\partial f}{\partial v}(x^*, y^*, v^*) + f(x^*, y^*, v^*).$$

When $R_0 > 1$, it is clear that $a_1, a_2, a_3$ and $a_4$ are positive. In addition,

$$\begin{vmatrix} a_1 & 1 \\ a_3 & a_2 \end{vmatrix} = a_1 a_2 - a_3 > 0. \qquad (30)$$

In the same way, we have

$$\begin{vmatrix} a_1 & 1 & 0 \\ a_3 & a_2 & a_1 \\ 0 & a_4 & a_1 \end{vmatrix} = a_3 \begin{vmatrix} a_1 & 1 \\ a_2 & a_3 \end{vmatrix} - a_1^2 a_4 > 0. \qquad (31)$$

From the Routh-Hurwitz theorem, we know that all roots of (29) have negative real parts. Thus, the chronic infection equilibrium $E^*$ is locally asymptotically stable for $R_0 > 1$. ■

## 5. Global Stability of Equilibria

In this section, we study the global stability of both equilibria. For the infection-free equilibrium $E_f$, we assume that $a \ge d$. Therefore, we have the following result.

**Theorem 5.1.** The infection-free equilibrium $E_f$ is globally asymptotically stable if $R_0 \le 1$.
**Proof.** Consider

$$ \Omega = \{ (x, y, v, z) \in \mathbb{R}_+^4 : z \ge \frac{s}{b} \}. $$

We see that any solution $(x(t), y(t), v(t), z(t))$ starting in $\Omega$ remains there forever. In fact, by Theorem 2.1 we have $(x(t), y(t), v(t), z(t)) \in \mathbb{R}_+^4$. it remains to prove that $z \ge \frac{s}{b}$ with $z_0 \ge \frac{s}{b}$. According to the fourth equation of (1), we have

$$ z(t) \ge \frac{s}{b} + \left(z_0 - \frac{s}{b}\right) e^{-bt}, \qquad (32) $$

which implies that $z \ge \frac{s}{b}$. Then $(x(t), y(t), v(t), z(t)) \in \Omega$. We construct the Lyapunov functional $L$ on $\Omega$ as follows:

$$ L(t) = y(t) + \frac{f(\frac{\lambda}{d}, 0, 0)}{\mu} v(t). $$

Calculating the time derivative of $L$ along the positive solution of (1), we get

$$ \frac{dL}{dt} = \left( f(x, y, v) - f(\frac{\lambda}{d}, 0, 0) \right) v + \left( a + \rho + p z \right) \left( \frac{k f(\frac{\lambda}{d}, 0, 0) + \mu g(x, y)}{\mu(a + \rho + p z)} - 1 \right) y. $$

It is not hard to see that $\limsup_{t \to \infty} x(t) \le \frac{\lambda}{d}$ and $\limsup_{t \to \infty} z(t) \ge \frac{s}{b}$. This yields that all omega limit points satisfy $x(t) \le \frac{\lambda}{d}$ and $z(t) \ge \frac{s}{b}$. So, it suffices to consider solutions for which $x(t) \le \frac{\lambda}{d}$ and $z(t) \ge \frac{s}{b}$. Using the expression of $R_0$ given in (15), we obtain

$$ \begin{aligned} \frac{dL}{dt} &\le \left( f(x, 0, 0) - f(\frac{\lambda}{d}, 0, 0) \right) v + (a + \rho + p z)(R_0 - 1)y \\ &\le (a + \rho + p z)(R_0 - 1)y. \end{aligned} $$

Since $R_0 \le 1$, we have $\frac{dL}{dt} \le 0$. In addition, it is easy to show that the largest compact invariant set in $\{(x, y, v, z) : \frac{dL}{dt}(t) = 0\}$ is the singleton $\{E_f\}$. By the LaSalle invariance principle, the infection-free equilibrium $E_f$ is globally asymptotically stable for $R_0 \le 1$. ■

Next, we will study the global dynamics of model (1) when $R_0 > 1$. Firstly, we need the following lemma.

**Lemma 5.2.** The model (1) is uniformly persistent if $R_0 > 1$.

**Proof.** The maximal invariant set $M$ on the boundary $\partial\Omega$ is the singleton $\{E_f\}$ and it is isolated. Further, Theorem 4.2 given in [8] ensures the equivalence between the uniform persistence of model (1) and the instability of the infection-free equilibrium $E_f$. From Theorem 4.1, we have $E_f$ is unstable if $R_0 > 1$. Therefore, model (1) is uniformly persistent if $R_0 > 1$. ■

Now, we focus on the global stability of $E^*$ under the assumption $R_0 > 1$ and the incidence function $f$ satisfies the following property:
$(H_4) \quad f(x, y, v) + v \frac{\partial f}{\partial v}(x, y, v) \ge 0 \text{ for all } x \ge 0, y \ge 0 \text{ and } v \ge 0.$

**Theorem 5.3.** Suppose $R_0 > 1$ and $(H_4)$ holds. Then the chronic infection equilibrium $E^*$ is globally asymptotically stable.

**Proof.** To study the global stability of $E^*$, we apply the geometrical approach given in [9]. So, we consider the following sub-system:

$$ egin{cases} \dot{x} = \lambda - dx - f(x, y, v)v - g(x, y)y + \rho y, \\ \dot{y} = f(x, y, v)v + g(x, y)y - (a + \rho)y - pyz, \\ \dot{v} = ky - \mu v. \end{cases} \tag{33} $$

The Jacobian matrix of system (33) is:

$$ J = \begin{pmatrix} -d - v \frac{\partial f}{\partial x} - y \frac{\partial g}{\partial x} & -v \frac{\partial f}{\partial y} - y \frac{\partial g}{\partial y} - g(x, y) + \rho & -v \frac{\partial f}{\partial v} - f(x, y, v) \\ v \frac{\partial f}{\partial x} + y \frac{\partial g}{\partial x} & v \frac{\partial f}{\partial y} + y \frac{\partial g}{\partial y} + g(x, y) - (a + \rho + vz) & v \frac{\partial f}{\partial v} + f(x, y, v) \\ 0 & k & -\mu \end{pmatrix} \tag{34} $$

and its second additive compound matrix is:

$$ J^{[2]} = \begin{pmatrix} j_{11} + j_{22} & j_{23} & -j_{13} \\ j_{32} & j_{11} + j_{33} & j_{12} \\ -j_{31} & j_{21} & j_{22} + j_{33} \end{pmatrix}, \tag{35} $$

where $j_{kl}$ is the $(k,l)$th entry of the matrix $J$. In this case, we choose $P = \text{diag}(1, \frac{y}{v}, \frac{y}{v})$. Then:

$$ P_f P^{-1} = \text{diag}(0, \frac{\dot{y}}{y} - \frac{\dot{v}}{v}, \frac{\dot{y}}{y} - \frac{\dot{v}}{v}), $$

where $P_f$ is obtained by replacing each entry $p_{ij}$ of $P$ by its derivative in the direction of solution of (33). Also, we have:

$$ B = P_f P^{-1} + P J^{[2]} P^{-1} = \begin{pmatrix} B_{11} & B_{12} \\ B_{21} & B_{22} \end{pmatrix}, \tag{36} $$

where:

$$ B_{11} = -(a + d + \rho + vz) - \frac{\partial f}{\partial x}v - \frac{\partial g}{\partial x}y + \frac{\partial f}{\partial y}v + \frac{\partial g}{\partial y}y + g(x, y), $$

$$ B_{12} = \left( \frac{v}{y} \left( \frac{\partial f}{\partial v} v + f \right) \quad \frac{v}{y} \left( \frac{\partial f}{\partial v} v + f \right) \right), $$

$$ B_{21} = \left( \frac{ky}{v} \right), $$

$$ B_{22} = \begin{pmatrix} \frac{\dot{y}}{y} - \frac{\dot{v}}{v} - \mu - d - \frac{\partial f}{\partial x}v - \frac{\partial g}{\partial x}y & \rho - \frac{\partial f}{\partial y}v - g(x, y) - \frac{\partial g}{\partial y}y, \\ \frac{\partial f}{\partial x}v + \frac{\partial g}{\partial x}y & \frac{\dot{y}}{y} - \frac{\dot{v}}{v} - a - \rho - vz - \mu + \frac{\partial f}{\partial y}v + g(x, y) + \frac{\partial g}{\partial y}y \end{pmatrix}. $$

We choose a norm in $\mathbb{R}^3$ as follows $|\omega_1, \omega_2, \omega_3| = \max\{|\omega_1|, |\omega_2| + |\omega_3|\}$ for $(\omega_1, \omega_2, \omega_3) \in \mathbb{R}^3$. Thus, the Lozinskii measure $\mu$ with respect to this norm $|\cdot|$ can be estimated as follows (see [10]):

$$ \mu(B) \le \sup\{g_1, g_2\}, \tag{37} $$
where $g_1 = \mu_1(B_{11}) + |B_{12}|)$ and $g_2 = |B_{21}| + \mu_1(B_{22})$.

Here, $\mu_1$ represents the Lozinskii measure with respect to $l_1$ vector norm, $|B_{12}|$ and $|B_{21}|$ are matrix norms with respect to $l_1$ norm. Furthermore, we have

$$ \mu_1(B_{11}) = -(a+d+\rho+pz) - \frac{\partial f}{\partial x}v - \frac{\partial g}{\partial x}y + \frac{\partial f}{\partial y}v + \frac{\partial g}{\partial y}y + g(x, y), $$

$$ \begin{aligned} |B_{12}| &= \frac{v}{y} \left( v \frac{\partial f}{\partial v} + f(x, y, v) \right) \\ &= \frac{\dot{y}}{y} + a + \rho + p z + \frac{v^2}{y} \frac{\partial f}{\partial v} - g(x, y), \end{aligned} $$

$$ |B_{21}| = k \frac{\dot{y}}{y} = \frac{\dot{v}}{v} + \mu, $$

$$ \begin{aligned} \mu_1(B_{22}) &= \max \left\{ \frac{\dot{y}}{y} - \frac{\dot{v}}{v} - \mu - d, \frac{\dot{y}}{y} - \frac{\dot{v}}{v} - \mu - a - p z \right\} \\ &\leq \frac{\dot{y}}{y} - \frac{\dot{v}}{v} - \mu - \delta. \end{aligned} $$

Thus,

$$ \begin{aligned} g_1 &= \frac{\dot{y}}{y} - d + \frac{v^2}{y} \frac{\partial f}{\partial v} - \frac{\partial f}{\partial x}v - \frac{\partial g}{\partial x}y + \frac{\partial f}{\partial y}v + \frac{\partial g}{\partial y}y \\ &\leq \frac{\dot{y}}{y} - \delta, \end{aligned} $$

and

$$ g_2 \leq \frac{\dot{y}}{y} - \delta. $$

Therefore,

$$ \mu(B) \leq \frac{\dot{y}}{y} - \delta. $$

From Lemma 5.2, the model (1) is uniformly persistent for $R_0 > 1$. Then there exists a compact absorbing set $K \subset \Omega$ [11]. Along each solution $(x(t), y(t), v(t))$ of (1) with $X_0 = (x(0), y(0), v(0)) \in K$, we have

$$ \frac{1}{t} \int_{0}^{t} \mu(B) ds \leq \frac{1}{t} \ln \left( \frac{y(t)}{y(0)} \right) - \delta, $$

which implies that

$$ \overline{q_2} = \limsup_{t \to \infty} \sup_{X_0 \in K} \frac{1}{t} \int_0^t \mu(B) ds \leq \frac{-\delta}{2} < 0. $$

Hence, the positive equilibrium $(x^*, y^*, v^*)$ of the sub-system (33) is globally asymptotically stable. From the last equation of model(1), we have

$$ \dot{z} = s + \frac{cyz}{\omega + y} - bz, \quad (38) $$

and its limit system is

$$ \dot{z} = s + \frac{cy^*z}{\omega + y^*} - bz. \quad (39) $$
By $b - \frac{cy^*}{\omega + y^*} = \frac{s}{z^*}$, we obtain

$$ \dot{z} = s \left(1 - \frac{z}{z^*}\right). \qquad (40) $$

Then

$$ \lim_{t \to \infty} z(t) = z^*. \qquad (41) $$

Therefore, the chronic infection equilibrium $E^*$ is globally asymptotically stable for $R_0 > 1$. ■

# 6. Discussion and conclusions

In this work, we have proposed an HBV infection model that takes into the account cell-to-cell transmission and CTL immune response. In the proposed model, the infection processes for the two modes of transmission are modeled by two general functions. We first proved the existence, positivity, and the boundedness of solutions of the problem which ensures that our model is well-posed. Under some assumptions about the general incidence functions, the global dynamics of the model are fully characterized by a threshold parameter called the basic reproduction number $R_0$. More precisely, the infection-free equilibrium $E_f$ is globally asymptotically stable if $R_0 \le 1$ which biologically means that the virus is cleared and the infection die out. When $R_0 > 1$, $E_f$ becomes unstable and the chronic infection equilibrium $E^*$ is globally asymptotically stable. In this case, the HBV persists in the liver. From the above analytical results, we deduce a strategy to control the HBV infection. This strategy is based on the reduction of $R_0$ and makes its value less then or equal to 1. From explicit expression of $R_0$ given in (15), the value of $R_0$ can be reduced by increasing the export of the thymus. Therefore, we conclude that cellular immunity mediated by CTL cells plays an important role in the clearance of HBV from the liver.

# Acknowledgments

We would like to thank the editor and the anonymous reviewers for their comments and suggestions that greatly improved the quality of the manuscript.

# 7. Bibliographie

[1] WHO, Hepatitis B, July 2018, Available at : https://www.who.int/en/news-room/fact-sheets/detail/hepatitis-b.

[2] W. MOTHES, N. M. SHERER, J. JIN, P. ZHONG, Virus cell-to-cell transmission, *Journal of Virology*, vol. 84, n° 17, 2010, 8360–8368.

[3] P. ZHONG, L. M. AGOSTO, J. B. MUNRO, W. MOTHES, Cell-to-cell transmission of viruses, *Current Opinion in Virology*, vol. 3, n° 1, 2013, 44–50.

[4] G. BOCHAROV, B. LUDEWIG, A. BERTOLETTI ET AL, Underwhelming the immune response : effect of slow virus growth on CD 8+ T Lymphocytes responses, *Journal of Virology*, vol. 78, n° 5, 2004, 2247–2254.
[5] K. HATTAF, N. YOUSFI, A generalized virus dynamics model with cell-to-cell transmission and cure rate, *Advances in Difference Equations*, n° 174, 2016.

[6] K. HATTAF, N. YOUSFI, A numerical method for a delayed viral infection model with general incidence rate, *J. King Saud Univ. Sci.*, vol. 28, n° 4, 2016, 368-374.

[7] X.-Y. WANG, K. HATTAF, H.-F. HUO, H. XIANG, Stability analysis of a delayed social epidemics model with general contact rate and its optimal control, *Journal of Industrial and Management Optimization*, vol. 12, n° 4, 2016, 1267-1285.

[8] H. I. FREEDMAN, S. RUAN, M. TANG, Uniform persistence and flows near a closed positively invariant set, *Journal of Dynamics and Differential Equations*, vol. 6, n° 4, 1994, 583-600.

[9] M. Y. LI, J. S. MULDOWNEY, A geometric approach to the global-stability problems, *SIAM Journal on Mathematical Analysis*, vol. 27, n° 4, 1996, 1070-1083.

[10] J. R. MARTIN, H. ROBERT, Logarithmic norms and projections applied to linear differential systems, *Journal of Mathematical Analysis and Applications*, vol. 45, n° 2, 1974, 432-454.

[11] G. BUTLER, P. WALTMAN, Persistence in dynamical systems, *Journal of Differential Equations*, vol. 63, n° 2, 1986, 255-263.