Stochastic Model for the Spread of the COVID-19 Virus

Majdi Elhiwi

Probability and Statistics Laboratory UR/11/ES-76, Faculty of Sciences of Sfax, Sfax, Tunisia.

DOI: 10.4236/am.2021.121003                

Abstract

The COVID-19 pandemic has become a great challenge to scientific, biological and medical research as well as to economic and social sciences. Hence, the objective of infectious disease modeling-based data analysis is to recover these dynamics of infectious disease spread and to estimate parameters that govern these dynamics. The random aspect of epidemics leads to the development of stochastic epidemiological models. We establish a stochastic combined model using numerical scheme Euler, Markov chain and Susceptible-Exposed-Infected-Recovery (SEIR) model. The combined SEIR model was used to predict how epidemics will develop and then to act accordingly. These COVID-19 data were analyzed from several countries such as Italy, Russia, USA and Iran.

Keywords

Stochastic Model, COVID-19, EULER Method, SEIR Model, Markov Chain

1. Introduction

The epidemic outbreak of COVID-19 virus is of great interest to researches because it is fast spreading and contagious. Mathematical analysis and modeling as with Kingman, 1982 [1] are fundamental to infectious disease epidemiology. So, in this paper, we will introduce the stochastic process of disease transmission. Using this stochastic process with its mathematical representation, we can predict the future and, most importantly, we can say if the situation is still under control or not. In this framework, the combined SEIR model allows us to extrapolate from current data the state of the epidemic as well as its progress. Also, it enables us to predict the future and to quantify the uncertainty in these predictions. Accordingly, the combined SIR model with stochastic parameters is of primary importance in planning control measures.

2. Mathematical Models of Infectious Disease Spread

Epidemiology models predict how N individuals among a population move among four groups: susceptible S, exposed E, infective I and removed R. In this context “removed” means individuals who are either recovered from the disease and immune to further infection, or dead as in Tang, 2019 [2]. In this framework, we will study the dynamics of the number of individuals in each state compared to the corona virus in several countries. The cumulative number of individuals for each different state at a given time t will be noted by the initial letter of the state in question, for example I will denote the cumulative number of individuals infected at time t.

2.1. Infectious Disease Representation

In epidemiology, the dynamics of infectious disease are often characterized by the nonlinear and the quick variation in sizes of different population groups (states). One of the most widely used models is the Susceptible-Infectious-Removed (SIR) model Bailey, 1975 [3] Anderson and May, 1992 [4]. Unlike traditional infectious disease epidemiology models, we suppose that infected individual is governed by a d-dimensional Wiener process ${\left(I_t\right)}_{t\in \left[\mathrm{0,}T\right]}$. The solution of the stochastic differential equation is given by the following stochastic integral

$I_t=I_0+{\displaystyle {\int }_0^t}b\left(u\mathrm{,}{I}_u\right)\text{d}u+{\displaystyle {\int }_0^t}\sigma \left(u\mathrm{,}{I}_u\right)\text{d}{B}_u$ (1)

where $b\left(u\mathrm{,}{I}_u\right)\in {ℝ}^d\times {ℝ}^d\to {ℝ}^d$ and $\sigma \left(u\mathrm{,}{I}_u\right)\in {ℝ}^d\times {ℝ}^d\to {ℝ}^{d\times d}$ denote respectively the drift of the disease transmission (infection rate) and the Brownian motion diffusion $\left(B_t\right)\mathrm{,}\forall t\in \left[\mathrm{0,}T\right]$.

2.2. Euler Scheme

In order to approximate the number of affected individuals $I^n=\left(I_t^n\mathrm{<mo>\; ;\; </mo>}t\in \left[\mathrm{0,}T\right]\right)$, we apply the Euler scheme. The principle of this numerical scheme is based on the following points:

1) divide the interval $\left[\mathrm{0,}T\right]$ : $0=t_0<t_1<\cdots <t_{n-1}<t_n=T$ and we denote

$\Pi =\underset{i=0,\cdots ,n-1}{\sup }{t}_{i+1}-t_i.$

2) approximate integrals in Equation (1) and we get

${\displaystyle {\int }_{t_i}^{t_{i+1}}}b\left(s\mathrm{,}{I}_s\right)\text{d}s\simeq b\left(t_i\mathrm{,}{I}_{t_i}\right)\left(t_{i+1}-t_i\right),$

${\displaystyle {\int }_{t_i}^{t_{i+1}}}\sigma \left(s\mathrm{,}{I}_s\right)\text{d}{B}_s\simeq \sigma \left(t_i\mathrm{,}{I}_{t_i}\right)\left(B_{t_{i+1}}-B_{t_i}\right)\mathrm{.}$

3) For $i=0,\cdots ,n-1$, calculate $I_{t_{i+1}}^n$ in terms of $I_{t_i}^n$ as following

$I_{t_{i+1}}^n=I_{t_i}^n+b\left(t_i,I_{t_i}^n\right)\left(t_{i+1}-t_i\right)+\sigma \left(t_i,I_{t_i}^n\right)\left(B_{t_{i+1}}-B_{t_i}\right).$

The sequence of random variables $\left(B_{t_{i+1}}-B_{t_i}\right)$ is that of independent random variables normally distributed.

3. Stochastic SEIR Model

3.1. SIR Model

The mostly widely studied class of epidemic models, and the one on which we focus in this paper, is the class of susceptible/infected-infectious/removed or SIR models. Under SIR (Susceptible-Infectious-Removed) model, the population consists of three groups/compartments categories.

$\boxed{S\mathrm{,}I-\mathrm{1,}R+1}\underset{\gamma I}{\underset{︸}{\leftarrow }}\boxed{S\mathrm{,}I\mathrm{,}R}\underset{\beta SI}{\underset{︸}{\to }}\boxed{S-\mathrm{1,}I+\mathrm{1,}R}$

In each group, the population transitions from the current state $\left(S\mathrm{,}I\mathrm{,}R\right)$ to $\left(S-\mathrm{1,}I+\mathrm{1,}R\right)$ with rate $\beta SI$ and to $\left(S\mathrm{,}I-\mathrm{1,}R+1\right)$ with rate $\gamma I$.

3.2. Combined SEIR Model

We use a Combined Susceptible-Exposed-Infected-Recovery (CSEIR) model to illustrate stochastic model representation. CSEIR model is an extension of the SIR model. The transition events between different states of the CSEIR model are depicted in Figure 1.

The current state is $\left(S\mathrm{,}E\mathrm{,}I\mathrm{,}R\right)$. The infection events appear in transitions to state $\left(S-\mathrm{1,}E+\mathrm{1,}I\mathrm{,}R\right)$ and $\left(S-\mathrm{1,}E+\mathrm{1,}I\mathrm{,}R\right)$ with respectively rates $\delta \beta S$ and $\beta SI$. We suppose that Exposed compartment is proportional to infected group by the factor $\delta$. Infected individual becomes infectious is within the transition to state $\left(S\mathrm{,}E-\mathrm{1,}I+\mathrm{1,}R\right)$ with rate $\alpha SE$. Last transition to state $\left(S\mathrm{,}E\mathrm{,}I-\mathrm{1,}R+1\right)$ explains a removal event with rate $\gamma I$. Next section focuses on the transition probabilities based on the following stochastic differential equations

$\{\begin{array}{l}\frac{\text{d}S}{\text{d}t}=-\frac{\beta S}{N}\left(I+\delta E\right)\hfill \\ \frac{\text{d}I}{\text{d}t}=\frac{\beta SI}{N}+\alpha E-\gamma I\hfill \\ \frac{\text{d}E}{\text{d}t}=\frac{\delta \beta SE}{N}-\alpha E\hfill \\ \frac{\text{d}R}{\text{d}t}=\gamma I.\hfill \end{array}$

Figure 1. Combined SEIR model.

The challenge is to calculate transition probability (see Gillespie, 1977 and 2000 [5] [6], Keeling and Ross, 2007 [7] ) with non constant rates by solving the above system of differential equations and using Euler scheme.

3.3. Transition Probabilities

In order to assess the random aspect of epidemics models, we emphasize on the probabilistic version of the CSEIR model. As mentioned our model divides population into four compartments. Each compartment is modeled as a continuous time Markov chain over discrete state space as in Figure 2. We discretize time to $t_n=n\Delta t=\left\{0,\Delta t,2\Delta t,\cdots \right\}$ where $\Delta t$ is small enough so that, in a time interval $\Delta t$, there is almost one change of states. We assume that recovered individuals never interact with individuals in other states. We get a homogeneous discrete Markov chain ${\left(S_n\mathrm{,}{E}_n\mathrm{,}{I}_n\mathrm{,}{R}_n\right)}_{n\in \left\{\mathrm{0,}\cdots \mathrm{,}+\infty \right\}}$ with the following transition probabilities:

$\begin{array}{l}\text{p}\left(s,e,i\right);\left(s^{\prime },e^{\prime },i^{\prime }\right)\\ =\{\begin{array}{l}\mathbb{P}\left[\left(S_{n+1},E_{n+1},I_{n+1}\right)=\left(s^{\prime },e^{\prime },i^{\prime }\right)|\left(S_n,E_n,I_n\right)=\left(s,e,i\right)\right]\\ \mathbb{P}\left[\left(S_{n+1},E_{n+1},I_{n+1}\right)=\left(s^{\prime }=s+k,e^{\prime }=e+l,i^{\prime }=i+j\right)|\left(S_n,E_n,I_n\right)=\left(s,e,i\right)\right].\end{array}\end{array}$

Since the total population size should remain constant, we assume that $k\in \left\{-\mathrm{1,0}\right\}$ and $l\mathrm{,}j\in \left\{-\mathrm{1,0,1}\right\}$. Under the assumption that no other instantaneous transitions are allowed, we get:

$p_{\left(s,e,i\right);\left(s^{\prime },e^{\prime },i^{\prime }\right)}=\{\begin{array}{ll}\frac{\beta si}{N}\Delta t\hfill & \left(k,l,j\right)=\left(-1,0,1\right);\hfill \\ \delta \beta se\Delta t\hfill & \left(k,l,j\right)=\left(-1,1,0\right);\hfill \\ \alpha e\Delta t\hfill & \left(k,l,j\right)=\left(0,-1,1\right);\hfill \\ \gamma i\Delta t\hfill & \left(k,l,j\right)=\left(0,0,-1\right)\hfill \\ 1-\left(\frac{\beta si}{N}+\delta \beta se+\alpha e+\gamma i\right)\Delta t\hfill & \left(k,l,j\right)=\left(0,0,0\right);\hfill \\ 0\hfill & \text{otherwise}\text{.}\hfill \end{array}$

Figure 2. Combined SEIR model under Markov jump process.

3.4. Simulation Study

In this section we have taken up the simulation study of the COVID 19 virus. To assess the performance of our combined SEIR model we have considered the number of infected and recovered cases of COVID 19 epidemic registered up during March and April 2020 in four countries Italy, Russia, USA and Iran. The results below of the simulation study based on the CSEIR model are incorporated graphically. Statistical data for the spread of COVID 19 in Italy, Russia, USA and Iran are collected from the website of World Health Organization provided by the Department of the Ministry of Health in each country [8]. Simulation of infected individuals using CSEIR is shown in Figures 3-6 with red color. At present, the simulation of infected curves in Italy and Iran has reached their peak. Our model predicted around 110.000 infected Italian by the 20^th April 2020 and 33.000 infected Iranian by the first week of April 2020 as in Figure 3 and Figure 4. Theoretically, it is possible that Italian and Iranian authorities could control the spread of new outbreaks by quickly identifying and isolating new patients. Since the curves of infected haven’t reached yet their peaks of epidemic by April 30^th 2020, the situation in USA and Russia is less critical. Simulations using the CSEIR model predict around one million infected in USA and more than 150.000 infected in Russia by the end of April 2020 as in Figure 5 and Figure 6. Consequently, the decision to impose a preventative lockdown in USA and Russia seems perfectly reasonable.

Simulation of recovered individuals using CSEIR is shown in Figures 7-9 and Figure 10 with green color. There are two stages in the simulation of recovered individuals in Italy and Iran. The first stage is before March 20^th, 2020 during the initial outbreak of COVID-19 in both Italy and Iran; the second stage is from March 20^th to April 30^th. Our model predicted around 30.000 recovered Italian as in Figure 7, 100.000 recovered Iranian as in Figure 8, 250.000 recovered American as in Figure 9 and 20.000 recovered Russian as in Figure 10 by April 30^th 2020. The decrease of the recovered curve depends on the strength of the healthcare system and the steps taken up to stop the spreading of this virus in each country.

Figure 3. Simulation of infected in Italy using CSEIR model.

Figure 4. Simulation of infected in Iran using CSEIR model.

Figure 5. Simulation of infected in USA using CSEIR model.

Figure 6. Simulation of infected in Russia using CSEIR model.

Figure 7. Simulation of recovered in Italy using CSEIR model.

Figure 8. Simulation of recovered in Iran using CSEIR model.

Figure 9. Simulation of recovered in USA using CSEIR model.

Figure 10. Simulation of recovered in Russia using CSEIR model.

4. Conclusion

In this paper, we established a combined stochastic model through analyzing the existing data of an epidemic situation. Here, we studied the spread of COVID-19 (Infected individuals) to see whether the situation is still under control or not. Basically, if we could down the spread of the virus enough, we will be able to keep infections at a manageable level so that hospitals won’t be overcrowded. Finally, Statisticians and data experts around the globe need a constant flow of new and disaggregated data to monitor the spread of the virus. Stepping up this challenge is saving lives.

Conflicts of Interest

The author declares no conflicts of interest regarding the publication of this paper.

References

[1]	Kingman, J.F.C. (1982) The Coalescent. Stochastic Processes and Their Applications, 13, 235-248. https://doi.org/10.1016/0304-4149(82)90011-4
[2]	Tang, M.W. (2019) Fitting Stochastic Epidemic Models to Multiple Data Types. University of Washington, Washington.
[3]	Bailey, N.T.J. (1975) The Mathematical Theory of Infectious Diseases and Its Applications. Hafner Press/MacMillian Pub. Co, Sumas, Washington, Etats-Unis.
[4]	Anderson, R.M. and May, R.M. (1992) Infectious Diseases of Humans: Dynamics and Control. Vol. 28, Wiley Online Library, Hoboken.
[5]	Gillespie, D.T. (1977) Exact Stochastic Simulation of Coupled Chemical Reactions. The Journal of Physical Chemistry, 81, 2340-2361. https://doi.org/10.1021/j100540a008
[6]	Gillespie, D.T. (2000) The Chemical Langevin Equation. The Journal of Chemical Physics, 113, 297-306. https://doi.org/10.1063/1.481811
[7]	Keeling, M.J. and Ross, J.V. (2007) On Methods for Studying Stochastic Disease Dynamics. Journal of the Royal Society Interface, 5, 171-181. https://doi.org/10.1098/rsif.2007.1106
[8]	https://www.worldometers.info/coronavirus/