A Nested Multi-Scale Model for COVID-19 Viral Infection

Bishal Chhetri^a,∗, D. K. K. Vamsi^a, Carani B Sanjeevi^b
^aDepartment of Mathematics and Computer Science, Sri Sathya Sai Institute of Higher Learning, Prasanthi Nilayam,
Puttaparthi, Anantapur District - 515134, Andhra Pradesh, India
^b Vice-Chancellor, Sri Sathya Sai Institute of Higher Learning - SSSIHL, India
bishalchhetri@sssihl.edu.in, dkkvamsi@sssihl.edu.in^∗,
sanjeevi.carani@sssihl.edu.in, sanjeevi.carani@ki.se
^∗ Corresponding Author

Abstract

In this study, we develop and analyze a nested multi-scale model for COVID -19 disease that integrates within-host scale and between-host scale sub-models. First, the well-posedness of the multi-scale model is discussed, followed by the stability analysis of the equilibrium points. The disease-free equilibrium point is shown to be globally asymptotically stable for $R_{0}<1$ . When $R_{0}$ exceeds unity, a unique infected equilibrium exists, and the system is found to undergo a forward (trans-critical) bifurcation at $R_{0}=1$ . Two parameter heat plots are also done to find the parameter combinations for which the equilibrium points are stable. The parameters $\beta,\pi$ and $\Lambda$ are found to be most sensitive to $R_{0}$ . The influence of within-host sub-model parameter on the between-host sub-model variables is numerically illustrated. The spread of infection in a community is shown to be influenced by within-host level sub-model parameters, such as the production of viral particles by infected cells $(\alpha)$ , the clearance rate of infected cells by the immune system $(x)$ , and the clearance rate of viral particles by the immune system $(y)$ . The comparative effectiveness of the three health interventions (antiviral drugs, immunomodulators, and generalized social distancing) for COVID-19 infection was examined using the effective reproductive number $R_{E}$ as an indicator of the effectiveness of the interventions. The results suggest that a combined strategy of antiviral drugs, immunomodulators and generalized social distancing would be the best strategy to implement to contain the spread of infection in the community. We believe that the results presented in this study will help physicians, medical professionals, and researchers to make informed decisions about COVID -19 disease prevention and treatment interventions.

keywords

Multi-Scale Models, Well Posedness, Elasticity, Comparative Effectiveness, Heat Plot

1 Introduction

COVID -19 is a contagious respiratory and vascular disease that has resulted in more than 209 million cases and 4.4 million deaths worldwide. It is caused by severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2). It was declared an international public health emergency on January 30.

Mathematical modeling of infectious diseases is one of the highly researched area in applied mathematics. Mathematical epidemiology has contributed to a better understanding of the dynamic behavior of infectious diseases, their impact, and possible future predictions of their spread. In general, the transmission of an infectious disease system is a consequence of processes acting on at least two different scales, ie, within-host and between-host scales [18, 15]. At within-host scale there are three ways in which the transmission of an infectious disease system is regulated that is, when the immune system impacts the dynamics of an infectious disease system in a population: first, by regulating pathogen dynamics at within-host scale; second, by affecting herd immunity at between-host scale; and finally, by altered immuno-demography through change in the immune status due to age, disease, or therapy. There are number of studies on infectious disease systems which has established that transmission of infectious diseases depends critically on within-host scale processes. In particular, these studies established that the transmission potential of an infectious host increases with increasing pathogen load in the infected host [18]. Studies on HIV and dengue virus transmission established a sigmoid functional relationship between host infectiousness and within-host scale pathogen load [26, 30]. Other studies that have also confirmed a similar functional relationship between host infectiousness and within-host scale pathogen load include[20] for malaria, and [21] for human T lymphotropic virus. Studies such as [12, 2] considered transmission rate at between-host scale as a linear function of viral load. Several within-host and between-host compartment models are developed to study the spread of COVID-19 infection at individual level as well as at community level. Some of the important mathematical modelling studies that deal with transmission and spread of COVID-19 at between-host level can be found in [7, 23, 25, 39, 32, 28, 40, 24, 38, 10, 41, 4, 33]. Within-host mathematical modelling studies that deal with the interplay between the viral dynamics and immune response of COVID-19 disease can be found in [8, 36].

A thorough understanding of infectious disease transmission requires knowledge of the processes at different levels of infectious disease and of the interplay between these levels. To get a clear idea of the dynamics of the disease, the different time scale models need to be integrated. Multi-scale models of infectious disease systems integrate the within-host scale and the between-host scale. There are five main different categories of multi-scale models that can be developed at the different levels of organization of an infectious disease system, which are: individual-based multiscale models (IMSMs), nested multiscale models (NMSMs), embedded multiscale models (EMSMs), hybrid multiscale models (HMSMs), and coupled multiscale models (CMSMs) [13]. The multi-scale models for influenza infections can be found in [19, 16].

Multi-scale modeling of COVID -19 disease is still in its infancy, and to our knowledge there is only one literature in which a multi-scale model has been developed for COVID -19 disease [29]. A multi-scale model would be extremely helpful in understanding the spread of COVID -19 infection and evaluate the efficacy of the interventions not only at individual level but also the at the population level. Therefore, in this article, we develop a nested multi-scale model that integrates within-host and between-host sub-models. As the importance of multi-scale modeling in disease dynamics is increasingly recognized, we believe that our study contributes to the growing knowledge on multi-scale modeling of COVID -19 disease. The work in this study is aimed at two types of audiences - disease modelers and public health planners. For disease modelers, this study provides an alternative approach to modeling acute viral infections. For public health planners, this study provides a model-based approach to evaluating the effectiveness of public health interventions.

The paper is organised as follows: In section 2, we develop a nested multi-scale model, study the stability and bifurcation analysis and numerically illustrate the theoretical results obtained. We also do the sensitivity of $R_{0}$ with the model parameters and heat plots. In section 3 we evaluate the comparative effectiveness of the three health interventions. The study concludes with section 4 where we summarize the findings of the work.

2 Multi-Scale Model

The Multi-Scale model that we develop describes COVID-19 disease dynamics across two different time scales, that is, within-host scale and between-host scale. The multi-scale model is based on the monitoring of seven variables, namely, susceptible epithelial cells $U$ , infected epithelial cells $U^{*}$ , viral load $V$ , and four between-host variables namely, susceptible human $S$ , infection human population $I$ , exposed human population $E$ , and recovered human population $R$ . We make the following assumptions for the multi-scale model.

a) The dynamics of the within-host scale variables are assumed to occur at fast time scale $s$ so that $U=U(s)$ , $U^{*}=U^{*}(s)$ and $V=V(s)$ while the dynamics of the between-host scale variables is assumed to occur at slow time scale $t$ so that $S=S(t)$ , $I=I(t)$ , $E=E(t)$ , and $R=R(t)$ .

b) The transmission rate $\beta$ and disease induced death rate $d$ at between-host scale are assumed to be a functions of viral load.

c) Immune response is captured in the model through the parameters $b_{1},b_{2},b_{3},b_{4},b_{5},b_{6}$ and $d_{1},d_{2},d_{3},d_{4},d_{5},d_{6}$ where these parameters denote the rate at which infected cell and virus are cleared by the release of cytokines and chemokines such as IL-6 TNF- $\alpha$ , INF- $\alpha$ , CCL5, CXCL8 , and CXCL10.

Based on the assumptions above the multi-scale model is described by the following system of differential equations.

$\displaystyle\frac{dU}{ds}$	$\displaystyle=$	$\displaystyle\omega\ -kU(s)V(s)-\mu_{c}U(s)$	(1)
$\displaystyle\frac{dU^{*}}{ds}$	$\displaystyle=$	$\displaystyle kU(s)V(s)\ -{\bigg{(}d_{1}+d_{2}+d_{3}+d_{4}+d_{5}+d_{6}\bigg{)}U^{}(s)}\ -\mu_{c}U^{}(s)$	(2)
$\displaystyle\frac{dV}{ds}$	$\displaystyle=$	$\displaystyle\alpha U^{*}(s)\ -\bigg{(}b_{1}+b_{2}+b_{3}+b_{4}+b_{5}+b_{6}\bigg{)}V(s)\ -\mu_{v}V(s)$	(3)
$\displaystyle\frac{dS}{dt}$	$\displaystyle=$	$\displaystyle\Lambda\ -\beta(V(s))S(t)I(t)-\mu S(t)$	(4)
$\displaystyle\frac{dE}{dt}$	$\displaystyle=$	$\displaystyle\beta(V(s))S(t)I(t)\ -(\mu+\pi+\gamma_{1})E(t)$	(5)
$\displaystyle\frac{dI}{dt}$	$\displaystyle=$	$\displaystyle\pi E(t)-(\mu+\gamma_{2})I(t)-d(V(s))I(t)$	(6)
$\displaystyle\frac{dR}{dt}$	$\displaystyle=$	$\displaystyle\gamma_{1}E(t)+\gamma_{2}I(t)-\mu R(t)$	(7)

The meaning of each of the variables and parameters of the model is given in table 1.

Table 1: Meanings of the Parameters

Parameters	Biological Meaning
$\omega$	natural birth rate of cells
$\Lambda$	birth rate of human population
$\alpha$	burst rate of the virus
$\mu$	natural death rate of human population
$\mu_{c}$	natural death rate of cells
$\mu_{v}$	natural death rate of virus
$\pi$	infection rate of exposed population
$k$	infection rate of susceptible cell
$\gamma_{1},\gamma_{2}$	recovary rate of the exposed and infected human population

Considering experimental observations of the impact of viral load on disease transmission [22] and disease induced deaths, the linking of within-host and between-host sub-models is done considering $\beta=\beta(V(s))$ and $d=d(V(s))$ . Though the exact functional relationship between viral load, transmission rate and disease induced death rate is not known [2], as in [12, 2] we considered a linear form of the coupling functions $\beta$ and $d$ : $\beta(V(s))=\beta V(s)$ and $d(V(s))=dV(s)$ . We also observe that the first three equation of the between-host SEIR sub-model is independent of recovered population $R(t)$ . Therefore without loss of generality we omit the last equation of the between-host SEIR sub-model. Taking $\beta$ and $d$ as a linear function of viral load and omitting the equation for $R(t)$ in between-host SEIR sub-model the final multi-scale model for COVID-19 disease dynamics is given by the following set of equations:

$\displaystyle\frac{dU}{ds}$	$\displaystyle=$	$\displaystyle\omega\ -kU(s)V(s)-\mu_{c}U(s)$	(8)
$\displaystyle\frac{dU^{*}}{ds}$	$\displaystyle=$	$\displaystyle kU(s)V(s)\ -{\bigg{(}d_{1}+d_{2}+d_{3}+d_{4}+d_{5}+d_{6}\bigg{)}U^{}(s)}\ -\mu_{c}U^{}(s)$	(9)
$\displaystyle\frac{dV}{ds}$	$\displaystyle=$	$\displaystyle\alpha U^{*}(s)\ -\bigg{(}b_{1}+b_{2}+b_{3}+b_{4}+b_{5}+b_{6}\bigg{)}V(s)\ -\mu_{v}V(s)$	(10)
$\displaystyle\frac{dS}{dt}$	$\displaystyle=$	$\displaystyle\Lambda\ -\beta V(s)S(t)I(t)-\mu S(t)$	(11)
$\displaystyle\frac{dE}{dt}$	$\displaystyle=$	$\displaystyle\beta V(s)S(t)I(t)\ -(\mu+\pi+\gamma_{1})E(t)$	(12)
$\displaystyle\frac{dI}{dt}$	$\displaystyle=$	$\displaystyle\pi E(t)-(\mu+\gamma_{2})I(t)-dV(s)I(t)$	(13)

2.1 Reduced Multi-Scale Model for COVID-19 Dynamics

The multi-scale model given by equations $(2.8)-(2.13)$ is not easy to analyze. The difficulty arises from the mismatch of time scales, since the within-host sub-model is in the form of a fast time scale $s$ , while the between-host sub-model is in the form of a slow time scale $t$ . To overcome these problems, we simplify the multi-scale model by using the within-host scale sub-model to define another quantity, which we use as a proxy for the infectivity of the individual host and which is called the area under the viral load curve [17, 14]. Consider the within-host scale sub-model from the multi-scale model $(2.8)-(2.13)$ .

$\displaystyle\frac{dU}{ds}$	$\displaystyle=$	$\displaystyle\omega\ -kU(s)V(s)-\mu_{c}U(s)$	(14)
$\displaystyle\frac{dU^{*}}{ds}$	$\displaystyle=$	$\displaystyle kU(s)V(s)\ -{xU^{}(s)}\ -\mu_{c}U^{}(s)$	(15)
$\displaystyle\frac{dV}{ds}$	$\displaystyle=$	$\displaystyle\alpha U^{*}(s)\ -yV(s)\ -\mu_{v}V(s)$	(16)

where

x=\bigg{(}d_{1}+d_{2}+d_{3}+d_{4}+d_{5}+d_{6}\bigg{)}

y=\bigg{(}b_{1}+b_{2}+b_{3}+b_{4}+b_{5}+b_{6}\bigg{)}

The viral load $V(s)$ at within-host scale provide the link between the within-host scale and the between-host scale. We use the within-host scale sub-model given by the model system $(2.14)-(2.16)$ to estimate the individual host infectiousness. To derive an expression for the area under the viral load curve we use the within-host sub-model $(2.14)-(2.16)$ . We denote the area under the viral load curve by $N_{h}$ . The quantity $N_{h}$ measures the total amount of COVID-19 virus produced by an infected human throughout the period of host infectivity, and thus can be considered a proxy for individual host infectivity [17]. Let $d_{1}$ and $d_{2}$ be the times at which the viral load takes on the value of the detection limit at the beginning and end of infection. Then $d_{2}-d_{1}$ can be considered as the duration of host infectivity. Using the method described in [14] the area under the viral load curve is given by,

N_{h}=\frac{\alpha\int_{d_{1}}^{d_{2}}U^{*}ds}{y+\mu_{v}}

(17)

Where, $U^{*}$ is the infected cells. For the chosen parameter values $N_{h}$ will have a fixed numerical value.

Now the within-host scale sub-model $(2.14)-(2.16)$ is reduced to a composed parameter $N_{h}$ given by Equation $(2.17)$ . This $N_{h}$ replaces $V(s)$ in the between-host scale sub-model as follows:

$\displaystyle\frac{dS}{dt}$	$\displaystyle=$	$\displaystyle\Lambda\ -\beta N_{h}S(t)I(t)-\mu S(t)$	(18)
$\displaystyle\frac{dE}{dt}$	$\displaystyle=$	$\displaystyle\beta N_{h}S(t)I(t)\ -(\mu+\pi+\gamma_{1})E(t)$	(19)
$\displaystyle\frac{dI}{dt}$	$\displaystyle=$	$\displaystyle\pi E(t)-(\mu+\gamma_{2})I(t)-dN_{h}I(t)$	(20)

The reduced simplified model $(2.18)-(2.20)$ is now at a single time scale $t$ and is much easier to analyse now. We will study the dynamics of COVID-19 disease using the reduced simplified model $(2.18)-(2.20)$ and also study the influence of within-host sub-model parameters on the spread of infection.

2.2 Well-Posedness of the Model

The existence, the positivity, and the boundedness of the solutions of the proposed model $(2.18)-(2.20)$ need to be proved to ensure that the model has a mathematical and biological meaning.

2.2.1 Positivity and Boundedness

Positivity:

Lemma 1.

Let $t_{0}>0$ and $S(t_{0})>0$ , $E(t_{0})>0$ , $I(t_{0})>0$ then the solution $S(t),E(t)$ and $I(t)$ of the system $(2.18)-(2.20)$ are positive for all $t\geq 0$ .

proof: From equation $(2.18)$ we have

\begin{split}\frac{dS}{dt}&=\Lambda-\beta N_{h}SI-\mu S\\[4.0pt] \frac{dS}{dt}&\geq-(\beta N_{h}SI-\mu)S\\ \frac{dS}{S}&=(\beta N_{h}I-\mu)dt\\ \end{split}

Integrating both sides from $t_{0}$ to $t$ we get

S(t)\geq S(t_{0})e^{-\int_{t_{0}}^{t}(\beta N_{h}I+\mu)}

Therefore $S(t)>0$ for all $t>0$ .

Also from equation $2.19$ we have,

\begin{split}\frac{dE}{dt}&=\beta N_{h}SI-(\mu+\pi+\gamma_{1})E\\[4.0pt] \frac{dE}{dt}&=\frac{\beta N_{h}SIE}{E}-(\mu+\pi+\gamma_{1})E\\ \frac{dE}{E}&=\bigg{(}\frac{\beta N_{h}SI}{E}-(\mu+\pi+\gamma_{1})\bigg{)}dt\\ \frac{dE}{E}&=f(S,E,I)dt\\ \end{split}

where

f(S,E,I)=\frac{\beta N_{h}SI}{E}-(\mu+\pi+\gamma_{1})

Integrating both sides from $t_{0}$ to $t$ we get

E(t)=E(t_{0})e^{\int_{t_{0}}^{t}f(S,E,I)}

Therefore $E(t)>0$ for all $t>0$ .

From the last equation $2.20$ we have

\begin{split}\frac{dI}{dt}&=\pi E-(\mu+\gamma_{2}+dN_{h})I\\[4.0pt] \frac{dI}{dt}&\geq-(\mu+\gamma_{2}+dN_{h})I\\ \frac{dI}{I}&=(\mu+\gamma_{2}+dN_{h})dt\\ \end{split}

Integrating both sides from $t_{0}$ to $t$ we get

I(t)\geq I(t_{0})e^{-\int_{t_{0}}^{t}(\mu+\gamma_{2}+dN_{h})}

Therefore $I(t)>0$ for all $t>0$ . Hence, we conclude that all the solutions of the the system $(2.18)-(2.20)$ remain positive for any time $t>0$ provided that the initial conditions are positive. This establishes the positivity of the solutions of the system $(2.18)-(2.20)$ .

Boundedness:
Let $N(t)=S(t)+E(t)+I(t)$
Now,

\begin{split}\frac{dN}{dt}&=\frac{dS}{dt}+\frac{dE}{dt}+\frac{dI}{dt}\\[4.0pt] &=\Lambda-\mu(S(t)+E(t)+I(t))-\gamma_{1}E-\gamma_{2}I-dN_{h}\\ &\leq\Lambda-\mu N(t)\end{split}

Therefore,

\frac{dN}{dt}+\mu N(t)\leq\Lambda

The integrating factor is given by $e^{\mu t}.$ Therefore, after integration we get,

N(t)\leq\frac{\Lambda}{\mu}

Thus we have shown that the solutions of the system $(2.18)-(2.20)$ is bounded.

Therefore, the biologically feasible region is given by the following set,

\Omega=\bigg{\{}\bigg{(}S(t),E(t),I(t)\bigg{)}\in\mathbb{R}^{3}_{+}:S(t)+E(t)+I(t)\leq\frac{\Lambda}{\mu},\ t\geq 0\bigg{\}}

We summarize the above discussion on boundedness by the following lemma.

Lemma 2.

The set

\Omega=\bigg{\{}\bigg{(}S(t),E(t),I(t)\bigg{)}\in\mathbb{R}^{3}_{+}:S(t)+E(t)+I(t)\leq\frac{\Lambda}{\mu},\ t\geq 0\bigg{\}}

is a positive invariant and an attracting set for system $(2.18)-(2.20)$ .

2.2.2 Existence and Uniqueness of Solution

For the general first order ODE of the form

\dot{x}=f(t,x),\hskip 56.9055ptx(t_{0})=x_{0}

(21)

One would have interest in knowing the answers to the following questions:
(i) Under what conditions solution exists for the $(2.21)$ ?
(ii) Under what conditions unique solution exists for the system $(2.21)$ ?

We use the following theorem to established the existence and uniqueness of solution for our SEI model $(2.18)-(2.20)$ .

Theorem 1.

Let D denote the domain:

|t-t_{0}|\leq a,||x-x_{0}||\leq b,x=(x_{1},x_{2},...,x_{n}),x_{0}=(x_{10},..,x_{n0})

and suppose that $f(t,x)$ satisfies the Lipschitz condition:

||f(t,x_{2})-f(t,x_{1})||\leq k||x_{2}-x_{1}||

(22)

and whenever the pairs $(t,x_{1})$ and $(t,x_{2})$ belong to the domain $D$ , where $k$ is used to represent a positive constant. Then, there exist a constant $\delta>0$ such that there exists a unique (exactly one) continuous vector solution $x(t)$ of the system $(2.21)$ in the interval $|t-t_{0}|\leq\delta$ . It is important to note that condition $(2.22)$ is satisfied by requirement that:

\frac{\partial f_{i}}{\partial x_{j}},i,j=1,2,..,n

be continuous and bounded in the domain D.

Theorem 2.

Existence of Solution
Let $D$ be the domain defined in above such that $(2.22)$ hold. Then there exist a solution of model system of equations $(2.18)-(2.20)$ which is bounded in the domain $D$ .
proof
Let:

$\displaystyle f_{1}$	$\displaystyle=$	$\displaystyle\Lambda\ -\beta N_{h}S(t)I(t)-\mu S(t)$	(23)
$\displaystyle f_{2}$	$\displaystyle=$	$\displaystyle\beta N_{h}S(t)I(t)\ -(\mu+\pi+\gamma_{1})E(t)$	(24)
$\displaystyle f_{3}$	$\displaystyle=$	$\displaystyle\pi E(t)-(\mu+\gamma_{2})I(t)-dN_{h}I(t)$	(25)

We will show that

\frac{\partial f_{i}}{\partial x_{j}},i,j=1,2,..,n

is continuous and bounded in the domain D.

From equation $(2.23)$ we have

$\displaystyle\frac{\partial f_{1}}{\partial S}$	$\displaystyle=$	$\displaystyle-\beta N_{h}I-\mu,\|\frac{\partial f_{1}}{\partial S}\|=\|-\beta N_{h}I-\mu\|<\infty$
$\displaystyle\frac{\partial f_{1}}{\partial E}$	$\displaystyle=0$	$\displaystyle,\|\frac{\partial f_{1}}{\partial E}\|<\infty$
$\displaystyle\frac{\partial f_{1}}{\partial I}$	$\displaystyle=$	$\displaystyle-\beta N_{h}S,\|\frac{\partial f_{1}}{\partial I}\|=\|-\beta N_{h}S\|<\infty$

Similarly from equation $(2.24)$ we have

$\displaystyle\frac{\partial f_{2}}{\partial S}$	$\displaystyle=$	$\displaystyle\beta N_{h}I,\|\frac{\partial f_{2}}{\partial S}\|=\|\beta N_{h}I\|<\infty$
$\displaystyle\frac{\partial f_{2}}{\partial E}$	$\displaystyle=$	$\displaystyle-\mu-\pi-\gamma_{1},\|\frac{\partial f_{2}}{\partial E}\|=\|-(\mu+\pi+\gamma_{1})\|<\infty$
$\displaystyle\frac{\partial f_{2}}{\partial I}$	$\displaystyle=$	$\displaystyle\beta N_{h}S,\|\frac{\partial f_{1}}{\partial I}\|=\|\beta N_{h}S\|<\infty$

Finally from $(2.25)$ we have

$\displaystyle\frac{\partial f_{3}}{\partial S}$	$\displaystyle=$	$\displaystyle 0,\|\frac{\partial f_{3}}{\partial S}\|<\infty$
$\displaystyle\frac{\partial f_{3}}{\partial E}$	$\displaystyle=$	$\displaystyle\pi,\|\frac{\partial f_{3}}{\partial E}\|=\|\pi\|<\infty$
$\displaystyle\frac{\partial f_{3}}{\partial I}$	$\displaystyle=$	$\displaystyle-(\mu+\gamma_{2}+dN_{h}),\|\frac{\partial f_{3}}{\partial I}\|=\|-(\mu+\gamma_{2}+dN_{h})\|<\infty$

Hence we have shown that all the partial derivatives are continuous and bounded. Therefore, Lipschitz condition $(2.22)$ is satisfied. Hence, by theorem 2.1 there exists a unique solution of system $(2.18)-(2.20)$ in the region $D$ .

2.3 Stability Analysis

The basic reproduction number denoted by $R_{0}$ that gives the average number of secondary cases per primary case is calculated using the next generation matrix method [11] and the expression for $R_{0}$ for the system $(2.18)-(2.20)$ is given by

\mathbf{R_{0}}=\mathbf{\frac{\beta N_{h}\pi\Lambda}{\mu(\mu+\pi+\gamma_{1})(\mu+\gamma_{2}+dN_{h})}}\\

(26)

System $(2.18)-(2.20)$ admits two equilibria namely, the infection free equilibrium $E_{0}=\bigg{(}\frac{\omega}{\mu},0,0\bigg{)}$ and the infected equilibrium $E_{1}=(S^{*},E^{*},I^{*})$ where,

S^{*}=\frac{\Lambda}{(\beta N_{h}I^{*}+\mu)}

E^{*}=\frac{(\mu+\gamma_{2}+dN_{h})I^{*}}{\pi}

I^{*}=\frac{\mu(R_{0}-1)}{\beta N_{h}}

Since negative population does not make sense, the existence condition for the infected equilibrium point $E_{1}$ is that $R_{0}>1$ .

2.3.1 Stability Analysis of $E_{0}$

We analyse the stability of equilibrium points $E_{0}$ . This is done based on the nature of the eigenvalues of the jacobian matrix evaluated at $E_{0}$ .

The jacobian matrix of the system $(2.18)-(2.20)$ at the infection free equilibrium $E_{0}$ is given by,

J_{E_{0}}=\begin{pmatrix}-\mu&0&\frac{-\beta N_{h}\Lambda}{\mu}\\ 0&-(\mu+\pi+\gamma_{1})&\frac{\beta N_{h}\Lambda}{\mu}\\ 0&\pi&-(\mu+\gamma_{2}+dN_{h})\end{pmatrix}

The characteristic equation of $J_{E_{0}}$ is given by,

(-\mu-\lambda)\bigg{(}\lambda^{2}+(2\mu+\pi+\gamma_{1}+\gamma_{2}+dN_{h})\lambda-(R_{0}-1)(\mu+\pi+\gamma_{1})(\mu+\gamma_{2}+dN_{h})\bigg{)}

(27)

One of the eigenvalue of characteristic equation $(\ref{111b})$ is $-\mu$ and the other two are the roots of the following quadratic equation.

\lambda^{2}+(2\mu+\pi+\gamma_{1}+\gamma_{2}+dN_{h})\lambda-(R_{0}-1)(\mu+\pi+\gamma_{1})(\mu+\gamma_{2}+dN_{h})

(28)

We see that when $R_{0}<1$ both the eigenvalues of equation $(2.28)$ is negative. Therefore, $E_{0}$ remains locally asymptotically stable whenever $R_{0}<1$ . When $R_{0}>1$ the quadratic equation $(2.28)$ has one positive and one negative root. Therefore the characteristic equation $(2.27)$ has two negative and one positive root. Hence $E_{0}$ becomes unstable in this case. We summarize the local asymptotic stability of infection free equilibrium $E_{0}$ in the following theorem.

Theorem 3.

The infection free equilibrium point $E_{0}$ of system $(2.18)-(2.20)$ is locally asymptotically stable provided $R_{0}<1$ . If $R_{0}$ crosses unity $E_{0}$ loses its stability and becomes unstable.

Global Stability of $E_{0}$

To establish the global stability of the infection free equilibrium $E_{0}$ we make use of the method discussed in Castillo-Chavez et al [5].

Theorem 4.

Consider the following general system,

\begin{split}\frac{dX}{dt}&=F(X,Y)\\[6.0pt] \frac{dY}{dt}&=G(X,Y)\qquad\end{split}

(29)

where $X$ denotes the uninfected population compartments and $Y$ denotes the infected population compartments including latent, infectious etc. Here the function $G$ is such that it satisfies $G(X,0)=0$ . Let $U_{0}=(X_{0},\bar{0})$ denote the equilibrium point of the above general system.

If the following two conditions are satisfied then the infection free equilibrium point $U_{0}$ is globally asymptotically stable for the above general system provided $R_{0}<1$

$A_{1}$ : For the subsystem $\frac{dX}{dt}=F(X,0)$ , $X_{0}$ is globally asymptotically stable.

$A_{2}$ : The function $G=G(X,Y)$ can be written as $G(X,Y)=AY-\widehat{G}(X,Y)$ , where $\widehat{G_{j}}(X,Y)\geq 0$ $\forall\;(X,Y)$ in the biologically feasible region $\Omega$ for j=1,2 and $A=D_{Y}G(X,Y)$ at $(X_{0},\bar{0})$ is a M-matrix(matrix with non-negative off diagonal element).

Theorem 5.

The infection free equilibrium point $E_{0}$ of the system $(2.18)-(2.20)$ is globally asymptotically stable whenever $R_{0}<1$
proof:
We will prove global stability of $E_{0}=(\frac{\omega}{\mu},0,0)$ of system $(2.18)-(2.20)$ by showing that system $(2.18)-(2.20)$ can be written as the above general form and both the conditions $A_{1}$ and $A_{2}$ are satisfied .
Comparing the above general system (29) to the system $(2.18)-(2.20)$ the functions $F$ and $G$ are given by

F(X,Y)=\Lambda-\beta N_{h}SI-\mu S

\hskip 36.98866ptG(X,Y)=\bigg{(}\beta N_{h}SI-(mu+\pi+\gamma_{1})E,\;\pi E-(\mu+\gamma_{2}+dN_{h})I\bigg{)}

where $X=S$ and $Y=(E,I)$
The disease free equilibrium point is $U_{0}=(X_{0},\bar{0}),$ where,

X_{0}=\frac{\Lambda}{\mu}\quad\text{and}\quad\bar{0}=(0,0)

From the stability analysis of $E_{0}$ , we know that $U_{0}$ is locally asymptotically stable iff $R_{0}<1$ . Clearly, we see that $G(X,\bar{0})=(0,\bar{0})$ . Now, we show that $X_{0}=(\frac{\Lambda}{\mu})$ is globally asymptotically stable for the subsystem

\frac{dS}{dt}=F(S,\bar{0})=\omega-\mu S

(30)

The integrating factor is $e^{\mu t}$ and therefore after performing integration on the above equation (30) we get,

S(t)e^{\mu t}=\frac{\Lambda e^{\mu t}}{\mu}+c

As $t\rightarrow\infty$ we get,

S(t)=\frac{\omega}{\mu}

which is independent of c. This independency implies that $X_{0}=\frac{\Lambda}{\mu_{1}}$ is globally asymptotically stable for the subsystem $\frac{dS}{dt}=\Lambda-\mu S$ . So, the assumption $A_{1}$ is satisfied.

Now, we will show that assumption $A_{2}$ holds. First, we will find the matrix $A$ . As per the theorem, $A=D_{Y}G(X,Y)$ at $X=X_{0}$ and $Y=\bar{0}$ . Now

D_{Y}G(X,Y)=\begin{bmatrix}-(\mu+\pi+\gamma_{1})&\beta N_{h}S\\[6.0pt] \pi&-(\mu+\gamma_{2}+dN_{h})\end{bmatrix}

At $X=X_{0}$ and $Y=\bar{0}$ , we obtain,

A=\begin{bmatrix}-(\mu+\pi+\gamma_{1})&\frac{\beta N_{h}\Lambda}{\mu}\\[6.0pt] \pi&-(\mu+\gamma_{2}+dN_{h})\end{bmatrix}

Clearly, matrix $A$ has non-negative off-diagonal elements. Hence, $A$ is a M-matrix. Using $\widehat{G}(X,Y)=AY-G(X,Y)$ , we get,

\widehat{G}(X,Y)\;\;=\;\;\begin{bmatrix}\widehat{G_{1}}(X,Y)\\[6.0pt] \widehat{G_{2}}(X,Y)\end{bmatrix}\;\;=\;\;\begin{bmatrix}\beta N_{h}I(\frac{\Lambda}{\mu}-S)\\[6.0pt] 0\end{bmatrix}

Hence $\widehat{G_{1}}(X,Y)=\beta N_{h}I(\frac{\Lambda}{\mu}-S)\geq 0$ because $S(t)\leq\frac{\Lambda}{\mu}$ and $\widehat{G_{2}}(X,Y)=0$
Thus both the assumptions $A_{1}$ and $A_{2}$ are satisfied and therefore infection free equilibrium point $E_{0}$ is globally asymptotically stable provided $R_{0}<1.$

2.3.2 Stability Analysis of $E_{1}$

The jacobian matrix of the system $(2.18)-(2.20)$ at $E_{1}$ is given by,

J=\begin{pmatrix}-(\beta N_{h}I^{*}+\mu)&0&-\beta N_{h}S\\ \beta N_{h}I^{*}&-(\mu+\pi+\gamma_{1})&\beta N_{h}S^{*}\\ 0&\pi&-(\mu+\gamma_{2}+dN_{h})\end{pmatrix}

The characteristic equation of the jacobian $J$ evaluated at $E_{1}$ is given by,

\lambda^{3}+A_{1}\lambda^{2}+B_{1}\lambda+C_{1}=0\\

(31)

where

A_{1}=3\mu+\pi+\beta N_{h}I^{*}+\gamma_{1}+\gamma_{2}+dN_{h}

B_{1}=(\mu+\beta N_{h}I^{*})(2\mu+\pi+\gamma_{1}+\gamma_{2}+dN_{h})+(\mu+\pi+\gamma_{1})(\mu+\gamma_{2}+dN_{h})-\beta N_{h}\pi S^{*}

C_{1}=(\mu+\beta N_{h}I^{*})\bigg{(}(\mu+\pi+\gamma_{1})(\mu+\gamma_{2}+dN_{h})-\beta N_{h}\pi S^{*}\bigg{)}+\beta^{2}N_{h}^{2}S^{*}I^{*}\pi

Clearly, $A_{1}>0$ . By Routh- Hurwitz criterion, all the roots of characteristic equation $(2.31)$ are negative iff $C_{1}>0$ and $A_{1}B_{1}-C_{1}>0$ .
Simplifying the expression for $C_{1}$ we get,

C_{1}=\mu(\mu+\pi+\gamma_{1})(\mu+\gamma_{2}+dN_{h})(R_{0}-1)

Therefore, $C_{1}>0$ iff $R_{0}>1$ . Hence we conclude that the infected equilibrium point $E_{1}$ exists and remains locally asymptotically stable provided $R_{0}>1.$ and $(A_{1}B_{1}-C_{1})>0$ . In the following theorem we summarize the above discussion on the stability of $E_{1}$ .

Theorem 6.

There exists a unique infected equilibrium point $E_{1}$ of the system $(2.18)-(2.20)$ if the following conditions are satisfied:
(i) $R_{0}>1$
(ii) $(A_{1}B_{1}-C_{1})>0$

2.4 Bifurcation Analysis

We now use the method given by Chavez and Song in [6] to do the bifurcation analysis.

Theorem 7.

Consider a system,

\frac{dX}{dt}=f(X,\phi)

where $X\in\mathbb{R}^{n}$ , $\phi\in\mathbb{R}$ is the bifurcation parameter and $f:\mathbb{R}^{n}\times\mathbb{R}\rightarrow\mathbb{R}^{n}$ where $f\in\mathbb{C}^{2}(\mathbb{R}^{n},\mathbb{R})$ . Let $\bar{0}$ be the equilibrium point of the system such that $f(\bar{0},\phi)=\bar{0},\forall\;\phi\in\mathbb{R}$ . Let the following conditions hold :

1.

For the matrix $A=D_{X}f(\bar{0},0)$ , zero is the simple eigenvalue and all other eigenvalues have negative real parts.
2.

Corresponding to zero eigenvalue, matrix A has non-negative right eigenvector, denoted as $u$ and non-negative left eigenvectors, denoted as $v$ .

Let $f_{k}$ be the $k^{th}$ component of $f$ . Let $a$ and $b$ be defined as follows -

a=\sum_{k,i,j=1}^{n}\Bigg{[}v_{k}w_{i}w_{j}\bigg{(}\frac{\partial^{2}f_{k}}{\partial x_{i}\partial x{j}}(\bar{0},0)\bigg{)}\Bigg{]}

b=\sum_{k,i=1}^{n}\Bigg{[}v_{k}w_{i}\bigg{(}\frac{\partial^{2}f_{k}}{\partial x_{i}\partial\phi}(\bar{0},0)\bigg{)}\Bigg{]}

Then local dynamics of the system near the equilibrium point $\bar{0}$ is totally determined by the signs of $a$ and $b$ . Here are the following conclusions :

1.

If $a>0$ and $b>0$ , then whenever $\phi<0$ with $\mid\phi\mid\ll 1$ , the equilibrium $\bar{0}$ is locally asymptotically stable, and moreover there exists a positive unstable equilibrium. However when $0<\phi\ll 1$ , $\bar{0}$ is an unstable equilibrium and there exists a negative and locally asymptotically stable equilibrium.
2.

If $a<0$ , $b<0$ , then whenever $\phi<0$ with $\mid\phi\ll 1$ , $\bar{0}$ is an unstable equilibrium whereas if $0<\phi\ll 1$ , $\bar{0}$ is locally asymptotically stable equilibrium and there exists a positive unstable equilibrium.
3.

If $a>0$ , $b<0$ , then whenever $\phi<0$ with $\mid\phi\mid\ll 1$ , $\bar{0}$ is an unstable equilibrium, and there exists a locally asymptotically stable negative equilibrium. However if $0<\phi\ll 1$ , $\bar{0}$ is stable, and a there appears a positive unstable equilibrium.
4.

If $a<0$ , $b>0$ , then whenever $\phi$ changes its value from negative to positive, the equilibrium $\bar{0}$ changes its stability from stable to unstable. Correspondingly a negative equilibrium, unstable in nature, becomes positive and locally asymptotically stable.

Applying the Theorem 2.7 to our system $(2.18)-(2.20)$ :

In our case, we have $x=(S,E,I)\in\mathbb{R}^{3}$ where $x_{1}=S$ , $x_{2}=E$ and $x_{3}=I$ . Let us consider $\beta$ (transmission rate of the infection) to be the bifurcation parameter.
We know that,

R_{0}=\frac{\beta N_{h}\Lambda\pi}{\mu(\mu+\pi+\gamma_{1})(\mu+\gamma_{2}+dN_{h})}

Therefore we have,

\beta=\frac{R_{0}\mu(\mu+\pi+\gamma_{1})(\mu+\gamma_{2}+dN_{h})}{\beta N_{h}\Lambda\pi}

Let $\beta=\beta^{*}$ at $R_{0}=1$ . So, we have,

\beta^{*}=\frac{\mu(\mu+\pi+\gamma_{1})(\mu+\gamma_{2}+dN_{h})}{\beta N_{h}\Lambda\pi}

With $x=(x_{1},x_{2},x_{3})=(S,I,V)$ system $(2.18)-(2.20)$ can be written as follows :

	$\displaystyle\frac{dx_{1}}{dt}$	$\displaystyle=\Lambda-\beta N_{h}x_{1}x_{3}-\mu x_{1}=f_{1}$
	$\displaystyle\frac{dx_{2}}{dt}$	$\displaystyle=\beta N_{h}x_{1}x_{3}-(\mu+\pi+\gamma_{1})x_{2}=f_{2}$
	$\displaystyle\frac{dx_{3}}{dt}$	$\displaystyle=\pi x_{2}-(\mu+\gamma_{2}+dN_{h})x_{3}\quad=f_{3}$

The disease free equilibrium point $E_{0}$ is given by,

x^{*}=\bigg{(}\frac{\Lambda}{\mu},0,0\bigg{)}=(x_{1}^{*},x_{2}^{*},x_{3}^{*})

Clearly, $f(x^{*},\beta)=0,\;\forall\;\beta\in\mathbb{R}$ , where $f=(f_{1},f_{2},f_{3})$ . Let $D_{x}f(x^{*},\beta^{*})$ denote the Jacobian matrix of the above system at the equilibrium point $x^{*}$ and $R_{0}=1$ . Now we see that,

D_{x}f(x^{*},\beta^{*})=\begin{bmatrix}-\mu&0&\frac{-\beta^{*}N_{h}\Lambda}{\mu}\\[6.0pt] 0&-(\mu+\pi+\gamma_{1})&\frac{\beta^{*}N_{h}\Lambda}{\mu}\\[6.0pt] 0&\pi&-(\mu+\gamma_{2}+dN_{h})\end{bmatrix}

The characteristic polynomial of the above matrix $D_{x}f(x^{*},\beta^{*})$ is given by,

(-\mu-\lambda)\bigg{[}((\mu+\pi+\gamma_{1})+\lambda)((\mu+\gamma_{2}+dN_{h})+\lambda)-\bigg{(}\frac{\alpha\beta^{*}N_{h}\Lambda\pi}{\mu}\bigg{)}\bigg{]}=0

(32)

Hence, we obtain the first eigenvalue of (32) as

\boldsymbol{\lambda_{1}=-\mu<0}

The other eigenvalues $\lambda_{2,3}$ of (32) are the solutions of the following equation,

\lambda^{2}+\bigg{(}2\mu+\pi+\gamma_{1}+\gamma_{2}+dN_{h}\bigg{)}\lambda+(\mu+\pi+\gamma_{1})(\mu+\gamma_{2}+dN_{h})-\frac{\beta^{*}N_{h}\pi\Lambda}{\mu}=0

(33)

substituting the expression for $\beta^{*}$ in (33) we get,

\lambda^{2}+\bigg{(}2\mu+\pi+\gamma_{1}+\gamma_{2}+dN_{h}\bigg{)}\lambda=0

(34)

The eigen values of (34) are $\lambda_{2}=0$ and $\lambda_{3}=-(2\mu+\pi+\gamma_{1}+\gamma_{2}+dN_{h})$

Hence, the matrix $D_{x}f(x^{*},\beta^{*})$ has zero as its simple eigenvalue and all other eigenvalues with negative real parts. Thus, the condition 1 of the theorem 2.7 is satisfied.
Next, for proving condition 2, we need to find the right and left eigenvectors of the zero eigenvalue ( $\lambda_{2}$ ). Let us denote the right and left eigen vectors by $\boldsymbol{w}$ and $\boldsymbol{v}$ respectively. To find $w$ , we use $(D_{x}f(x^{*},\beta^{*})-\lambda_{2}I_{d})w=0$ , which implies that

\ \begin{bmatrix}-\mu&0&\frac{-\beta^{*}N_{h}\Lambda}{\mu}\\[6.0pt] 0&-(\mu+\pi+\gamma_{1})&\frac{\beta^{*}N_{h}\Lambda}{\mu}\\[6.0pt] 0&\pi&-(\mu+\gamma_{2}+dN_{h})\par\end{bmatrix}\begin{bmatrix}w_{1}\\[6.0pt] w_{2}\\[6.0pt] w_{3}\\[6.0pt] \end{bmatrix}\;\;=\;\;\begin{bmatrix}0\\[6.0pt] 0\\[6.0pt] 0\end{bmatrix}

where $w=(w_{1},w_{2},w_{3})^{T}$ . As a result, we obtain the system of simultaneous equations as follows :

-\mu w_{1}-\frac{\beta^{*}N_{h}\Lambda}{\mu}w_{3}=0

(35)

-(\mu+\pi+\gamma_{1})w_{2}+\frac{\beta^{*}N_{h}\Lambda}{\mu}w_{3}=0

(36)

\pi w_{2}-(\mu+\gamma_{2}+dN_{h})w_{3}=0

(37)

By choosing $w_{3}=\mu$ in the above simultaneous equation (35)-(37) we obtain

w_{2}=\frac{\beta^{*}N_{h}\Lambda}{(\mu+\pi+\gamma_{1})}\;\;\;\text{and}\;\;\;{w_{1}=-\frac{\beta^{*}N_{h}\Lambda}{\mu}}

Therefore, the right eigen vector of zero eigenvalue is given by

\boldsymbol{w=\bigg{(}-\frac{\beta^{*}N_{h}\Lambda}{\mu},\;\frac{\beta^{*}N_{h}\Lambda}{(\mu+\pi+\gamma_{1})},\;\mu\bigg{)}}

Similarly, to find the left eigenvector $v$ , we use $v(D_{x}f(x^{*},\beta^{*})-\lambda_{2}I_{d})=0$ , which implies that

\begin{bmatrix}v_{1}&v_{2}&v_{3}\end{bmatrix}\begin{bmatrix}-\mu&0&\frac{-\beta^{*}N_{h}\Lambda}{\mu}\\[6.0pt] 0&-(\mu+\pi+\gamma_{1})&\frac{\beta^{*}N_{h}\Lambda}{\mu}\\[6.0pt] 0&\pi&-(\mu+\gamma_{2}+dN_{h})\par\end{bmatrix}\;\;=\;\;\begin{bmatrix}0&0&0\end{bmatrix}

where $v=(v_{1},v_{2},v_{3})$ . The simultaneous equations obtained thereby are as follows :

-\mu v_{1}=0

(38)

-(\mu+\pi+\gamma_{1})v_{2}+\pi v_{3}=0

(39)

\frac{-\beta^{*}N_{h}\Lambda}{\mu}v_{1}+\frac{\beta^{*}N_{h}\Lambda}{\mu}v_{2}-(\mu+\gamma_{2}+dN_{h})v_{3}=0

(40)

Therefore solving the above simultaneous equation (38 - 40) we obtain $v_{1}=0$ .

By choosing $v_{2}=1$ we get

v_{3}=\frac{\beta^{*}N_{h}\Lambda}{(\mu(\mu+\gamma_{2}+dN_{h}))}

Hence, the left eigen vector is given by

\boldsymbol{v=\bigg{(}0,\;1,\;\frac{\beta^{*}N_{h}\Lambda}{(\mu(\mu+\gamma_{2}+dN_{h}))}\;\bigg{)}}

Now, we need to find $a$ and $b$ . As per the Theorem 2.7, $a$ and $b$ are given by

\boldsymbol{a=\sum_{k,i,j=1}^{3}\Bigg{[}v_{k}u_{i}u_{j}\bigg{(}\frac{\partial^{2}f_{k}}{\partial x_{i}\partial x_{j}}(x^{*},\beta^{*})\bigg{)}\Bigg{]}}

\boldsymbol{b=\sum_{k,i=1}^{3}\Bigg{[}v_{k}u_{i}\bigg{(}\frac{\partial^{2}f_{k}}{\partial x_{i}\partial\beta}(x^{*},\beta^{*})\bigg{)}\Bigg{]}}

Expanding the summation in the expression for $a$ , it reduces to

a=w_{1}w_{3}\frac{\partial^{2}f_{2}}{\partial x_{1}\partial x_{3}}\thinspace+\thinspace\thinspace w_{3}w_{1}\frac{\partial^{2}f_{2}}{\partial x_{3}\partial x_{1}}

where partial derivatives are found at $(x^{*},\beta^{*})$ . Now

\frac{\partial^{2}f_{2}}{\partial x_{1}\partial x_{3}}(x^{*},\beta^{*})=\beta^{*}N_{h}\qquad\frac{\partial^{2}f_{2}}{\partial x_{3}\partial x_{1}}(x^{*},\beta^{*})=\beta^{*}N_{h}

Substituting these partial derivatives along with $w$ in the expression of $a$ , we get,

\boldsymbol{a=-2\beta^{*2}N_{h}^{2}\Lambda<0}

Next, expanding the summation in the expression for $b$ , we get,

b=v_{2}w_{2}\bigg{(}\frac{\partial^{2}f_{2}}{\partial x_{1}\partial\beta}(x^{*},\beta^{*})\bigg{)}+v_{2}w_{3}\bigg{(}\frac{\partial^{2}f_{2}}{\partial x_{3}\partial\beta}(x^{*},\beta^{*})\bigg{)}

Now

\frac{\partial^{2}f_{2}}{\partial x_{3}\partial\beta}(x^{*},\beta^{*})=\frac{N_{h}\lambda}{\mu}

\frac{\partial^{2}f_{2}}{\partial x_{1}\partial\beta}(x^{*},\beta^{*})=0

substituting in the expression of $b$ we get,

\boldsymbol{b=N_{h}\Lambda>0}

Hence $a<0$ and $b>0$ . We notice that condition (iv) of the theorem 2.7 is satisfied. Hence, we conclude that the system undergoes bifurcation at $\beta=\beta^{*}$ implying $R_{0}=1$ .

Thus, we conclude that when $R_{0}<1$ , there exists a unique disease free equilibrium which is globally asymptotically stable and negative infected equilibrium which is unstable . Since negative values of population is not practical, therefore we ignore it in this case. Further, as $R_{0}$ crosses unity from below, the disease free equilibrium point loses its stable nature and become unstable, the bifurcation point being at $\beta=\beta^{*}$ implying $R_{0}=1$ and there appears a positive locally asymptotically stable infected equilibrium point. There is an exchange of stability between disease free equilibrium and infected equilibrium at $R_{0}=1$ . Hence, a forward bifurcation (trans-critical bifurcation) takes place at the break point $\beta=\beta^{*}$ . We summarize the above discussion on bifurcation by the following theorem.

Theorem 8.

As $R_{0}$ crosses unity the disease free equilibrium changes its stability from stable to unstable and there exists a locally asymptotically infected equilibrium when $R_{0}>1$ i.e . direction of bifurcation is forward (transcritical) at $R_{0}>1$ .

2.5 Sensitivity and Elasticity

The Basic Reproduction number denoted by $R_{0}$ is one of the most important quantity in any infectious disease models. The expression for $R_{0}$ for the reduced multi-scale model $(2.18)-(2.20)$ is given by,

R_{0}=\frac{\beta N_{h}\Lambda\pi}{\mu(\mu+\pi+\gamma_{1})(\mu+\gamma_{2}+dN_{h})}

To determine best control measures, knowledge of the relative importance of the different factors responsible for transmission is useful. Initially disease transmission is related to $R_{0}$ and sensitivity predicts which parameters have a high impact on $R_{0}$ . The sensitivity index of $R_{0}$ with respect to a parameter $\mu$ is $\frac{\partial R_{0}}{\partial\mu}$ : Another measure is the elasticity index (normalized sensitivity index) that measures the relative change of $R_{0}$ with respect to $\mu$ , denoted by $\phi_{\mu}^{R_{0}}$ , and defined as

\phi_{\mu}^{R_{0}}=\frac{\partial R_{0}}{\partial\mu}\frac{\mu}{R_{0}}

The sign of the elasticity index tells whether $R_{0}$ increases (positive sign) or decreases (negative sign) with the parameter; whereas the magnitude determines the relative importance of the parameter [3, 35]. These indices can guide control by indicating the most important parameters to target, although feasibility and cost play a role in practical control strategy. If $R_{0}$ is known explicitly, then the elasticity index for each parameter can be computed explicitly, and evaluated for a given set of parameters.

The elasticity index of $R_{0}$ with the model parameters is given by,

\phi_{\beta}^{R_{0}}=\frac{\partial R_{0}}{\partial\beta}\frac{\beta}{R_{0}}=1

\phi_{\Lambda}^{R_{0}}=\frac{\partial R_{0}}{\partial\Lambda}\frac{\Lambda}{R_{0}}=1

\phi_{\pi}^{R_{0}}=\frac{\partial R_{0}}{\partial\pi}\frac{\pi}{R_{0}}=1

\phi_{\mu}^{R_{0}}=\frac{\partial R_{0}}{\partial\mu}\frac{\mu}{R_{0}}=\frac{-\beta N_{h}\Lambda\pi\bigg{(}3\mu^{2}+2\mu(\pi+\gamma_{1}+\gamma_{2}+dN_{h})+(\gamma_{2}+d)(\pi+\gamma_{1})\bigg{)}}{\bigg{(}\mu^{3}+\mu^{2}(\pi+\gamma_{1}+\gamma_{2}+dN_{h})+\mu(\gamma_{2}+d)(\pi+\gamma_{1})\bigg{)}}\frac{\mu}{R_{0}}

\phi_{\gamma_{1}}^{R_{0}}=\frac{\partial R_{0}}{\partial\gamma_{1}}\frac{\gamma_{1}}{R_{0}}=-\frac{\gamma_{1}}{(\mu+\pi+\gamma_{1})}

\phi_{\gamma_{2}}^{R_{0}}=\frac{\partial R_{0}}{\partial\gamma_{2}}\frac{\gamma_{2}}{R_{0}}=-\frac{\gamma_{2}}{(\mu+\gamma_{1}+dN_{h})}

\phi_{N_{h}}^{R_{0}}=\frac{\partial R_{0}}{\partial N_{h}}\frac{N_{h}}{R_{0}}=-\frac{\mu+\gamma_{2}}{(\mu+\gamma_{2}+dN_{h})}

\phi_{d}^{R_{0}}=\frac{\partial R_{0}}{\partial d}\frac{d}{R_{0}}=-\frac{dN_{h}}{(\mu+\gamma_{2}+dN_{h})}

The elasticity indices of the parameters of $R_{0}$ are given in Table 2. The elastic index of parameters $\beta,\pi,\Lambda$ and $N_{h}$ are positive and the remaining are negative. This implies that the increase in the values of these parameters increases $R_{0}$ , whereas increase in the values of parameters $\mu,\gamma_{1}$ and $\gamma_{2}$ decreases $R_{0}$ . For parameter $\beta$ , $\phi_{\beta}^{R_{0}}=1$ implies an increase (decrease) of $\beta$ by $y\%$ increases (decreases) $R_{0}$ by the same percentage. From table 2 we see that the basic reproduction number is most sensitive to the the parameters $\beta,\pi$ and $\Lambda$ . The implication of this is that an increase in the transmission rate increases the spread of the disease in the community.

Table 2: Elasticity Indices of

R_{0}

Parameters	Elastic Index	Value
$\beta$	$\phi_{\beta}^{R_{0}}$	1
$\pi$	$\phi_{\pi}^{R_{0}}$	1
$\Lambda$	$\phi_{\Lambda}^{R_{0}}$	1
$\mu$	$\phi_{\mu}^{R_{0}}$	-0.2785
$\gamma_{1}$	$\phi_{\gamma_{1}}^{R_{0}}$	-0.3196
$\gamma_{2}$	$\phi_{\gamma_{2}}^{R_{0}}$	-0.00082
$N_{h}$	$\phi_{N_{h}}^{R_{0}}$	0.0018
$d$	$\phi_{d}^{R_{0}}$	-0.99

2.6 Numerical Illustrations

2.6.1 Numerical Illustrations of the stability of equilibrium points

Now we numerically illustrate the stability of the equilibrium points admitted by the system $(2.18)-(2.20)$ . The simulation is done using matlab software and ode solver ode45 is used to solve the system of equations. The parameter values of the within-host sub-model used in simulation is taken from [9]. These within-host parameter values are used in calculating the area under the viral load curve $N_{h}$ . The parameter values of the reduced multi-scale mode $(2.18)-(2.20)$ are taken from [27, 32]. All the parameter values are listed in table 4. These parameter values are used in illustrating the stability of the equilibrium points admitted by the system $(2.18)-(2.20)$ . We also illustrate the influence of key within-host scale sub-model parameters on the between-host scale sub-model variables.

Table 3: parameter Values

Symbols	Values	Source
$\omega$	$2$	[9]
$k$	$0.05$	[9]
$\mu_{c}$	$0.1$	[9]
$\mu_{v}$	$0.1$	[9]
$\alpha$	$0.24$	[27]
$d_{1},d_{2},d_{3},d_{4},d_{5},d_{6}$	$0.027,0.22,0.1,0.428,0.01,0.01$	[9]
$b_{1},b_{2},b_{3},b_{4},b_{5},b_{6}$	$0.1,0.1,0.08,0.11,0.01,0.07$	[9]
$\Lambda$	$\mu N(0)$	[27]
$\beta$	0.0115	[32]
$\mu$	0.062	[27]
$\pi$	0.09	[32]
$d$	0.0018	[32]
$\gamma_{1},\gamma_{2}$	0.05, 0.0714	[32]

The initial values used in the simulation is given in table 5.

Table 4: Initial Values of the Variables

Variable	Initial Values	Source
$U(s)$	$3.2*10^{5}$	[9]
$U^{*}(s)$	$0$	[9]
$V(s)$	$5.2$	[9]
$S(t)$	$1000$	Assumed
$E(t)$	$100$	Assumed
$I(t)$	$50$	Assumed

With the parameter values from table 3, the area under the viral load curve $N_{h}$ is found to be $3.3759\times 10^{4}$ . The value of basic reproduction number $\boldsymbol{R_{0}}$ with $\beta=0.00115$ , $\mu=0.72$ and other parameters from table 3 is calculated to be $\boldsymbol{R_{0}}=0.84$ . Since $R_{0}<1$ , by theorem 2.5 the disease free equilibrium point, $\boldsymbol{E_{0}=(98.99,0,0)}$ is globally asymptotically stable for the system $(2.18)-(2.20)$ . The plots in figure 1 for different initial conditions depict the global stability of the infection free equilibrium $E_{0}$ of the system $(2.18)-(2.20)$ .

Refer to caption — Figure 1: Figure depicting the stability of $E_{0}$ of the system $(2.18)-(2.20)$ . The first frame in the figure depicts the local asymptotic stability and the second frame depicts the global stability.

We know that the infected equilibrium $E_{1}$ exists only if $R_{0}>1$ and it is also locally asymptotically stable whenever $R_{0}>1$ . For the parameter values in the following table 3 the value of $R_{0}$ was calculated to be 135.7936 and $E_{1}$ to be $(8.4687,316.8,165.03)$ . Figure 2 demonstrates that $E_{1}$ is locally asymptotically stable whenever $R_{0}>1$ . In figure 3 the occurrence of the trans-critical bifurcation at $R_{0}=1$ is shown.

2.6.2 Heat Plots

Here we vary two parameters of the model $(2.18)-(2.20)$ at a time in a certain interval and plot the the value of $R_{0}$ as heat plots. Heat plot has two different colours: one corresponding to the region with $R_{0}<1$ and the other corresponding to $R_{0}>1$ . This plot helps to find the combination of parameter values for which $R_{0}<1$ and $R_{0}>1$ . The blue colour region in these plots corresponds to the region where $R_{0}<1$ and therefore, from theorem 2.5, the disease free equilibrium is globally stable in this region. The other region with green colour corresponds to the region where $R_{0}>1$ and In these region, the disease free equilibrium is unstable and there exists a unique infected equilibrium point whose stability is determined using theorem 2.6. The stability of infected equilibria is determined using theorem 2.6.

Parameter $\mu$ and $d$
Heat plot varying $\mu$ in the interval $(0.5,0.94)$ and $d$ in the interval $(0.001,0.44)$ is given in figure 4

Parameter $\pi$ and $d$
Heat plot varying $\pi$ in the interval $(0.0005,0.005)$ and $d$ in the interval $(0.05,0.5)$ is given in figure 5.

Parameter $\mu$ and $\pi$
Heat plot varying $\pi$ in the interval $(0.004,0.049)$ and $\mu$ in the interval $(1.05,1.5)$ is given in figure 6.

Parameter $\beta$ and $d$
Heat plot varying $\beta$ in the interval $(0.0005,0.005)$ and $d$ in the interval $(0.05,0.5)$ is given in figure 7.

2.6.3 Influence of within-host sub-model parameters on the between-host sub-model

In this section we study the influence of the within-host sub-model parameters on the between-host sub-model variables. The reduced multi-scale model $(2.18)-(2.20)$ is categorized as a nested multi-scale model according to the categorization of multi-scale models infectious disease systems [13]. Therefore, the multi-scale model $(2.18)-(2.20)$ is unidirectionally coupled. In this only the within-host scale sub-model influences the between-host scale sub-model without any reciprocal feedback. Here we illustrate the influence of the key within-host scale sub-model parameters such as $\alpha,y$ and $x$ on the between-host scale sub-model variable $I$ . At the within-host scale sub-model, $\alpha,y$ and $x$ describes the production of virus by infected cells, the clearance of the infected cells by the immune system, and the clearance rate of free virus particles by the immune system. The parameter values for the simulation are taken from table 3.

In figure 8, the effect of variation of burst rate of virus on infected population is plotted. The infected population $I(t)$ is plotted for three different values of burst rate $\alpha=0.24$ , $\alpha=0.5$ , and $\alpha=0.7$ . These numerical illustration show that the between-host scale variable $I(t)$ is influenced by the within-host scale parameter $\alpha$ . We see from figure 8 that as the infected cell burst size increases, transmission of the disease in the community also increases. Therefore, antiviral drugs such as remdesivir, arbidol, and HCQ that reduce the average infected cell viral production rate at within-host scale will possibly reduce the transmission of COVID-19 disease at between-host scale.

In figure 9, the effect of variation of rate of clearance of infected cells by the immune system on infected population is plotted. The infected population $I(t)$ is plotted for three different values $x=0.5$ , $x=0.795$ , and $x=0.85$ . These numerical illustration show that the between-host scale variable $I(t)$ is influenced by the within-host scale parameter $x$ . We see from figure 9 that as rate of clearance of infected cells by the immune system increases, the infection at population level decreases. Therefore, drugs that kill infected cells could have community level benefits of reducing COVID-19 transmission.

In figure 10, the effect of variation of rate of clearance of virus particles by the immune system on infected population is plotted. The infected population $I(t)$ is plotted for three different values $y=0.56$ , $y=0.7$ , and $y=0.8$ . These numerical illustration show that the between-host scale variable $I(t)$ is influenced by the within-host scale parameter $y$ . We see from figure 10 that as rate of clearance of virus by the immune system increases, the infection in the community decreases. Therefore, in addition to the benefits at individual level, treatments that increase the clearance rate of free virus particles in an infected individual could also have potential benefits in reducing the transmission at between-host level.

3 Comparative Effectiveness of Health Interventions

In this section, in similar lines to the study done in [14, 29] we extend the work on the multi-scale model discussed in Section 2 by including health care interventions. In the context of infectious diseases, disease dynamics can cause a large difference between the performance of a health intervention at the individual level (within-host scale), which is easy to determine, and its performance at the population level (between-host scale), which is difficult to determine [14]. As a result, in situations where the effectiveness of a health intervention cannot be determined, health interventions with proven effectiveness may be recommended over those with potentially higher comparable effectiveness [14]. We use three health care interventions that are as follows:

(1) Antiviral Drugs:

(a) Drugs such as Remdesivir inhibit RNA-dependent RNA polymerase and drugs Lopinavir and Ritonavir inhibit the viral protease there by reducing viral replication [34]. Considering these antiviral drugs are administered the burst rate of the within-host scale sub-model $\alpha$ now gets modified to $\alpha(1-\epsilon)$ where $\epsilon$ is the efficacy of antiviral drugs and $0<\epsilon<1.$

(b) HCQ acts by preventing the SARS-CoV-2 virus from binding to cell membranes and also CQ/HCQ could inhibit viral entry by acting as inhibitors of the biosynthesis of sialic acids, critical actors of virus-cell ligand recognition [31]. The administration of this drug decreases the infection rate $k$ of susceptible cells.This parameter gets modified to become $k(1-\gamma)$ where $\gamma$ is the efficacy of these drug and $0<\gamma<1$ .

(2) Immunomodulators: Interferons are broad spectrum antivirals, exhibiting both direct inhibitory effect on viral replication and supporting an immune response to clear virus infection [37]. Due to the administration of immunomodulators such as INF immunity power of an individual increases as a result of which the clearance rate of infected cells and virus particles by immune response increases. Therefore the parameter $x$ and $y$ now gets modified to $x(1+\delta)$ and $y(1+\delta)$ $\delta$ is the efficacy of immunomodulators and $0<\gamma<1$ .

(3) Generalized social distancing: This intervention involves introducing measures such as restriction of mass gatherings, wearing of mask, temporarily closing schools etc. Assuming that generalized social distancing is introduced to control COVID-19 epidemic, then the rate of contact with community $\beta$ gets modified to $\beta(1-\rho)$ where $\rho$ is the efficacy of generalized social distancing and $0<\rho<1$ .

Considering all the modifications in the parameters, the multi-scale model incorporating the effects of all the above health interventions becomes

$\displaystyle\frac{dS}{dt}$	$\displaystyle=$	$\displaystyle\Lambda\ -\beta(1-\rho)N_{m}S(t)I(t)-\mu S(t)$	(41)
$\displaystyle\frac{dE}{dt}$	$\displaystyle=$	$\displaystyle\beta N_{m}S(t)I(t)\ -(\mu+\pi+\gamma_{1})E(t)$	(42)
$\displaystyle\frac{dI}{dt}$	$\displaystyle=$	$\displaystyle\pi E(t)-(\mu+\gamma_{2})I(t)-dN_{m}I(t)$	(43)

where $N_{m}$ is modified $N_{h}$ given by

N_{m}=\frac{\alpha(1-\epsilon)\int_{d_{1}}^{d_{2}}U^{*}ds}{y(1+\delta)+\mu_{v}}

(44)

The effective reproduction number after incorporating the health interventions is given by,

R_{E}=\frac{\beta(1-\rho)N_{m}\Lambda\pi}{\mu(\mu+\pi+\gamma_{1})(\mu+\gamma_{2}+dN_{m})}

Here the comparative effectiveness of the above three health interventions are evaluated using $R_{E}$ as indicators of intervention effectiveness. In order to evaluate these we calculate the percentage reduction of $\mathcal{R}_{0}$ for single and multiple combination of the interventions at different efficacy levels such as (a) low efficacy of $0.3$ , (b) medium efficacy of $0.6$ , and (c)high efficacy of $0.9$ . The different efficacy levels of the health interventions are chosen from [29].

Percentage reduction of $\mathcal{R}_{0}$ is given by

\bigg{[}\frac{\mathcal{R}_{0}-\mathcal{R}_{E_{j}}}{\mathcal{R}_{0}}\bigg{]}\times 100

where $R_{Ej}$ means effective reproductive number when intervention/combination of interventions with efficacy/efficacies $j$ .

We now consider $8$ different combinations of these three health interventions corresponding to efficacy values 0.3, 0.6, and 0.9 obtained using the effective reproductive number $R_{E}$ as the indicator of intervention effectiveness. Then for each efficacy level, we rank the percentage reductions on $\mathcal{R}_{0}$ in ascending order from 1 to 8 corresponding to the different combinations of three health interventions considered in this study. The comparative effectiveness is calculated and measured on a scale from $1$ to $8$ with $1$ denoting the lowest comparative effectiveness and $8$ denoting the highest comparative effectiveness. In Table 5 the abbreviations a) CEL stands for ”Comparative Effectiveness at Low efficacy,” which is 0.3, b) CEM stands for ”Comparative Effectiveness at Medium efficacy,” which is 0.6, c) CEH stands for ”Comparative Effectiveness at High efficacy” which is 0.9.

Table 5: Comparative effectiveness for

\mathcal{R}_{0}

No.	Indicator	%age	CEL	%age	CEM	%age	CEH
1	$\mathcal{R}_{0}$	0	1	0	1	0	1
2	$\mathcal{R_{E_{\rho}}}$	30	5	60	5	90	5
3	$\mathcal{R_{E_{\delta}}}$	0.106	3	0.24	2	0.45	4
4	$\mathcal{R_{E_{\epsilon}}}$	0.075	2	0.26	3	0.38	3
5	$\mathcal{R_{E_{\rho\delta}}}$	30.07	7	60.12	7	90.45	7
6	$\mathcal{R_{E_{\rho\epsilon}}}$	30.05	6	60.1	6	90.23	6
7	$\mathcal{R_{E_{\epsilon\delta}}}$	0.22	4	0.35	4	0.67	2
8	$\mathcal{R_{E_{\rho\delta\epsilon}}}$	30.16	8	60.34	8	90.5	8

From table 5, we deduce the following results regarding comparative effectiveness of three different interventions considered:

(a) When a single intervention strategy is implemented, intervention involving introducing measures such as restriction of mass gatherings, wearing of mask, temporarily closing schools etc. show significant decrease in $R_{0}$ compared to the implementation of antiviral drugs and immunomodulators at all efficacy levels.

(b) When considering two interventions, we observe that the generalized social distancing along with Immunomodulators that boost the immune response would be highly effective in limiting the spread of infection in the community.

(c) A combined strategy involving treatment with anitviral drugs, immunomodulators, and generalized social distancing seems to perform the best among all the combinations considered at all efficacy levels.

From this section, we conclude that a combined strategy is the best strategy to contain the spread of infection. However, the implementation of a combined strategy may not always be cost-effective and feasible. In a single intervention strategy, generalized social distancing is shown to lower the value of $R_{0}$ better than individual use of antiviral drugs and immunomodulators. Therefore, in a situation where resources are limited and costly, interventions such as restricting mass gatherings, wearing masks, temporarily closing schools, etc. would be highly effective in containing infection and also incur less cost.

4 Discussions and Conclusions

COVID-19 is a contagious respiratory and vascular disease that has resulted in more than 195 million cases and 4.18 million deaths worldwide. It is caused by severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2). On 30 january it was declared as a public health emergency of international concern [1]. Mathematical models have proved to provide useful information about the dynamics of the infectious diseases. To understand the dynamics of the COVID-19 disease several compartment models has been developed [32, 28, 40, 24, 38, 10, 41, 7, 4, 33].

An in-depth understanding of the transmission of infectious disease systems requires knowledge of the processes at the various scales of infectious diseases and how these scales interact. In this study, we develop a nested multi-scale model for COVID-19 disease that integrates the within-host scale and the between-host scale sub-models. The transmission rate and COVID-19 induced death rate at between-host scale are assumed to be a linear function of viral load. Because of the difficulties of working at two different time scales, we approximate individual-level host infectiousness by some surrogate measurable quantity called area under viral load and reduce the multi-scale model at two different times scales to a model with area under the viral load curve acting as a proxy of individual level host infectiousness.

Initially, the well posedness of the reduced multi-scale model is discussed followed by the stability analysis of the equilibrium points admitted by the reduced multi-scale model. The disease free equilibrium point of the model is found to remain globally asymptotically stable whenever the value of basic reproduction number is less than unity. As the value of basic reproduction number crossed unity a unique infected equilibrium point exists and remains stable depending on the sign of $(AB-C)$ . The model is shown to undergo a forward (trans-critical) bifurcation at $R_{0}=1$ . To predict the sensitivity of the model parameters on $R_{0}$ , elastic index that measures the relative change of $R_{0}$ with respect to parameters is calculated for each parameters in the definition of $R_{0}$ . The parameters $\beta$ , $\pi$ and $\Lambda$ were found to be most sensitive towards $R_{0}$ . The theoretical results are supported with numerical illustrations. To separate out the region of stability and instability of the equilibrium points, two parameter heat plot is done by varying the parameter values in certain range. The influence of the key within-host scale sub-model parameters such as $\alpha,y$ and $x$ on between-host scale are also numerically illustrated. It is found that the spread of infection in a community is influenced by the production of virus particles by an infected cells ( $\alpha$ ), clearance rate of infected cell $(x)$ and clearance rate of virus particles by immune system $(y)$ .

We also use the reduced multi-scale model developed to study the comparative effectiveness of the three health interventions (antiviral drugs, immunomodulators and generalized social distancing) for COVID-19 viral infection using $R_{E}$ as the indicator of intervention effectiveness. The result suggested that a combined strategy involving treatment with anitviral drugs, immunomodulators, and generalized social distancing would be the best strategy to limit the spread of infection in the community. However, implementing a combined strategy may not be always cost-effective or feasible. In a single intervention strategy, general social distancing has been shown to lower the value of $R_{0}$ better than individual use of antiviral drugs and immunomodulators. Therefore, in a situation where resources are limited and costly, interventions such as restricting mass gatherings, wearing masks, temporarily closing schools, etc., would be highly effective in containing infection and also more cost-effective. We believe that the results presented in this study will help physicians, doctors and researchers in making informed decisions about COVID-19 disease prevention and treatment interventions.

DEDICATION

The authors from SSSIHL dedicate this paper to the founder chancellor of SSSIHL, Bhagawan Sri Sathya Sai Baba. The first author dedicates this paper to his loving father Purna Chhetri and the second author to his loving elder brother D. A. C. Prakash who still lives in his heart.

References

[1] https://www.who.int/publications/m/item/covid-19-public-health-emergency-of-international-concern-(pheic)-global-research-and-innovation-forum.
[2] Alexis Erich S Almocera, Esteban A Hernandez-Vargas, et al., Multiscale model within-host and between-host for viral infectious diseases, Journal of mathematical biology 77 (2018), no. 4, 1035–1057.
[3] Hailay Weldegiorgis Berhe, Oluwole Daniel Makinde, and David Mwangi Theuri, Parameter estimation and sensitivity analysis of dysentery diarrhea epidemic model, Journal of Applied Mathematics 2019 (2019).
[4] Sudhanshu Kumar Biswas, Jayanta Kumar Ghosh, Susmita Sarkar, and Uttam Ghosh, Covid-19 pandemic in india: a mathematical model study, Nonlinear dynamics 102 (2020), no. 1, 537–553.
[5] Zhilan Feng Carlos Castillo-Chavez and Wenzhang Huang, On the computation of reproduction number and its role in global stability., Institute for Mathematics and Its Applications 125 (2002), no. 2, 229–250.
[6] Carlos Castillo-Chavez and Baojun Song, Dynamical models of tuberculosis and their applications., Mathematical Biosciences and Engineering 1 (2004), no. 2, 361–404.
[7] Tian-Mu Chen, Jia Rui, Qiu-Peng Wang, Ze-Yu Zhao, Jing-An Cui, and Ling Yin, A mathematical model for simulating the phase-based transmissibility of a novel coronavirus, Infectious diseases of poverty 9 (2020), no. 1, 1–8.
[8] Bishal Chhetri, Vijay M Bhagat, DKK Vamsi, VS Ananth, Roshan Mandale, Swapna Muthusamy, Carani B Sanjeevi, et al., Within-host mathematical modeling on crucial inflammatory mediators and drug interventions in covid-19 identifies combination therapy to be most effective and optimal, Alexandria Engineering Journal 60 (2021), no. 2, 2491–2512.
[9] Bishal Chhetri, Vijay M Bhagat, DKK Vamsi, VS Ananth, Bhanu Prakash, Roshan Mandale, Swapna Muthusamy, and Carani B Sanjeevi, Within-host mathematical modeling on crucial inflammatory mediators and drug interventions in covid-19 identifies combination therapy to be most effective and optimal, Alexandria Engineering Journal (2020).
[10] Mohammadali Dashtbali and Mehdi Mirzaie, A compartmental model that predicts the effect of social distancing and vaccination on controlling covid-19, Scientific Reports 11 (2021), no. 1, 1–11.
[11] Odo Diekmann, JAP Heesterbeek, and Michael G Roberts, The construction of next-generation matrices for compartmental epidemic models, Journal of the Royal Society Interface 7 (2010), no. 47, 873–885.
[12] Zhilan Feng, Jorge Velasco-Hernandez, Brenda Tapia-Santos, and Maria Conceição A Leite, A model for coupling within-host and between-host dynamics in an infectious disease, Nonlinear Dynamics 68 (2012), no. 3, 401–411.
[13] Winston Garira, A complete categorization of multiscale models of infectious disease systems, Journal of biological dynamics 11 (2017), no. 1, 378–435.
[14] Winston Garira and Dephney Mathebula, Development and application of multiscale models of acute viral infections in intervention research, Mathematical Methods in the Applied Sciences 43 (2020), no. 6, 3280–3306.
[15] Julia R Gog, Lorenzo Pellis, James LN Wood, Angela R McLean, Nimalan Arinaminpathy, and James O Lloyd-Smith, Seven challenges in modeling pathogen dynamics within-host and across scales, Epidemics 10 (2015), 45–48.
[16] Dongmin Guo, King C Li, Timothy R Peters, Beverly M Snively, Katherine A Poehling, and Xiaobo Zhou, Multi-scale modeling for the transmission of influenza and the evaluation of interventions toward it, Scientific reports 5 (2015), no. 1, 1–9.
[17] Christoforos Hadjichrysanthou, Emilie Cauët, Emma Lawrence, Carolin Vegvari, Frank De Wolf, and Roy M Anderson, Understanding the within-host dynamics of influenza a virus: from theory to clinical implications, Journal of The Royal Society Interface 13 (2016), no. 119, 20160289.
[18] Andreas Handel and Pejman Rohani, Crossing the scale from within-host infection dynamics to between-host transmission fitness: a discussion of current assumptions and knowledge, Philosophical Transactions of the Royal Society B: Biological Sciences 370 (2015), no. 1675, 20140302.
[19] Frank S Heldt, Timo Frensing, Antje Pflugmacher, Robin Gröpler, Britta Peschel, and Udo Reichl, Multiscale modeling of influenza a virus infection supports the development of direct-acting antivirals, PLoS computational biology 9 (2013), no. 11, e1003372.
[20] Geoffrey M Jeffery and Don E Eyles, Infectivity to mosquitoes of plasmodium falciparum as related to gametocyte density and duration of infection1, The American journal of tropical medicine and hygiene 4 (1955), no. 5, 781–789.
[21] Jonathan E Kaplan, Rima F Khabbaz, Edward L Murphy, Sigurd Hermansen, Chester Roberts, Renu Lal, Walid Heneine, David Wright, Lauri Matijas, Ruth Thomson, et al., Male-to-female transmission of human t-cell lymphotropic virus types i and ii: association with viral load, JAIDS Journal of Acquired Immune Deficiency Syndromes 12 (1996), no. 2, 193–201.
[22] Hitoshi Kawasuji, Yusuke Takegoshi, Makito Kaneda, Akitoshi Ueno, Yuki Miyajima, Koyomi Kawago, Yasutaka Fukui, Yoshihiro Yoshida, Miyuki Kimura, Hiroshi Yamada, et al., Transmissibility of covid-19 depends on the viral load around onset in adult and symptomatic patients, PloS one 15 (2020), no. 12, e0243597.
[23] Adam J Kucharski, Timothy W Russell, Charlie Diamond, Yang Liu, John Edmunds, Sebastian Funk, Rosalind M Eggo, Fiona Sun, Mark Jit, James D Munday, et al., Early dynamics of transmission and control of covid-19: a mathematical modelling study, The lancet infectious diseases (2020).
[24] Alexandros Leontitsis, Abiola Senok, Alawi Alsheikh-Ali, Younus Al Nasser, Tom Loney, and Aamena Alshamsi, Seahir: A specialized compartmental model for covid-19, International Journal of Environmental Research and Public Health 18 (2021), no. 5, 2667.
[25] Qianying Lin, Shi Zhao, Daozhou Gao, Yijun Lou, Shu Yang, Salihu S Musa, Maggie H Wang, Yongli Cai, Weiming Wang, Lin Yang, et al., A conceptual model for the coronavirus disease 2019 (covid-19) outbreak in wuhan, china with individual reaction and governmental action, International journal of infectious diseases 93 (2020), 211–216.
[26] John W Mellors, Charles R Rinaldo, Phalguni Gupta, Roseanne M White, John A Todd, and Lawrence A Kingsley, Prognosis in hiv-1 infection predicted by the quantity of virus in plasma, Science 272 (1996), no. 5265, 1167–1170.
[27] Samuel Mwalili, Mark Kimathi, Viona Ojiambo, Duncan Gathungu, and Rachel Mbogo, Seir model for covid-19 dynamics incorporating the environment and social distancing, BMC Research Notes 13 (2020), no. 1, 1–5.
[28] Faïçal Ndaïrou, Iván Area, Juan J Nieto, and Delfim FM Torres, Mathematical modeling of covid-19 transmission dynamics with a case study of wuhan, Chaos, Solitons & Fractals 135 (2020), 109846.
[29] D Bhanu Prakash, DKK Vamsi, D Bangaru Rajesh, and Carani B Sanjeevi, Control intervention strategies for within-host, between-host and their efficacy in the treatment, spread of covid-19: A multi scale modeling approach, Computational and Mathematical Biophysics 8 (2020), no. 1, 198–210.
[30] Thomas C Quinn, Maria J Wawer, Nelson Sewankambo, David Serwadda, Chuanjun Li, Fred Wabwire-Mangen, Mary O Meehan, Thomas Lutalo, and Ronald H Gray, Viral load and heterosexual transmission of human immunodeficiency virus type 1, New England journal of medicine 342 (2000), no. 13, 921–929.
[31] Eugenia Quiros Roldan, Giorgio Biasiotto, Paola Magro, and Isabella Zanella, The possible mechanisms of action of 4-aminoquinolines (chloroquine/hydroxychloroquine) against sars-cov-2 infection (covid-19): A role for iron homeostasis?, Pharmacological research (2020), 104904.
[32] Piu Samui, Jayanta Mondal, and Subhas Khajanchi, A mathematical model for covid-19 transmission dynamics with a case study of india, Chaos, Solitons & Fractals 140 (2020), 110173.
[33] Kankan Sarkar, Subhas Khajanchi, and Juan J Nieto, Modeling and forecasting the covid-19 pandemic in india, Chaos, Solitons & Fractals 139 (2020), 110049.
[34] Yung-Fang Tu, Chian-Shiu Chien, Aliaksandr A Yarmishyn, Yi-Ying Lin, Yung-Hung Luo, Yi-Tsung Lin, Wei-Yi Lai, De-Ming Yang, Shih-Jie Chou, Yi-Ping Yang, et al., A review of sars-cov-2 and the ongoing clinical trials, International Journal of Molecular Sciences 21 (2020), no. 7, 2657.
[35] Pauline Van den Driessche, Reproduction numbers of infectious disease models, Infectious Disease Modelling 2 (2017), no. 3, 288–303.
[36] Esteban Abelardo Hernandez Vargas and Jorge X Velasco-Hernandez, In-host modelling of covid-19 kinetics in humans, medrxiv (2020).
[37] Ben X Wang and Eleanor N Fish, Global virus outbreaks: Interferons as 1st responders, Seminars in immunology, vol. 43, Elsevier, 2019, p. 101300.
[38] Tianbing Wang, Yanqiu Wu, Johnson Yiu-Nam Lau, Yingqi Yu, Liyu Liu, Jing Li, Kang Zhang, Weiwei Tong, and Baoguo Jiang, A four-compartment model for the covid-19 infection—implications on infection kinetics, control measures, and lockdown exit strategies, Precision Clinical Medicine 3 (2020), no. 2, 104–112.
[39] Chayu Yang and Jin Wang, A mathematical model for the novel coronavirus epidemic in wuhan, china, Mathematical Biosciences and Engineering 17 (2020), no. 3, 2708–2724.
[40] Anwar Zeb, Ebraheem Alzahrani, Vedat Suat Erturk, and Gul Zaman, Mathematical model for coronavirus disease 2019 (covid-19) containing isolation class, BioMed research international 2020 (2020).
[41] Ze-Yu Zhao, Yuan-Zhao Zhu, Jing-Wen Xu, Shi-Xiong Hu, Qing-Qing Hu, Zhao Lei, Jia Rui, Xing-Chun Liu, Yao Wang, Meng Yang, et al., A five-compartment model of age-specific transmissibility of sars-cov-2, Infectious diseases of poverty 9 (2020), no. 1, 1–15.

A Nested Multi-Scale Model for COVID-19 Viral Infection

Abstract

keywords

1 Introduction

2 Multi-Scale Model

2.1 Reduced Multi-Scale Model for COVID-19 Dynamics

2.2 Well-Posedness of the Model

2.2.1 Positivity and Boundedness

Lemma 1.

Lemma 2.

2.2.2 Existence and Uniqueness of Solution

Theorem 1.

Theorem 2.

2.3 Stability Analysis

2.3.1 Stability Analysis of E0E_{0}

Theorem 3.

Theorem 4.

Theorem 5.

2.3.2 Stability Analysis of E1E_{1}

Theorem 6.

2.4 Bifurcation Analysis

Theorem 7.

Theorem 8.

2.5 Sensitivity and Elasticity

2.6 Numerical Illustrations

2.6.1 Numerical Illustrations of the stability of equilibrium points

2.6.2 Heat Plots

2.6.3 Influence of within-host sub-model parameters on the between-host sub-model

3 Comparative Effectiveness of Health Interventions

4 Discussions and Conclusions

References

2.3.1 Stability Analysis of $E_{0}$

2.3.2 Stability Analysis of $E_{1}$