Dynamical properties of two-diffusion SIR epidemic model with Markovian switching

1 Introduction and motivation

The mathematical approach in epidemiology is related to the study of the spread of infectious diseases in a population. These diseases can be transmitted from animals to humans or from humans to humans through physical contact with an infected person or animal, as well as through the contact with person or animal that carry the virus, so that they can transmit it to someone else, even though they are not infected themselves. The infection can also be transmitted through air, water, soil, etc. Given that experiments in the field of epidemiology are not possible due to numerous ethical and technical reasons, it turned out that mathematical models are a very useful tool that can be used to predict the outbreak of epidemics, but also to measure the effects of various preventive measures on their control and suppression. Different types of the diseases are modelled by different epidemiological models. In mathematical literature, there are many studies that consider deterministic models that describe the spread of a disease in a large population (see [1–3], for instance). These models are mainly based on the so-called mass action law, so it is assumed that people within a population mix homogeneously, which in fact means that they come into contact independently of each other and everyone has an equal chance to encounter any other individual. The classic susceptible-infectious-recovered (SIR) epidemiological model is one of the most commonly implemented models where the total population at time t denoted by N ( t ) is divided into three compartments: “susceptible” S ( t ) , “infected” I ( t ) , and “recovered” R ( t ) . It is used to model diseases after which the recovered patients acquire permanent immunity. It was initially proposed and studied by Kermack and McKendrick in a set of three articles from 1927, 1932, and 1933. Because of their importance to the field of theoretical epidemiology, these articles were republished [4–6]. The model is of the form

with initial condition S ( 0 ) = S 0 , I ( 0 ) = I 0 , R ( 0 ) = R 0 , while β is the transmission rate and γ natural recovery rate.

In [1], the authors constructed the SIR epidemiological model using a nonlinear saturated incidence rate given by β I S 1 + k I , unlike the incidence rate in (1), β I S . This is great improvement because the incidence rate of such form does not allow unlimited number of contacts, and the term 1 1 + k I represents precautionary measures for individuals with a compromised immune system, making them more susceptible to infection. In addition to precautionary measures, treatment is another method to prevent and control the spread of infectious diseases. Zhou and Fan [1] considered a nonlinear saturated treatment function expressed as α I ω + I and investigated how limited medical resources and their supply efficiency affect the transmission of infectious diseases. The model is represented by a system of ordinary differential equations of the form

with the initial conditions S ( 0 ) = S 0 , I ( 0 ) = I 0 , R ( 0 ) = R 0 , while biological interpretation of the parameters is described in Table 1.

Since people within the population do not meet completely by chance, but each person has a certain circle of people with whom he/she is in daily contact, the assumption of homogeneous mixing in deterministic models is not realistic. That is why, stochastic models represent the more appropriate way of modelling epidemics than their deterministic analogues. In the literature, there are many studies that consider dynamics of the stochastic epidemiological models [7–13]. Due to the fact that almost all diseases have an incubation period, many authors naturally assumed that epidemiological models with time delay are more suitable to analyse the spread of the infectious diseases (see [14–16], among others).

Based on model (2), Milunović and Krstić [12] the authors constructed stochastic SIR epidemiological model of the form

with the initial conditions S ( 0 ) = S 0 , I ( 0 ) = I 0 , R ( 0 ) = R 0 , where σ 1 and σ 2 are the real constants, and B 1 ( t ) and B 2 ( t ) are two independent standard Brownian motions defined on complete probability space ( Ω , ℱ , P ) with the filtration { ℱ t } t ≥ 0 , satisfying the usual conditions (it is right continuous and increasing, while ℱ 0 contains all P-null sets), while other parameters are defined as in Table 1. For model (3), the authors showed the existence and uniqueness of the global positive solution and that the solution remains in the set

provided that ( S 0 , I 0 , R 0 ) ∈ Γ . Moreover, they proved that the solution admits the unique stationary distribution, which means that the system is stochastically weak stabile, i.e. the disease will prevail in the population in long time, and then, they gave the sufficient conditions under which the disease-free equilibrium of system (3) is p th moment and almost surely exponentially stable, as well as the conditions under which extinction of the disease occurs.

Population dynamics is often affected by various types of external environmental noises, such as the changes in weather, temperature, and climate. Since these fluctuations may cause dramatic changes in some important parameters of infectious diseases (e.g., transmission rate, death rate, recovery rate), they are generally described as switching between two or more regimes of the environment. In the literature, the switching can be generally illustrated as a continuous-time Markov chain, which drives the changes of the main parameters of population and epidemiological models with state switchings of homogeneous Markov chain [17–26]. With that motivation, we introduce coloured or telegraphic noise in (3), which will be represented by right-continuous Markov chain { ϑ ( t ) , t ≥ 0 } , defined on complete probability space ( Ω , ℱ , P ) with the filtration { ℱ t } t ≥ 0 , independent of Brownian motions B 1 ( t ) , B 2 ( t ) , and taking values in a finite space S = { 1 , 2 , … , N } . Hence, model (3) becomes

with the initial conditions S ( 0 ) = S 0 , I ( 0 ) = I 0 , R ( 0 ) = R 0 , ϑ ( 0 ) ∈ S , and for any ι ∈ S , coefficients Λ ( ι ) , β ( ι ) , k ( ι ) , μ ( ι ) , ε ( ι ) , γ ( ι ) , α ( ι ) , ω ( ι ) , σ 1 ( ι ) , and σ 2 ( ι ) are positive. The infinitesimal generator matrix Q = [ λ ι κ ] N × N is defined by

where λ ι κ ≥ 0 is the transition rate from state ι to κ if ι ≠ κ while λ ι ι = − ∑ κ ≠ ι λ ι κ . We assume that λ ι κ > 0 for ι , κ = 1 , 2 , … , N with ι ≠ κ . This assumption provides that the Markov chain ϑ ( t ) is irreducible, i.e. there exists a unique stationary distribution π = ( π ι ) 1 × N such that π Q = 0 , ∑ ι = 1 N π ι = 1 and π ι > 0 for any ι ∈ S .

For discussion convenience, define ρ ˆ = max ι ∈ S ρ ( ι ) and ρ ˇ = min ι ∈ S ρ ( ι ) .

As it was already mentioned, there are many articles in the mathematical literature dealing with the dynamical behaviour of SIR epidemiological models with or without Markovian switching. Since we present a detailed analysis of a stochastic SIR epidemiological model with treatment, non-linear incidence rate, two diffusions, and Markovian switching, we regard this study as a generalization of the articles that can be found in literature. More precisely, we can consider articles [18,25] as a special case of our model, because in the study [18], there is no treatment function and second diffusion, i.e. α = 0 and σ 2 = 0 . This is a shortcoming of [18,25] as it ignores the impact of the environment on limited medical resources, especially considering the coronavirus disease 2019 (COVID-19) pandemic, where the unpredictability of this parameter has been shown. Also, the recent catastrophic earthquake in Turkey and Syria, which destroyed many hospitals and killed many medical workers, proves the importance of the parameter σ 2 in our model. Furthermore, compared to the extinction results obtained in [18], we provide the weaker conditions for model parameters, which is improvement of the existing results.

For all epidemiological models, the basic reproduction number is a very important quantity, which describes the dynamics of the disease. It represents the average number of secondary infections that single infected person may generate in a susceptible population. For model (2), the basic reproduction number is computed in [1] by the next-generation method, and it is given by

In general, it is more difficult to determine the basic reproduction number for stochastic models. For model (4), we will show that the quantity ℛ 0 s defined by

where

and model parameters satisfy the condition ℛ 0 s > 0 , i.e.

determines the extinction and persistence in mean of the disease. It is easy to see that ℛ 0 s becomes ℛ 0 when there is no white and coloured noise.

The rest of this article is organized as follows. In the next section, we prove that model (4) has a unique global positive solution for any initial positive value. In Section 3, we give sufficient conditions under which the disease is persistent in mean. Results on extinction are presented in Section 4. In Section 5, we illustrate the theoretical results from Sections 3 and 4 by using real data for Ebola outbreak in Sierra Leone and COVID-19 pandemic in Pakistan. The obtained simulations show that the model (4) predicts the course of the disease well, and can be used for all diseases and disease variants in which permanent immunity is acquired. We close the study with a conclusion section.

2 Existence and uniqueness of global positive solution

From a biological point of view, it is important to see whether the solution of (4) is global and positive. Consider the set

Based on Theorem 3.1 in [12], we present the following theorem.

3 Persistence in mean

In this section, we consider persistence in mean for model (4). From the epidemiological point of view, it is not desirable property because it indicates that the disease prevails in population. The definition of persistence in mean may be found in [8], for example.

To facilitate the discussion, we define

where N is given by (7). It can be seen that ℛ 1 s is an approximation of the basic reproduction number ℛ 0 s defined by (6).