Abstract
In this paper, three stochastic mathematical models are developed for the spread of the coronavirus disease (COVID-19). These models take into account the known special characteristics of this disease such as the existence of infectious undetected cases and the different social and infectiousness conditions of infected people. In particular, they include a novel approach that considers the social structure, the fraction of detected cases over the real total infected cases, the influx of undetected infected people from outside the borders, as well as contact-tracing and quarantine period for travellers. Two of these models are discrete time–discrete state space models (one is simplified and the other is complete) while the third one is a continuous time–continuous state space stochastic integro-differential model obtained by a formal passing to the limit from the proposed simplified discrete model. From a numerical point of view, the particular case of Lebanon has been studied and its reported data have been used to estimate the complete discrete model parameters, which can be of interest in estimating the spread of COVID-19 in other countries. The obtained simulation results have shown a good agreement with the reported data. Moreover, a parameters’ analysis is presented in order to better understand the role of some of the parameters. This may help policy makers in deciding on different social distancing measures.
1. Introduction
After the identification of a novel strain of coronavirus in China in December 2019, later labelled COVID-19, the virus spread in China and other countries of the world. The World Health Organization declared the outbreak as a public health emergency of international concern (Roda et al., 2020) by January 2020 and the situation was designated as a pandemic by March 2020 (Raue et al., 2009). With the uncertainties about the COVID-19 virus, predictive mathematical models prove to be fundamental, from a strategic and healthcare management perspective, in order to estimate the severity of the epidemic, to forecast its time course and to minimize its socioeconomic consequences. Various models have been proposed in the literature for disease spreads. They can be categorized into agent-based models (ABMs) (Ajelli et al., 2010; Bonabeau, 2002) and compartmental models (Brauer et al., 2012; Diekmann et al., 2012; Hethcote, 2000). Deterministic compartmental models are built on differential equations and assume that the population is perfectly mixed with people moving between compartments such as susceptible (S), infected (I) and recovered (R) (Hurley et al., 2006; Jin et al., 2007; Kermack & McKendrick, 1927; Ministry of Public Health, 2020). These models revealed the threshold nature of epidemics and triumphed in explaining ‘ herd immunity’. However, they fail to capture complex social networks and the behaviour of individuals who may adapt depending on disease prevalence. On the other hand, ABMs can capture irrational behaviour and are used to simulate the interactions of autonomous agents that can be either individuals or collective entities (Epstein, 2009). In addition to their pertinence to visualization, ABM approaches in epidemiology can simulate such complex dynamic systems with less oversimplification of the rich variation among individuals. The statistical variance (in ABM) is more evident than in deterministic compartmental models, whose smooth curves often misleadingly express more certainty than justified, due to the randomization at each run (Hunter et al., 2018, 2017; Kai et al., 2020; World Health Organization, 2020a). However, the ABM approach needs realistic data, typically obtained from a census, and important assumptions and data collection in order to set their structural parameters (Epstein et al., 2007) but this is not always possible at the early stages of the outbreak. Several models have been developed for the COVID-19 pandemic. Lin et al. (2020) extended an SEIR model to account for public perception of risk and the number of cumulative cases. Anastassopoulou et al. (2020) included dead individuals in a discrete-time SIRD model and provided estimates of the main epidemiological parameters for Hubei (China). Casella (2020) derived a simplified control-oriented model comparing the outcomes of different policies and Wu et al. (2020) inferred clinical severity estimates using transmission dynamics. Giordano et al. (2020) proposed a mean-field epidemiological model extending the classical SIR model for the COVID-19 epidemic in Italy. Tracy et al. (2018) also proposed an SEIRP including compartments for asymptomatic exposed individuals and passed-away population while Goel & Sharma (2020) suggested a mobility-based SIR model for COVID-19 pandemic. Stochastic transmission models have been proposed in Hellewell et al. (2020) and Kucharski et al. (2020) and ABMs have been used for a computational simulation of the pandemic in Australia by Chang et al. (2020) and for recommending universal masking by Kai et al. (2020). Moreover, in World Health Organization (2020b), Varotsos and Krapivin developed the COVID-19 decision making system to study disease transmission. Since the list of models developed so far for COVID-19 is much longer, and the aforementioned list is non-exhaustive, we also refer the reader to Comunian et al. (2020), Calafiore et al. (2020), Calvetti et al. (2020) and Scheiner et al. (2020) for some insight on parameter calibration and inverse modelling. In the present work, we propose two stochastic discrete time–discrete state space models (one simplified and one complete) for the spread of COVID-19, in view of the daily reporting by countries of the infection indicators as well as the randomness in the transmission of COVID-19. We also obtain by a limiting process a continuous time–continuous state space mathematical model. The models take into account the characteristics of the virus by attributing to them probability distributions. We mention, for instance, the incubation period, infectiousness and testing sensitivity. Since a realistic model should take into account that there are many heterogeneities in societies that affect disease transmission, the proposed models also incorporate probability distributions for social and individual conditions such as the following: family size within the same household, number of people contacted per day and their subdivision into known versus unknown contacts. We mention in the sequel that in Britton et al. (2020), the authors showed how population heterogeneity affects herd immunity by adding to an SEIR model two features namely age structure and social activity levels of individuals. Regarding discrete modelling of COVID-19 epidemic, we refer to Boulmezaoud (2020) where a discrete-time deterministic model is proposed for forecasting the temporal evolution of the epidemic from day to day, to Li et al. (2021, 2020a) where an SEIQR difference-equation model taking into account the transmission with discrete time imported case is analysed and He et al. (2020) where a discrete stochastic compartmental model is proposed. An important discrepancy between He et al. (2020) and our work is that we use two discrete time variables one for the time of infection and the second one is for the daily update. This is crucial in modelling the COVID-19 epidemic due to the prevalence of pre-symptomatic transmission. Indeed, the distinction of two time variables permits the tracking, in the future, of people who get infected on a particular day and how long they stay actively infecting others before they get detected either by appearance of symptoms or by contact tracing combined with polymerase chain reaction (PCR) testing. We also distinguish the transmission between family members, known contacts and unknown contacts, hence including the social network in the model. The proposed continuous model is obtained by formal derivation from the proposed simplified discrete model. It consists of an integro-differential stochastic model. For instance, the SIR model may be regarded as a particular case of the obtained continuous model under several simplifying assumptions (see Remark 4.2). The remainder of the paper is laid out in the following way. Section 2 describes the main parameters that must be accounted for in modelling epidemics such as COVID-19. Section 3 presents the discrete stochastic model in two steps. First, we describe a simplified version that does not include contact-tracing and quarantine period and, second, we present the complete model. In Section 4, we formally derive a continuous stochastic integro-differential model from the simple model. Numerical simulations for the validation, the prediction and the parameter analysis of the complete discrete stochastic model are presented in Section 5. Finally, we conclude with some remarks in Section 6.
2. Description of key parameters of the model
To understand the dynamics of COVID-19 epidemic, several key parameters must be integrated in the model. In particular, the chronological characteristics of the virus, the context of contact between individuals and the influx of infected undetected individuals into the population are crucial for a model to be able to forecast the evolution of the epidemic from day to day. We shed light herein on the rationale for selecting keys parameters as follows:
- -
the times of infection and detection: in COVID-19 epidemic, many infected people may be infecting others in the community without being detected. So, it is important to reduce the delay between the infection time and the detection time (Klinkenberg et al., 2006). In the proposed model, we use, for the flow of individuals from one category to another, two time variables to distinguish infection and detection times and to keep track in the future of the people infected on a particular day. This will be clarified in the next section. The importance of such delayed effects on the modelling of COVID-19 fatality trends was illustrated by Scheiner et al. in Varotsos & Krapivin (2020), where upon the introduction of an infection-to-death delay rule to the SEIR model, the data recorded in many countries were captured better than by the traditional death kinetics law.
- -
the probability of infection (or the transmission coefficient of the virus): this parameter depends on several factors including mask wearing, adherence to hygiene measures and social distancing in public. This probability differs also depending on the nature of contact between the infected person and his/her contacts where we differentiate family members in the same household from contacts outside the house such as at work or in public transportation. The probability of infecting family members sharing the same house is assumed to be higher during the period extending from the instant of infection to the time of detection.
- -
the incubation period (time from exposure to illness onset): the general incubation period for COVID-19 ranges from 2 to 14 days depending on each individual. According to several recent studies, the average incubation time is about 5 days with the 95th percentile of the distribution at 12.5 days and the log-normal distribution provides the best fit to the data for incubation period estimates (Li et al., 2020b; Linton et al., 2020) (see Fig. 1). We use herein a discrete random variable corresponding to the incubation period obtained from the suggested log-normal distribution in Li et al. (2020b) and Linton et al. (2020).
Fig. 1.PDF for the incubation period (log-normal distribution).
- -
the serial interval (duration from onset of symptoms in an infector (a primary-case patient) to the onset of symptoms in an infectee (a secondary-case patient)): the mean serial interval for COVID-19 was estimated as 7.5 days in Li et al. (2020b) and 6.9 days in Lavezzo et al. (2020) but more recent studies (Du et al., 2020; Nåsell, 1996) suggest that it is around 4 days. Being shorter than the mean incubation period, pre-symptomatic transmission is likely and may be more frequent than symptomatic transmission (Giordano et al., 2020; Nåsell, 1996). In Nåsell (1996) and Du et al. (2020), a normal distribution was suggested for the serial interval estimates. The knowledge of the distribution of the serial interval will help us to estimate the proportion of infected people who will be detected according to the appearance of symptoms opposed to those who will be isolated according to contact tracing.
- -
the connections: these are the contacts of a person that we classify into family members, known contacts such as coworkers or neighbours and unknown contacts such as people encountered in public transportation or in social events. To account for family members, we consider a discrete random variable obtained from a log-normal probability distribution function for the average family size or household size (see Fig. 2 for the case of Lebanon). For known contacts, we consider an exponential distribution for the daily new encounters. So if a person was in contact with the individual on the first day, he/she is not counted in the new encounters of the following days. In other words, the first encounter is considered for the possible infection of the person. Finally, the unknown contacts are chosen to be uniformly distributed.
Fig. 2.PDF for the family size in Lebanon.
- -
contact-tracing: identifying the source of infection and tracing and isolating its contacts are crucial for breaking the chain of transmission and for the control of the epidemic (Hellewell et al., 2020). Two parameters influence the effectiveness of contact-tracing and isolation: the transmissibility of the pathogen (measured by the basic reproductive number) and the proportion of presymptomatic transmission (Fraser et al., 2004). Very high levels of contact tracing are required in the case of presymptomatic infectiousness (Hellewell et al., 2020). In the present model, we account for contact tracing efficiency as detailed in Section 3.
- -
the effect of border opening: in order to account for the effect of airport and border opening, we split infected individuals into travellers and locals. The probability distribution function for travellers is based on data collected from local health authorities. Furthermore, we distinguish between travellers having a positive screening result and those with a false negative result. Those with a positive test result are assumed to be isolated and are no longer infectious. The number of travellers with false negative test result is deduced from the number of travellers with positive test result based on the ratio of false negative screening results in the country. Moreover, a quarantine period is imposed on travellers and it changes according to policies adopted by local authorities. However, the rate of compliance to this measure varies from one individual to another. On the other hand, the local infected residents are subdivided into different categories: the symptomatic individuals with positive test result versus the individuals screened due to their contact with an infected person before onset of symptoms.
Some of the aforementioned parameters are calibrated to fit the real data by numerical optimization (e.g. probability of infection) or obtained from different studies and publications (e.g. incubation period, serial interval).
3. Discrete mathematical models
Daily infection being reported by countries worldwide, the unit of time adopted in the model is a single day, denoted herein by |$n$|.
3.1 Simplified model
- -
|$P:\mathbb {N}\rightarrow \mathbb {R}^+$|, the cumulative number of travellers with positive PCR result,
- -
|$N:\mathbb {N}\rightarrow \mathbb {R}^+$|, the cumulative number of infected travellers with false negative PCR,
- -
|$F:\mathbb {N}\rightarrow \mathbb {R}^+$|, the cumulative number of infected family members, and
- -
|$C:\mathbb {N}\rightarrow \mathbb {R}^+$|, the cumulative number of infected contacts.
Let |$\varDelta I:\mathbb {N}\rightarrow \mathbb {R}^+$| denote the daily number of new infected individuals and |$\varDelta P:\mathbb {N}\rightarrow \mathbb {R}^+$| the daily number of new travellers with positive PCR result which is a source term to the model and it is equal to zero when borders are closed. Moreover, we introduce the following functions: |$\varDelta N(\cdot ,\cdot )$|, |$\varDelta F(\cdot ,\cdot )$| and |$\varDelta C(\cdot ,\cdot )$| defined from |$\mathbb {N}\times \mathbb {N}$| into |$\mathbb {R}^+$|, where |$\varDelta X(n,n)$| denotes the daily new number of infected individuals from category |$X\in \{N,F,C\}$| that have been infected on day |$n$| and |$\varDelta X(k,n)$|, for |$k>n$|, denotes the remaining number out of |$\varDelta X(n,n)$| that are still active on day |$k$|. By convention, we write |$\varDelta X(0,0) = 0$|.
Schematic of the infection process between the categories |$N,F$| and |$C$|.
where |$\alpha $| and |$\beta $| are the infection probability coefficients within and outside the same household, respectively, and that will be estimated according to the real data, |$p_S(n) = \dfrac {S(n-1)}{N_{total}}$| is the proportion of susceptible individuals within the total population, |$ dF$| is the random variable of the family size whose probability distribution for the case of Lebanon is given in Figure (2) and |$dC$| is a time-dependent random variable for the number of new contacts met per day starting from the infection time. We assume that |$dC$|, as a function of time, has an exponential form that depends on both social habits and measures imposed by authorities. On the first day after the infection time, |$dC$| is maximal and it decreases to zero as time increases.
In a first step, we assume that the removal of infected individuals takes place according to the appearance of symptoms which leads to isolation or to the cure of asymptomatic cases so that they are no longer infectious. Thus, there is no contact tracing in this simplified model. Also, we do not take explicitly the quarantine period imposed on the travellers by the government.
For a better understanding of the discrete model, in the following Lemma, we solely express |$\varDelta C(n,n)$| (given in Equation (8)) in terms of |$\varDelta F(k,k)$|, |$\varDelta C(k,k)$| and |$\varDelta N(k,k)$| that can be viewed as functions of one time variable |$k$|. The resulting system can be considered well fitting into the discrete time series models. Moreover, this way of expressions turns out to be helpful in the derivation of a continuous model.
Systems (18) and (19) show that the daily new number of infected family members |$\varDelta F(n,n)$| and contacts |$\varDelta C(n,n)$| has a similar formulation to a nonlinear autoregressive model with a delayed input source term |$\varDelta N(n,n)$| (nonlinear ARX model) but whose coefficients are random variables. These coefficients vanish after some interval of time related to the maximum incubation period.
3.2 Complete model
Schematic of the infection process for each category: |$N,F,C$| and |$U$|.
For the cases that are detected after the appearance of symptoms or after healing, we use the same formulation that we developed in the simplified model. Indeed, the numbers |$\varDelta F_S(n,n)$|, |$\varDelta C_S(n,n)$|, |$\varDelta N(n,n)$| and |$\varDelta U(n,n)$| of new infected people on day |$n$| will change during the next days according to Equation (10).
According to Equations (35) and (36), and using Equations 2427, we can now express |$Z(n)$| in terms of |$Z(k)$| for |$k=1,\cdots ,n-1$|. Again, as in the simplified model, we obtain a sort of nonlinear ARX model for the numbers of daily new cases in all of the categories |$F_S,F_T,C_S,C_T$| and |$U$| where |$N$| and |$P$| are source terms or delayed inputs.
4. Continuous stochastic model
In the previous section, the discrete time variable was denoted |$n$| and the unit time step was one day. Now, we denote by |$t$| the continuous time variable and we use the same notations for the functions |$F,C,P,N,S,I,\cdots $| that are considered functions of the time variable |$t$|. The functions |$P$| and |$N$| are known inputs to the model.
In this section, we will present a formal passing to the limit in the discrete simplified model in order to derive an equivalent continuous model.
If in addition to the assumptions of Remark 4.1, we do not distinguish between categories |$F$| and |$C$|, and we do not take into consideration the travellers, then |$F$|, |$N$| and |$P$| are null. In this case, the number of total infected cases is |$I(t) = C(t)$| and the number of active infected cases is |$I_a(t) = I(t)-R(t)$|, and |$R^{\prime}(t) = RC^{\prime}(t)$|. Also, we have |$S(t)+I_a(t)+R(t) = S(t) + C(t) = N_{total}$|.
5. Parameters calibration and numerical simulations
In this section, we use the complete discrete model first in a general context with different values of the parameters. In all simulations, we use, for the family size variable, the distribution given in Fig. 2 and an exponential distribution for the daily new known contacts. The model being stochastic, we use 1000 runs for each simulation and we consider the average as a representative for these 1000 runs. In some figures, we display the average of the 1000 runs along with the corresponding standard deviation. The code is implemented with MATLAB and is publicly available on the link: https://github.com/AymanMourad/Stochastic-Corona-Model.
5.1 A sample population
We consider, herein, a fully susceptible population of |$6,000,000$| individuals (i.e. |$S=N_{total}=6,000,000$|) and we assume that on day 1, 10 infected travellers were introduced (i.e. |$\varDelta N(1,1)=10$|) and then the borders are closed so that no more infectious individuals are introduced into the population. We use the complete model with the parameters as given in Table 1. To show the effect of contact tracing, we used in one simulation |$p=0.1$| and in the other |$p=0.3$|, so that in the first there is |$90\%$| efficiency of contact tracing while in the second it is |$70\%$|. In Fig. 5, we see that with no measures imposed and no isolation of the infected travellers, the number of new infected individuals reaches the peak on day 68, when |$p=0.3$| and the total infected population on day 100 counts 2,435,000 individuals. On the other hand, when |$p=0.1$|, the peak is reached on day 80 and the total count reaches 1,710,000 individuals after 100 days. A plateau is reached in the total number of infections by that time suggesting that herd immunity has started going into effect (for the parameters given in Table 1). Also in another set of simulations, we fixed |$p=0.1$| and we varied the values of |$\beta $| and |$\gamma $|. The results are shown in Fig. 5.
Parameters of the sample simulation.
| Simulation | |$\textbf {p}$| | |${\boldsymbol \alpha }$| | |${\boldsymbol \beta }$| | |${\boldsymbol \gamma }$| | |$\textbf {q}_{max}$| |
|---|---|---|---|---|---|
| I | 0.1 | 0.15 | 0.045 | 0.045 | 1 |
| II | 0.3 | 0.15 | 0.045 | 0.045 | 1 |
| III | 0.1 | 0.15 | 0.042 | 0.042 | 1 |
| IV | 0.1 | 0.15 | 0.042 | 0.038 | 1 |
| V | 0.1 | 0.15 | 0.041 | 0.042 | 1 |
| Simulation | |$\textbf {p}$| | |${\boldsymbol \alpha }$| | |${\boldsymbol \beta }$| | |${\boldsymbol \gamma }$| | |$\textbf {q}_{max}$| |
|---|---|---|---|---|---|
| I | 0.1 | 0.15 | 0.045 | 0.045 | 1 |
| II | 0.3 | 0.15 | 0.045 | 0.045 | 1 |
| III | 0.1 | 0.15 | 0.042 | 0.042 | 1 |
| IV | 0.1 | 0.15 | 0.042 | 0.038 | 1 |
| V | 0.1 | 0.15 | 0.041 | 0.042 | 1 |
Parameters of the sample simulation.
| Simulation | |$\textbf {p}$| | |${\boldsymbol \alpha }$| | |${\boldsymbol \beta }$| | |${\boldsymbol \gamma }$| | |$\textbf {q}_{max}$| |
|---|---|---|---|---|---|
| I | 0.1 | 0.15 | 0.045 | 0.045 | 1 |
| II | 0.3 | 0.15 | 0.045 | 0.045 | 1 |
| III | 0.1 | 0.15 | 0.042 | 0.042 | 1 |
| IV | 0.1 | 0.15 | 0.042 | 0.038 | 1 |
| V | 0.1 | 0.15 | 0.041 | 0.042 | 1 |
| Simulation | |$\textbf {p}$| | |${\boldsymbol \alpha }$| | |${\boldsymbol \beta }$| | |${\boldsymbol \gamma }$| | |$\textbf {q}_{max}$| |
|---|---|---|---|---|---|
| I | 0.1 | 0.15 | 0.045 | 0.045 | 1 |
| II | 0.3 | 0.15 | 0.045 | 0.045 | 1 |
| III | 0.1 | 0.15 | 0.042 | 0.042 | 1 |
| IV | 0.1 | 0.15 | 0.042 | 0.038 | 1 |
| V | 0.1 | 0.15 | 0.041 | 0.042 | 1 |
Sample simulation results: left (I and II); right (III, IV, V). The daily infected individuals correspond to |$\varDelta I(d)$|, the daily detected |$\varDelta R(d)$| and the total detected |$R(d)$|. The last row shows the realizations of 1000 runs of scenario (II) and the histogram of the total number of removed cases at the final time of each realization of the simulations corresponding to scenario (V).
Figure 5 clearly shows that the higher the value of |$p$| (i.e contact tracing will be more effective), the lower the number of daily new cases. On the other hand, we see that when |$\beta $| or |$\gamma $| are reduced, the number of infected people will reduce. But we notice that reducing |$\beta $| a little will have an important impact on the number of infected cases while this is not the case when reducing |$\gamma $| a little. This remark makes sense because from the reported data, it is known that a high proportion of the infected people are infected from known contacts such as family members, friends, coworkers and so on. This issue will be reconsidered later when studying the real case of Lebanon. When looking at the results of the 1000 runs for each scenario, we noticed that, unlike the result of Barbour & Reinert (2013), the obtained curves are not time shifts of a deterministic path. This emphasizes the role of the stochastic variables introduced in our model. For instance, in the last row of Fig. 5, we see that the limit value (plateau phase) of the total number of removed cases changes from one run to another. Moreover, the distribution of these limit values can be depicted from the histogram in the bottom right figure of Fig. 5 as having a normal distribution. This backs up the use of the average as a representative of the runs.
5.2 The case of Lebanon
The first confirmed case of COVID-19 in Lebanon was detected on 21 February 2020 and by 29 February, there was a total of 7 confirmed cases going back from travel trips. Consequently, educational institutions were closed, and gradual measures were introduced until the declaration of general mobilization and state of emergency by mid-March where public gatherings were banned, cultural venues closed and social distancing measures imposed in public. Furthermore, the airport, the land borders and the seaports were closed as of 19 March. On 5 April, new arrivals were intermittently allowed through the airport. By the end of May, authorities in Lebanon have been gradually easing restrictions. Public transportation has resumed, with social-distancing measures. Government institutions and certain private companies, including various shops and stores, were permitted to return to normal operations from June 1 and the lockdown was lifted on 7 June. Starting 1 July, flights were resumed and the airport started operating up to 10 |$\%$| of its capacity (Nishiura et al., 2020). Travelers with negative PCR results upon arrival were isolated for a maximum of 3 days. However, people were not abiding by the preventive measures, and travellers were not respecting the isolation period (Houssari, 2020). The virus rebounded to reach a total number of 4205 by 29 July and the authorities reinstated lock-down from 30 July to 3 August and from 6 August to 10 August. However, due to the explosion in Beirut, the sanitary situation deteriorated in the country and the total number of cases reached 16,870 by 31 August. In this context, it is obvious that there is an urgent need for evidence-based decisions.
Calibration and parameters’ estimation
In order to obtain a realistic model for the case of Lebanon, first we used part of the real data to estimate the parameters that are involved in the model. Then we validated the model for the second part of the real data. Indeed, using numerical optimization techniques, we calibrated the model and computed the parameters in order to reproduce results that are close to real data for the period before August 31, published on the website of the ministry of public health of Lebanon (https://coronanews-lb.com/). During the period preceding 30 June, strict measures were imposed on travellers visiting Lebanon, and the number of social events (such as weddings) was still small. On the other hand, as of 1 July, measures were relaxed and people were less respectful to the measures so the average number of new people met per day was increased, the quarantine period for travellers was restricted to 3 days and the percentage of travellers respecting the quarantine duration decreased.
For the incubation period, we used a discrete distribution obtained from Li et al. (2020b) and Linton et al. (2020) as in Fig. 1. On the other hand, due to the lack of census in Lebanon, we used the reports Masri (2008) and Yaacoub & Badre (2012) to deduce a distribution for family size in the same household (see Fig. 2). In the absence of official reports, the distribution of new people met per day was also deduced from our familiarity with the Lebanese context where public transportation is not prevalent and the number of family members met per day is high. The estimated parameters are displayed in Table 2. Those parameters change depending on the period of time that corresponds to different measures in the country. For the whole period, the probability of infection for the same household was estimated as |$\alpha =0.15$| and the probability |$p$| related to contact tracing is |$p=0.1$|. Although during model calibration, different combinations of model parameters can show good fit to the data. This is known in the literature as non-identifiability of parameters (Lintusaari et al., 2016; Rouabah et al., 2020; Sameni, 2020). In our calculations, the parameters estimation was consistent with government’s measures and people response. For instance, at the start of the epidemic, when the measures were strictly imposed, the infection coefficients were small (see Table 2). Later, people’s responsiveness to the measures decreased and many were not respecting them anymore. So measures were only announced without forcing them by the authorities, for this scenario the infection coefficients were high. The result of the simulation along with the real data is shown in Fig. 6.
Parameters for Lebanon.
| Date | |${\boldsymbol \beta }$| | |${\boldsymbol \gamma }$| | |$\textbf {q}_{max}$| |
|---|---|---|---|
| 21 Feb to 8 Mar | 0.031 | 0.020 | 14 |
| 9 Mar to 15 Mar | 0.09 | 0.03 | 1 |
| 16 Mar to 20 Apr | 0.031 | 0.020 | 14 |
| 21 Apr to 30 June | 0.035 | 0.027 | 10 |
| 1 July to 31 July | 0.038 | 0.037 | 8 |
| 1 Aug to 14 Aug | 0.039 | 0.035 | 10 |
| 15 Aug to 31 Aug | 0.037 | 0.030 | 10 |
| Date | |${\boldsymbol \beta }$| | |${\boldsymbol \gamma }$| | |$\textbf {q}_{max}$| |
|---|---|---|---|
| 21 Feb to 8 Mar | 0.031 | 0.020 | 14 |
| 9 Mar to 15 Mar | 0.09 | 0.03 | 1 |
| 16 Mar to 20 Apr | 0.031 | 0.020 | 14 |
| 21 Apr to 30 June | 0.035 | 0.027 | 10 |
| 1 July to 31 July | 0.038 | 0.037 | 8 |
| 1 Aug to 14 Aug | 0.039 | 0.035 | 10 |
| 15 Aug to 31 Aug | 0.037 | 0.030 | 10 |
Parameters for Lebanon.
| Date | |${\boldsymbol \beta }$| | |${\boldsymbol \gamma }$| | |$\textbf {q}_{max}$| |
|---|---|---|---|
| 21 Feb to 8 Mar | 0.031 | 0.020 | 14 |
| 9 Mar to 15 Mar | 0.09 | 0.03 | 1 |
| 16 Mar to 20 Apr | 0.031 | 0.020 | 14 |
| 21 Apr to 30 June | 0.035 | 0.027 | 10 |
| 1 July to 31 July | 0.038 | 0.037 | 8 |
| 1 Aug to 14 Aug | 0.039 | 0.035 | 10 |
| 15 Aug to 31 Aug | 0.037 | 0.030 | 10 |
| Date | |${\boldsymbol \beta }$| | |${\boldsymbol \gamma }$| | |$\textbf {q}_{max}$| |
|---|---|---|---|
| 21 Feb to 8 Mar | 0.031 | 0.020 | 14 |
| 9 Mar to 15 Mar | 0.09 | 0.03 | 1 |
| 16 Mar to 20 Apr | 0.031 | 0.020 | 14 |
| 21 Apr to 30 June | 0.035 | 0.027 | 10 |
| 1 July to 31 July | 0.038 | 0.037 | 8 |
| 1 Aug to 14 Aug | 0.039 | 0.035 | 10 |
| 15 Aug to 31 Aug | 0.037 | 0.030 | 10 |
Real data of Lebanon and average simulation results with parameters as in Table 2.
Validation and prediction
In order to validate the estimated parameters in the previous paragraph, we compare the cumulative number of cases and daily new cases in the period extending from 1 September to 22 November, as forecast by the ministry of public health, to the results obtained in the simulation. For this period of time, new measures of partial lock-down in some cities were intermittently applied and a total lock-down was announced as of 14 November. Accordingly, the values of |$\beta $| and |$\gamma $| were reduced by |$1\%$| for this period compared to the period of the last week of August. The results of the simulations that are shown in Fig. 7 greatly overlap with the real data. The discrepancy is attributed to the fact that the total number of daily infected people is not fully reported as shown in Fig. 7. Indeed, many asymptomatic people are never detected because contact tracing is not |$100\%$| efficient while others may develop minor to mild symptoms but do not get tested for several reasons (some of them are financial). This highlights the importance of wide testing in containing the epidemic, also. Indeed, looking at the curve of the active infected people shown in Fig. 7, we can see the number of infected cases that are still infectious in the society. So, with a broad testing financed by the government, the number of active infected individuals will be reduced. Hence, the spread of the pathogen can be better controlled.
From top to bottom: Total number of cases till November 20 compared to real data - Daily new reported cases (|$\varDelta R$|) till April 30 with real data till November 20, in addition to the active infected individuals (|$I_a$|) obtained by simulation - Prediction of the total number of cases till April 30 with the parameters |$\beta $| and |$\gamma $| reduced by |$1\%$| for the whole period from September 1 till April 30.
Parameter analysis
For a deeper insight on the impacts of the main parameters involved in our models, we present in this paragraph a numerical parameter analysis. Indeed, we study by numerical simulations the effect of the change of the parameters for the case study of Lebanon. Mainly, we focus here on the role of the two main parameters |$\beta $| and |$\gamma $| that model the infection coefficients for known and unknown contacts, respectively, since these are the parameters that are mostly tied to measures. For instance, the parameter |$\alpha $| that models the infection coefficient between household members is conceived to be constant and unaffected by measures imposed by governments whereas |$\beta $| and |$\gamma $| can be highly reduced when lockdown is imposed in the country. Accordingly, different scenarios associated to the modifications of |$\beta $| and |$\gamma $| are simulated. First, in Scenario I, we reduce each coefficient by |$1\%$| with respect to the values used in the period of the last two weeks of August as shown in Table 2 of Section 5.2. Then in Scenario II we reduced |$\beta $| and |$\gamma $| by |$2\%$|, in Scenario III we reduced |$\beta $| by |$1\%$| and |$\gamma $| by |$2\%$| and in Scenario IV we reduced |$\beta $| by |$1\%$| and |$\gamma $| by |$10\%$|. The numerical results are illustrated in Fig. 8. They show that a slight modification of the parameter |$\beta $| has an important change on the number of infected people; however, this is not the case for the parameter |$\gamma $| where a big change is required in order to get a significant change in the result. We understand from these simulations that, for the case of Lebanon, the parameter |$\beta $| is the most effective in the control of the spread of the epidemic. Since |$\beta $| is the infection coefficient for known contacts, this result suggests that the epidemic can be efficiently controlled when strict measures are imposed in the context of workplaces, universities, schools, family gatherings or shops in villages where people know each other. On the other hand, if such measures are not taken in those places, authorities must impose much stricter measures on public transportation, malls, and restaurants in cities (to name a few places where people get in contact with strangers). It is clear by comparing Scenarios II and IV that to get a similar effect to a reduction of both |$\beta $| and |$\gamma $| by |$2\%$|, |$\gamma $| must be reduced by |$10\%$| if |$\beta $| is only reduced by |$1\%$|. We can also better understand this observation when we look at the number of infected people in each of the categories |$F$|, |$C$| and |$U$|. For instance, in Fig. 9, we show the corresponding results for the Scenarios (I) and (II). The curves of the number of infected cases in the category |$C$| are the highest; however, they are the lowest for the category |$U$|. Moreover, we see how reducing the parameter |$\beta $| affects the number of infected cases in the category |$F$| that corresponds to the household members. These observations hold for the case of Lebanon and these are not necessarily valid for the case of other countries where the category |$U$| may have an important role in transmitting the virus mainly when for instance public transportation is essential in the society.
Average simulation results for 4 scenarios with |$\beta $| reduced by 0.01 or 0.02 while |$\gamma $| is reduced either by 0.01 or by 0.02 or 0.1.
Cumulative number of infected individuals from categories F, C and U corresponding to Scenarios I and II from November 16 till April 30.
Parameters for the simulations of the epidemic in Lebanon during the period extending from 20 November to 30 April and average of total number of cases by 30 April 2021.
| Scenario | |${\boldsymbol \beta }$| | |${\boldsymbol \gamma }$| | Average |$\pm $| Std |
|---|---|---|---|
| I | 0.03663 | 0.0297 | 289,642 |$\pm $| 62,037 |
| II | 0.03626 | 0.0294 | 221,727 |$\pm $| 54,737 |
| III | 0.03663 | 0.0294 | 283,628 |$\pm $| 61,596 |
| IV | 0.03663 | 0.027 | 238,281 |$\pm $| 57,189 |
| Scenario | |${\boldsymbol \beta }$| | |${\boldsymbol \gamma }$| | Average |$\pm $| Std |
|---|---|---|---|
| I | 0.03663 | 0.0297 | 289,642 |$\pm $| 62,037 |
| II | 0.03626 | 0.0294 | 221,727 |$\pm $| 54,737 |
| III | 0.03663 | 0.0294 | 283,628 |$\pm $| 61,596 |
| IV | 0.03663 | 0.027 | 238,281 |$\pm $| 57,189 |
Parameters for the simulations of the epidemic in Lebanon during the period extending from 20 November to 30 April and average of total number of cases by 30 April 2021.
| Scenario | |${\boldsymbol \beta }$| | |${\boldsymbol \gamma }$| | Average |$\pm $| Std |
|---|---|---|---|
| I | 0.03663 | 0.0297 | 289,642 |$\pm $| 62,037 |
| II | 0.03626 | 0.0294 | 221,727 |$\pm $| 54,737 |
| III | 0.03663 | 0.0294 | 283,628 |$\pm $| 61,596 |
| IV | 0.03663 | 0.027 | 238,281 |$\pm $| 57,189 |
| Scenario | |${\boldsymbol \beta }$| | |${\boldsymbol \gamma }$| | Average |$\pm $| Std |
|---|---|---|---|
| I | 0.03663 | 0.0297 | 289,642 |$\pm $| 62,037 |
| II | 0.03626 | 0.0294 | 221,727 |$\pm $| 54,737 |
| III | 0.03663 | 0.0294 | 283,628 |$\pm $| 61,596 |
| IV | 0.03663 | 0.027 | 238,281 |$\pm $| 57,189 |
6. Conclusion
In this work, we proposed three stochastic models (two discrete and one continuous) for the spread of COVID-19. These models may be generalized to any other epidemic by considering appropriate characteristics of the pathogen. In these models, we include several social indicators related to family structure and habits by splitting the infected individuals into several categories and we take into account contact-tracing and quarantine period. We used two discrete time variables to permit tracing infected people from the time of infection until the time of detection. We were able to express the daily new cases at a given day in terms of the daily new cases of previous days. The resulting model is similar to a nonlinear autoregressive discrete time series model with delayed input corresponding to the travellers’ inflow to the country. By formal passage to the limit, we obtained a continuous stochastic integro-differential system. Under several simplifying assumptions, the proposed continuous model reduces to an SIR-like system. Furthermore, we validated the efficiency of the complete discrete model using the available data for the case of Lebanon. Indeed, we were able to reproduce the pattern of this data by estimating the parameters involved in the model using optimization techniques. As for the prediction, we considered four scenarios for the future period of time where each scenario corresponds to a set of parameters obtained for a previous time slot under particular mitigation measures imposed in the country. Although the proposed model was applied to the COVID-19 epidemic in Lebanon, the same approach may be applied to other countries by changing the distributions of household size, known contacts and unknown contacts. Moreover, one may add the spatial distribution of the population to the continuous integro-differential model and take into account the effect of mobility between different regions on the spread of the disease. In summary, we consider that this approach constitutes one step in building models that take into account societal factors and complex human behaviour without an extensive process of data collection.
References
