key: cord-0795260-s736xibz
authors: Đorđević, J.; Papić, I.; Šuvak, N.
title: A two diffusion stochastic model for the spread of the new corona virus SARS-CoV-2
date: 2021-04-30
journal: Chaos Solitons Fractals
DOI: 10.1016/j.chaos.2021.110991
sha: 5b6516fd0b6d66fd817d5196ef91a0905e24f8a3
doc_id: 795260
cord_uid: s736xibz

We propose a refined version of the stochastic SEIR model for epidemic of the new corona virus SARS-Cov-2, causing the COVID-19 disease, taking into account the spread of the virus due to the regular infected individuals (transmission coefficient [Formula: see text]), hospitalized individuals (transmission coefficient [Formula: see text] [Formula: see text]) and superspreaders (transmission coefficient [Formula: see text]). The model is constructed from the corresponding ordinary differential model by introducing two independent environmental white noises in transmission coefficients for above mentioned classes - one noise for infected and hospitalized individuals and the other for superspreaders. Therefore, the model is defined as a system of stochastic differential equations driven by two independent standard Brownian motions. Existence and uniqueness of the global positive solution is proven, and conditions under which extinction and persistence in mean hold are given. The theoretical results are illustrated via numerical simulations.

The COVID-19 disease, caused by the novel coronavirus SARS-CoV-2 (severe acute respiratory syndrome coronavirus type 2), today is in the focus of research in many scientific disciplines. In comparison to most of the other coronaviruses, the SARS-CoV-2 virus, discoveblack in December 2019 and believed to be of animal origin (see e.g. [1] ), causes the severe respiratory infection with a wide range of different symptoms, consequences and, in some cases, even death. In comparison to influenza in which one infected person spreads the virus to 1 − 2 others, in the case of SARS-CoV-2 virus it is more likely that one infected person infects 2 − 4 other persons, according to the Norwegian Institute of Public Health [2] . The virulency of this virus is highly dependent on its robustness which is, with respect to the different environmental conditions, studied in [3] . According to this study, after 14 days of incubation in the virus transport medium at 4 • C, the virus was stable. However, if the incubation temperature was increased to 70 • C, the virus was inactivated within 5 minutes. Regarding different surfaces inoculated with the virus at room temperature of 22 • C, it was shown that the printing and tissue paper was free of the infectious virus after 3 hours of incubation, while the virus is much more stable on smooth surfaces -after 4 days there was no infectious virus on glass surfaces and banknotes, while the stainless steel and plastic were free of infectious virus 7 days after the exposure. Furthermore, this study claims that various disinfectants at room temperature of 22 • C disable infectious virus within 5 minutes and that virus is highly stable in pH 3 − 10 environments.

However, the main transmission way is the close person-to-person contact during which the virus transmits via aerosols and respiratory droplets (see [4] ). According to [5] , this virus has the highest infectiousness potential just before or within the first five days after the appearance of symptoms. Although the positive reverse transcription polymerase chain reaction (RT-PCR) test is not equivalent to infectiousness of a person, it is important to point out that these tests detect the viral RNA in the upper respiratory tract of an infected person for a mean of 17 days (for more details on duration of incubation period and mean incubation time, based on different sources, see e.g. [6] and [7] ). In some cases the symptoms of the infection doesn't occur and these asymptomatic carriers silently spread the infection (see [4] ). However, according to [5] , pre-symptomatic (1 − 2 days before the symptom onset) and symptomatic individuals play a greater role in the spread of the infection than the asymptomatic carriers. For more information on asymptomatic carriers and the estimations of their proportion within the infected population we refer to [8] and [9] . Another category of virus carriers that plays an important role in spreading the disease are the so-called superspreaders. Superspreader status of an infected individual includes one or more of the following factors: high viral load due the immunity issues, underlying diseases, existing infectious co-factors and/or elevated social activity. Early modeling of COVID-19 data suggest a possibility that a small proportion of infected individuals could be responsible for high proportion of transmissions, indicating the importance of the role of the superspreaders in SARS-CoV-2 spread (see [10] and references therein). For more epidemiological details and references on this virus and the related disease we refer to [11] , [12] and [13] .

This pandemic triggeblack the adaptation of some classical and development of many new mathematical models for disease spread, and promoted the application of such models in biomedical research. One of the novel approaches, relying on recent successful modeling of different phenomena in online social dynamics, is the agent-based modeling approach, which is used for simulation of infection transmission driven by the social activity dynamics and comprising e.g. the characteristics of individuals as well as the virus survival time and its mutations, see e.g. [14] , [15] and [16] . Furthermore, we point out the application of machine learning algorithms within the standard time-series models (see e.g. [17] , [18] , [19] , [20] ) and use of deterministic SIR, SEIR and similar classical compartmental epidemiological models and many of their refinements (for different refinements of the SEIR model we refer to [21] , [22] , [23] and [24] ). As the spread of the SARS-CoV-2 virus is highly dependent of many parameters which are stochastic in nature (e.g. transmission rate and contact rate), this motivated the transition of deterministic models to their stochastic counterparts (see e.g. [25] ). For examples of stochastic models, quite different than the models in the present paper, we refer to the binomial distribution based discrete stochastic model in [26] , the stochastic version of the SEIR model in [27] and the refined SEIR model in [28] including stratification by age, mobility, social contacts and the impact of the asymptomatic cases. For another model incorporating spatial heterogeneity and stochasticity we refer to [29] .

To get more close to the approach used in this paper, we point out that due to environmental uncertainty and the behavior of infectious agents (e.g. viruses) in such environment, the classical literature suggests that the purely differential population dynamics models benefit from introducing the stochastic noise, see e.g. [30] , [31] , [32] , [33] and [34] . There are many ways of introducing environmental uncertainty into continuous-time epidemiological models, but most of them are concentrated on modeling the transmission rate β, which is in classical deterministic models either constant or given as a deterministic function of time t (e.g. piecewise constant function as in [35] ). Here we outline several main lines of introducing stochastic noise into the continuous-time epidemiological model. First line follows the idea of introducing uncertainty by modeling the transmission rate β by an Itô diffusion process

with drift µ(·) and diffusion σ(·), possibly depending on β(t) (see e.g. [36] where the Ornstein-Uhlenbeck process is incorporated in the classical SIS compartmental model). Introducing other Itô diffusions as models for transmission rate results in analytically complex stochastic systems depending on Itô integrals of unknown distribution which then has to be either approximated or estimated and therefore allow numerical rather than stochastic analysis. For the problem of estimation of the distribution of certain Itô integral observed by the authors in the case of modeling the transmission rates with certain heavy-tailed diffusions we refer to [37] and [38] . Furthermore, the transmission rate and possibly other model parameters can be modeled by a more general SDEs, comprising Lévy driven jumps beside the driving Brownian motion (see e.g. [39] and [40] ). Second line, referring to the more simple approach for introducing the stochastic noise into continuous-time compartmental models, used in e.g. [41] , [42] , [43] and [44] , is the perturbation of the constant transmission rate β by the additive noise

where (B(t), t ≥ 0) is the standard Brownian motion with a certain intensity. This approach, due to additive nature of the noise, results in system of Itô SDEs which are more suitable for explicit stochastic analysis than the systems incorporating Itô diffusions as models for transmission rate:

• The existence and uniqueness of the positive local solution follows from the classical stochastic analysis book [45] , while the extension to its global version relies on proving that the explosion time for the solution is almost surely (a.s.) infinite. The later problem, as well as the problem dealing with the derivation of the explicit conditions for the persistence of the disease, are solved by applying the transformation of the stochastic system by the appropriately defined Lyapunov function (see e.g. [46] and [47] ), the classical multidimensional Itô formula (see e.g. [45] ) and elementary stochastic estimates (inequalities).

• Explicit conditions for extinction of the disease, i.e. the condition on the intensities of the stochastic noises, are derived by applying the Itô formula (see [45] ) and stochastic inequalities on a suitable chosen equation from the system.

• All results and underlying theoretical conditions are illustrated by simulation for the appropriately chosen parameter values, see Figure ( 2) for the problem of persistence and Figure ( 3) for the problem of extinction of the disease.

Third line relies on the assumption of varying infectivity, where either the transmision coefficient or the contact rate are taken to be a simple piecewise defined deterministic functions or some rightcontinuous random functions, while in the system itself the distribution of the incubation time and the distribution of the period of infectiousness are used as a kind of a "distributional delay", see e.g. [35] , [48] and [49] . Delayed stochastic models are appealing since they comprise the memory of the process. Beside [35] , for different approaches of introducing delays in the epidemiological models we refer to [50] and [51] . This paper is structublack as follows. After extensive Introduction, In Section 2 (Mathematical models and methodology) we define deterministic and stochastic refinement of the SEIR model and give a short overview of the methodology used for its analysis. In section 3 (Main results), the existence and uniqueness of the positive global solution of this system are proved (subsection 3.1), the sufficient conditions for persistence and extinction of the disease are given in Theorems 2 and 3 (subsections 3.2 and 3.3, respectively) while subsection 3.4 is dedicated to simulations for reasonably chosen parameter values. Section 4 is dedicated to discussion and conclusion, comprising the guidelines for future research.

Torres et al. in [21] defined the deterministic model of SEIR type in which they introduced the classes of asymptomatic individuals A, hospitalized individuals H and superspreaders P (SEIPHAR model).

The model depends on several parameters, including the constant transmission coefficients related to symptomatic infected individuals (β), hospitalized individuals (lβ) and to superspreaders (β ). The importance of refinement of the infected class by introducing the classes of asymptomatic carriers, hospitalized individuals and superspreaders is visible in the corresponding figures within the simulation study at the end of the Section 3.

Naturally, the assumption of constant transmission coefficients is rather restrictive. Therefore, motivated by the nature of transmission rates itself as well as by [34] and other references from the Introduction discussing the stochastic nature of the model coefficients, we introduce into the model two independent stochastic noises by the additive perturbation of constant transmission coefficients β and β (for more details please see (3), (4) and (5)).

Beside this novelty, we generalize the particular deterministic model from [21] by introducing more refined parametrization. In particular, we introduced different recovery rates for asymptomatic carriers, symptomatic individuals and superspreaders.

Furthermore, we would like to emphasize the advantages of the stochastic model over its deterministic counterpart in short-term pblackictions as well as to show that for the long-term behaviour the deterministic and the stochastic models greatly coincide, which is clearly visible from the simulation study in subsection 3.4.

We observe the SEIPHAR compartmental model for epidemics of the new corona virus SARS-CoV-2, introduced in the recent paper [21] . Under this model, the human population is divided into seven mutually exclusive compartments:

• S -susceptible individuals,

• E -individuals exposed to the virus SARS-CoV-2, but not yet infectious to others (they may become infectious after a certain incubation period),

• I -symptomatic infectious individuals,

• P -superspreaders,

• A -asymptomatic infectious individuals,

• H -hospitalized individuals,

• R -recoveblack individuals.

Let us suppose that the total population size at time t is given by

We introduce the difference in the intensity of the human-to-human transmission (per day per person) between classes I, P and H. The transmission coefficient describing the human-to-human transmission due the regular infected individuals is denoted by β, while the transmission due to the superspreaders (who are more likely to infect others, compablack with a typical infected person) is described by the coefficient β . The relative infectiousness of the hospitalized individuals is described by the coefficient lβ, with the natural assumption β < β , where l is an estimated constant. An individual leaves the exposed class and becomes infectious (symptomatic, super-spreader or asymptomatic) at rate κ. The proportions of exposed individuals progressing to symptomatic, super-spreader or asymptomatic class are ρ 1 , ρ 2 and (1−ρ 1 −ρ 2 ), respectively. Symptomatic individuals and superspreaders become hospitalized at the average rate γ a , hospitalized individuals are recoveblack at rate γ r , while non-hospitalized individuals recover with rate γ i . The disease induced death rates are δ i , δ p and δ h for infected individuals, superspreaders and hospitalized patients, respectively, while death-rate in other classes is µ. The description of dynamics of asymptomatic infectious individuals within this model is limited, due to the lack of real measurements regarding this class. For more information on the role of asymptomatic carriers we refer to the recent paper by [52] . For estimated parameter values in similar model, based on the data from e.g. [53] and similar national resources, we refer to [22] . Furthermore, Table 1 in this paper contains sensible parameter values chosen for the simulation study.

Remark 2.1 It is sensible to assume that these parameters depend on the age a and location m of individuals and therefore they could be treated as functions of these two variables, allowing the application of this model to some age and/or location restricted subpopulation, for example κ = κ(a, m) (for similar approach see [28] ).

The dynamics of the SARS-CoV-2 virus spread in the framework of the model constructed in [21] is described by the following system of ordinary differential equations (ODE):

It follows that the number of disease induced deaths at time t is given by

(2)

In this subsection we introduce the stochastic version of the model (1). The model is constructed by introducing the perturbation in the form of the environmental white noise in transmission coefficients β and β :

where

are independent standard Brownian motions with intensities σ 1 and σ 2 , respectively. However, stochastic model under consideration is not only stochastic generalisation of the deterministic model (1) via introduced noises, but it is also qualitative generalisation of the model since we introduced parameter Λ (the expected number of new susceptible individuals, observed in appropriate units, e.g. per year per million people), general population death rate µ and finer assessment in terms of recovery rates. The resulting system of stochastic differential equations (SDEs) is of the following form:

where k 1 , k 2 are constants such that k 2 < k 1 < 1. Remark 2.2 Observe that the last two equations, regarding the dynamics of asymptomatic infectious patients and recoveblack individuals, are independent of the first five equations. Therefore, it is enough to analyze the problem of existence and uniqueness of the global positive solution of the system of first five SDEs from system (5), since it applies to the complete system (5) as well.

Schematically, the transitions in the framework of the stochastic model are shown in Figure 1 (note that under β and β in this scheme we understand their perturbed versions). Figure 1 : Scheme of the SDE system (5)

For the stochastic model (5) we observe three usual, but important problems, which are analyzed by standard tools from stochastic calculus.

First we prove existence and uniqueness of the positive global solution of stochastic system (5). This problem for a general system of Itô SDEs is widely coveblack in the classical and more recent mathematical literature. For example, in [30] it is proved that a presence of even a tiny amount of environmental noise in population dynamics system can suppress its explosion, i.e. even if solution of the deterministic system explodes to infinity in a finite time, under certain conditions with probability one the solution of the corresponding stochastic system does not (for more details please see [30] and references therein). Underling idea of the proof, also used in our paper, is that the explosion time is a.s. infinite, which provides the extension of the positive local to the the positive global solution. This statement for the model (5) is proven in Subsection 3.1.

Other two problems addressed in the paper are derivation of explicit conditions, comprising the intensities of Borwnian motions B 1 and B 2 , under which the disease is persistent or goes into extinction. We deal with this problems in Subsections 3.2 and 3.3, respectively. Motivation for the methodology used comes from the paper [34] . In particular, persistence is observed in the mean sense, which understands the expected value of the process of all infected individuals over the interval 0 to t (integral over time of the sum of processes regarding to all infected subgroups, divided by the length of the time interval). Conditions under which this integral will be greater than a certain positive constant as time tends to infinity are given in Theorem 2, meaning that under these conditions the disease remains persistent in the population. On the contrary, for the extinction of the disease, we generate the conditions under which the disease vanishes, i.e. the conditions under which the value of the sum of processes of all infected subgroups tends to zero as time tends to infinity. The sufficient conditions for the extinction of the COVID-19 disease are given in Theorem 3.

To emphasize more novelties of this paper that justify its position in the recent vast number of publications regarding the COVID-19 pandemic, we refer to two recent papers [43] and [44] , dealing with similar, but not the same refinements of the SEIR model of the same stochastic nature as the model (5).

The paper [43] deals with the SLIRD (susceptible -infected -latent -recoveblack -dead) model for which the problems of existence and uniqueness of the positive solution as well as the problems of persistence and the extinction of the disease are addressed. However, the advantage of the model (5) is the observation of additional classes of superspreaders, asymptomatic and hospitalized patients.

Very recent paper [44] deals with the SEIAQHR (susceptible -exposed -infected -asymptomatic -quarantined -hospitalized -recoveblack) model, where the class I contains infected individuals with symptoms and the class A infected individuals that are asymptomatic but can become symptomatic over time. Although very similar to our model (5) , this model treats asymptomatic carriers in a different way and doesn't cover the superspreader class. Regarding the mathematical analysis, the problem of existence and uniqueness of the positive solution and the problem of exponential stability of disease free equilibrium are addressed. However, the problem of persistence of the disease is not observed.

Suppose that independent Brownian motions B 1 = {B 1 (t), t ≥ 0} and B 2 = {B 2 (t), t ≥ 0}, that govern the SDE system (5) , are defined on the complete probability space (Ω, F, F, P) with filtration

where F Bi t , i = 1, 2, are σ-algebras from natural filtrations of Brownian motions B 1 and B 2 and it contains all P-null sets. The problem of existence and uniqueness of a positive global solution is treated in the following theorem, where we use the notation R 7 + for the space {(x 1 , x 2 , x 3 , x 4 , x 5 , x 6 , x 7 ) : x i > 0, ∀i = 1, . . . , 7}. Theorem 1 For any initial value (S(0), E(0), I(0), P (0), H(0), A(0), R(0)) ∈ R 7 + there exists a unique positive solution (S(t), E(t), I(t), P (t), H(t), A(t), R(t)) of the SDE system (5) for every t ≥ 0, which almost surely remains positive, i.e.

Proof 1 The proof of this theorem follows the idea from [41] and [42] . First we prove that N (t) → Λ/µ as t → ∞. By summing all seven equations from the system (5) we obtain the following

More precisely, we obtain the equation

with solution of the following form 

Now we proceed to proof of the existence and uniqueness of the positive global solution. Namely, since equations 1 − 5 in system (5) are independent of equations 6 and 7, it is enough to prove the existence of a positive global solution of the system consisting of the first five equations:

Since the coefficients of system (9) are locally Lipschitz continuous, according to [45] , for any initial value (S(0), E(0), I(0), P (0), H(0)) ∈ R 5 + = {(x 1 , x 2 , x 3 , x 4 , x 5 ) : x i > 0, ∀i = 1, . . . , 5}, this system has a unique local solution on [0, τ 0 , where τ 0 denotes the explosion time. It is necessary to prove that the solution is global, i.e. that τ 0 = +∞ P-a.s.

Let k 0 > 0 be a constant large enough so that all initial values S(0), E(0), I(0), P (0), H(0) belong to the interval [1/k 0 , k 0 ]. Then, for each k > k 0 we define the stopping time

with the usual agreement that inf ∅ = ∞. Note that τ k increases as k is increasing and denote

Now it follows that τ ∞ ≤ τ 0 P-a.s. To complete the proof, we need to prove that τ ∞ = ∞ P-a.s. Since inf ∅ = ∞ and τ k ≤ τ 0 , the proof of τ ∞ = ∞ P-a.s. completes the proof of this theorem. Namely, if τ ∞ = ∞ P-a.s., then τ 0 = ∞ P-a.s., which means that (S(t), E(t), I(t), P (t), H(t)) P-a.s. remains in R 5 + for all t > 0.

Let us suppose the opposite, i.e. let us suppose that there exists a pair of constants T ≥ 0 and ε ∈ (0, 1) such that

Hence, there exists k 1 ≥ k 0 such that Note that since log (x) ≤ x − 1 for every x ≥ 0, V is a non-negative function. By applying the multidimensional Itô formula (see e.g. [45] ) to the function V , we obtain the following dV (S, E, I, P, H)

Let us introduce the following notation in the first two equations from the system (9)

From (11) it follows that dV (S, E, I, P, H) = L(S, E, I, P, H)dt+

where function L : R 5 + → R + is given by the following expression:

Next, we need an upper bound for the function L, for each t > 0. Since all random variables and all coefficients in (13) are non-negative, by applying elementary inequalities we obtain the following upper bound of the function L L(S, E, I, P, H) ≤ (S(t) + 1) (λ + βK 1 (t) + β K 2 (t) + µ) + σ 2 1 K 2 1 (t) + σ 2 2 K 2 (t) 2 +

Furthermore, since N (t) ≥ S(t) + E(t) + I(t) + P (t) + H(t), it follows that for every t > 0,

which, together with the upper bound (14) , yields the following estimate for the function L L(S, E, I, P, H) ≤ (N + 1) (λ + (1 + l)β + β + κ + δ i + δ p + 2γ a + 2γ i + 2µ) + +2N ((1 + l)β + β + κ(ρ 1 + ρ 2 )) + N 2 + 1 (1 + l) 2 σ 2 1 + σ 2 2 . := N .

Now we proceed with the analysis of the expectation of V for random time (τ k ∧ T )

Then under the introduced assumptions, it follows that P(A k ) ≥ ε. Even more, for every ω ∈ A k at least one of random variables S, E, I, P, H is either less or equal than 1/k or it is greater or equal than k. It follows that

Now it follows that

where IA k is the indicator function of the event A k . If we now let k → ∞, it follows that

which is a contradiction. This means that the assumption is not correct, i.e. it follows that τ ∞ = ∞ P-a.s. 

From the objective point of view, the virus remains persistent in population if there is at least one symptomatic infectious, asymptomatic infectious, hospitalized individual or super-spreader. From the mathematical point of view we observe the persistence in mean and, in accordance to previously stated objective explanation, we say that the system (5) 

We introduce the notation

Theorem 2 Let initial value (S(0), E(0), I(0), P (0), A(0), H(0), R(0)) ∈ R 7 + , such that the solution of the system (5) is in Γ , where µ, β, β and l satisfy the relation

and where c is a small fixed constant such that inf t≥0 E(t)/N (t) ≥ c.

If we assume that noises satisfy the condition

than the solution (S(t), E(t), I(t), P (t), A(t), H(t), R(t)) has the property

Proof 2 In order to obtain [I(s) + P (s) + H(s) + A(s)] we derive a specific weighted sum of equation 1 with equations 3 − 6 from the system (5)

where (1 − ρ 1 − ρ 2 ) > 0 and weight coefficients are

, K p = γ r + γ a + δ p (γ a + k 2 γ i + δ p )(γ r + δ p ) .

Let us introduce the following notation

By applying the Itô formula on log (V (t)), we obtain the following

Since 0 < K i , K p ≤ 1, and 0 < 1

under assumptions of the theorem and using that µ/Λ < 1/N (t) ≤ 1/V (t) ≤ 1, β < β and l < 1, we obtain the following lower bound for d (log V (t))

By integrating the last inequality from 0 to t and by dividing the whole expression by t, we obtain the following result

where

are continuous local martingales with values 0 at time t = 0 and for which

and therefore M 1 (t) t → 0 and M 2 (t) t → 0 P − a.s. as t → ∞.

Therefore,

which is greater than zero for

which completes the proof.

In this section the conditions for complete eradication of the disease (virus) from the population are proven. Theorem 3 If noises satisfy that 1 2

than for any initial value (S(0), E(0), I(0), P (0), A(0), H(0), R(0)) ∈ R 7 + , such that the solution of the system (5) is in Γ , it follows that

Proof 3 Theorem 1 guaranties existence of unique positive global solution of system (5) . To derive the conditions for extinction of the disease we apply the Itô formula to log E(t) and obtain the following upper bound

where we have used the following upper bound

Integration of expression (25) from 0 to t and division by t yield the following result 

From the lower bound (27) it follows that

and therefore lim t→∞ E(t) = 0 P − a.s.

To verify that I(t), P (t), H(t), A(t), R(t) → 0 P-a.s. as t → ∞ we solve equations 3 − 7 from the system (5) explicitly:

,

,

,

.

Since E(t) → 0 P-a.s. as t → ∞, from previous solutions it follows that I(t) → 0, P (t) → 0 and A(t) → 0 P-a.s. as t → ∞. Furthermore, it follows that H(t) → 0 and R(t) → 0 P-a.s. as t → ∞.

Recall that

and since E(t)

which completes the proof.

In this section we provide simulation results which illustrate the obtained theoretical results provided in Theorem 2 and Theorem 3. All simulations were done assuming daily increments (dt = 1/365) with 182500 simulated points for each stochastic process (simulation of approximately 500 years).

The proven results are illustrated for the parameter and starting values given in Table 1 . Regarding persistence proven in Theorem 2, assumption on the noises is fulfilled if we take σ 2 1 + σ 2 2 = 0.0061, since cκ ρ 1 γ r + γ a + δ p (γ a + k 1 γ i + δ i )(γ r + δ p ) + ρ 2 γ r + γ a + δ p (γ a + k 2 γ i + δ p )(γ r + δ p )

On the other hand, extinction of the disease can be seen in Figure 3 , since all processes converge to 0, while only susceptibles tend to Λ/µ = 1.5 as t → ∞. 

Mathematical models of the spread of the infectious disease are the crucial tool for tracking the course of the epidemic as well as for forecasting its persistence or extinction in the population. The most widely used models in epidemiology are compartmental models. Within model presented in this paper the population is divided into several disjoint classes that can mutually communicate. The communication between classes is described by model parameters. For example, the transmission coefficient, which incorporates the force of infection, is one of them. The nature of this parameter is highly influenced by the environmental uncertainty and therefore it is much more suitable to assume the stochastic nature of this parameter than to take its value as a constant or define its change via some deterministic function depending on time.

In this paper we followed this exact idea. We introduced finer parametrization and perturbed the constant transmission rates regarding the symptomatic infected individuals (β), hospitalized individuals (lβ) and superspreaders (β ) in the system of ordinary differential equations first published in [21] . The perturbation of transmission coefficients was performed by adding of mutually independent Brownian noises dB 1 and dB 2 to β and β , respectively. The resulting system (5) is the two-diffusion stochastic model, consisting od two Itô SDEs (regarding classes S and E) and five ODEs (regarding classes I, P, H, A, R). By using the classical stochastic techniques from Itô's calculus explained in the Introduction as well as in the introductory part of subsection Problems and Methodology, we proved the existence and uniqueness of the solution of the system (5) In comparison to other models in recent publications, by the best of our knowledge, our model (5) is the only stochastic model that includes the analysis of four disjoint classes of infected individuals -symptomatic infected individuals, asymptomatic carriers, superspreaders and hospitalized individuals. Furthermore, this is the compartmental model for which, under the existence of the unique positive global solution, both the persistence in mean and the extinction of the disease are analyzed both theoretically and practically within the simulation study.

Stochastic model (5) can be generalized to the model in which all the coefficients can be functions of several variables, e.g. delayed functions of time t catching the incubation period and the duration of the infectiousness, of location m covering the different regional population densities and different intensities of movements within the regions, of age a etc. Furthermore, parameters of the model could be defined to depend of any observable value which is relevant for the spread of the disease.

Beside this "stratification" of model coefficients, our future plan is to introduce heavy-tailed diffusions as models for transmission rates and maybe for some other parameters, which yields much more complicated stochastic model depending od Itô integrals with unknown distribution which then has to be either approximated or estimated. Furthermore, these diffusions could be non-autonomous or driven by correlated Brownian motions, they could be driven by the noise of different nature (e.g. Lévy, fractional and time-changed Brownian motion) producing respectively jump processes and long-memory processes with trajectories allowing jumps and resting periods, which could be important for catching the superspreader events.

Tracing the origins of SARS-CoV-2 in coronavirus phylogenies: A review

Stability of SARS-CoV-2 in different environmental conditions

The transmission modes and sources of COVID-19: A systematic review

Virology, transmission, and pathogenesis of SARS-CoV-2

Asymptomatic transmission of COVID-19

Quantifying asymptomatic infection and transmission of COVID-19 in New York City using observed cases, serology, and testing capacity

Do superspreaders generate new superspreaders? A hypothesis to explain the propagation pattern of COVID-19

A comprehensive review of COVID-19 characteristics

A systematic review of COVID-19 epidemiology based on current evidence

Epidemiology of COVID-19: A systematic review and meta-analysis of clinical characteristics, risk factors, and outcomes

Modeling latent infection transmissions through biosocial stochastic dynamics

COVID-ABS: An agent-based model of COVID-19 epidemic to simulate health and economic effects of social distancing interventions

An agent-based modeling of COVID-19: Validation, analysis, and recommendations

Measuring and preventing COVID-19 using the SIR model and machine learning in smart health care

Regression analysis of COVID-19 using machine learning algorithms

Short-term forecasting COVID-19 cumulative confirmed cases: Perspectives for Brazil

A comparison of deep machine learning algorithms in COVID-19 disease diagnosis

Mathematical modeling of COVID-19 transmission dynamics with a case study of Wuhan

A new compartmental epidemiological model for COVID-19 with a case study of Portugal

Simulating the spread of COVID-19 via a spatially-resolved susceptible-exposed-infected-recovered-deceased (SEIRD) model with heterogeneous diffusion

A mathematical model of coronavirus disease (COVID-19) containing asymptomatic and symptomatic classes

Modeling and forecasting the spread of COVID-19 with stochastic and deterministic approaches: Africa and Europe

A discrete stochastic model of the COVID-19 outbreak: Forecast and control

Modeling the second wave of COVID-19 infections in France and Italy via a stochastic SEIR model

Modelling COVID-19 contagion: Risk assessment and targeted mitigation policies

Sequential data assimilation of the stochastic SEIR epidemic model for regional COVID-19 dynamics

Environmental Brownian noise suppresses explosions in population dynamics

An Introduction to Stochastic Epidemc Models

An Introduction to Stochastic Processes with Applications to Biology

Probabilistic Models of Population Evolution

A stochastic differential equation SIS epidemic model

Estimating the state of the COVID-19 epidemic in France using a model with memory

A stochastic differential equation SIS epidemic model incorporating Ornstein-Uhlenbeck process

On the integral of geometric Brownian motion

On the fundamental solution of the Kolmogorov-Shiryaev equation

Stochastic SIR model with jumps

Dynamics of a stochastic SIS model with double epidemic diseases driven by Lévy jumps

A stochastic SICA epidemic model for HIV transmission

A stochastic analysis of the impact of fluctuations in the environment on pre-exposure prophylaxis for HIV infection

Mathematical perspective of COVID-19 pandemic: Disease extinction criteria in deterministic and stochastic models

Mathematical analysis of a stochastic model for spread of coronavirus

Stochastic Differential Equations

Lyapunov functions and global stability for SIR, SIRS, and SIS epidemiological models

Nonlinear Systems Stability Analysis, Lyapunov-Based Approach

Epidemic models with varying infectivity

Functional central limit theorems for epidemic models with varying infectivity

Modeling and forecasting of COVID-19 spreading by delayed stochastic differential equations

A stochastic time-delayed model for the effectiveness of Moroccan COVID-19 deconfinement strategy

A model describing COVID-19 community transmission taking into account asymptomatic carriers and risk mitigation

Sample CRediT author statement • Jasmina Ðorđević

Furthermore, we thank the Editor, the Associate Editor and two anonymous Referees for their constructive suggestions that led to substantial improvement of the paper.

This is illustrated in Figure 2 . Moreover, Figures 2 (2h, 2j) show that necessary conditions are valid as well, i.e. proportion of the exposed individuals in population does not drop below c = 0.01. Persistence of stochastic process modeling infection (I, P, H, A) can be seen in Figures 2c -2f , while the validity sufficient condition (20) of the Theorem 2 is illustrated in Figure 2j . Formally, persistence in mean can be verified in Figure 2k where [I(t) + P (t) + H(t) + A(t)] never drops below the desiblack constant c κρ1 γr + γa + δp (γa + k1γi + δi)(γr + δp) + κρ2 γr + γa + δp (γa + k2γi + δp)(γr + δp) + κ(1 − ρ1 − ρ2) γi + µ − σ 2 1 + σ 2 2 c = 0.00243793.