Threshold dynamics of a stochastic SIHR epidemic model of COVID-19 with general population-size dependent contact rate

: In this paper, we propose a stochastic SIHR epidemic model of COVID-19. A basic reproduction number R s 0 is deﬁned to determine the extinction or persistence of the disease. If R s 0 < 1, the disease will be extinct. If R s 0 > 1, the disease will be strongly stochastically permanent. Based on realistic parameters of COVID-19, we numerically analyze the e ﬀ ect of key parameters such as transmission rate, conﬁrmation rate and noise intensity on the dynamics of disease transmission and obtain sensitivity indices of some parameters on R s 0 by sensitivity analysis. It is found that: 1) The threshold level of deterministic model is overestimated in case of neglecting the e ﬀ ect of environmental noise; 2) The decrease of transmission rate and the increase of conﬁrmed rate are beneﬁcial to control the spread of COVID-19. Moreover, our sensitivity analysis indicates that the parameters β , σ and δ have signiﬁcantly e ﬀ ects on R s 0 .


Introduction
The COVID-19 has been significantly impacting our lives since the emerge of the first case in early December 2019 in Wuhan, China. As of 14 January 2022, there have been 318,648,834 confirmed cases, including 5,518,343 deaths in the world [1]. Multiple mutant strains of COVID-19 have emerged, such as Alpha, Beta, Gamma, Delta, Omicron and so on. Among these strains, Omicron, the newly discovered strain, is known as the fast speed of transmission and the strong ability of infection. Therefore it is extremely urgent to study and control the transmission of COVID-19.
To better understand the transmission and develop efficient control strategies, researchers have employed mathematical models to analyze the dynamic behavior and control the outbreak of COVID-19 [2][3][4][5][6][7][8][9][10]. Recently, Allegretti et al. [3] considered a modified SIR model of COVID-19 and found that a high fraction of avoided contacts leads to the stability of the disease free equilibrium. Naik et al. [4] proposed a COVID-19 epidemic model and indicated that reducing transmission rate of the coronavirus is the most essential strategy to prevent the virus further spread. Okuonghae et al. [5] formulated a mathematical model and examined the impact of various non-pharmaceutical control measures on the population dynamics of COVID-19 by using the available data from Lagos, Nigeria. Fatma et al. [6] studied the interactions between COVID-19 and diabetes by using real data from Turkey and numerically visualized the population dynamics of COVID-19. Tang et al. [7] proposed a Filippov SIR model to investigate the impacts of three control strategies (media coverage, vaccination and treatment) and choose the switching policy properly to reduce the infected size. Humphrey et al. [8] developed an SEIRL model and found that testing and tracing asymptomatic individuals frequently can help in controlling new cases. Jin et al. [9] proposed a generalised SEIR model to seek optimal strategies for disease control, finding that reducing the transmission rates and increasing contact tracing are possible to hinder the fast spread of COVID-19. It is found that in the early stages of COVID-19, isolating confirmed cases was considered as a more effective control measure due to the inability to quickly produce highly effective vaccines [11]. Hence, Jiao et al. [11] proposed an SIHR model incorporating confirmed cases with general population-size dependent contact rate as follows: Here the infectious cases are divided into two sub-populations: non-confirmed cases (I) and confirmed cases (H). The non-confirmed cases are infected individuals who have not been tested by medical institutions. Once the nucleic acid tests are positive, they would become confirmed cases and be isolated. In model (1), I(t) and H(t) are the non-confirmed infected individuals and the confirmed individuals at time t, respectively; and S (t), R(t) and N(t) denote respectively the susceptible individuals, the recovered individuals and the total population. A is the recruitment rate of the population and it is assumed that all the newcomers are susceptible. β f (N)S I is the general population-size dependent incidence, in which the parameter β is the transmission rate from the infectious class to the susceptible class, and f (N) is a function of N and it comes in many forms, see [11][12][13]. µ and γ denote the natural death rate and the natural recovery rate, respectively. δ is the confirmation rate from the infected population to confirmed cases. µ 1 , µ 2 are respectively the extra disease-related death rate constants in compartments I and H. m is the transform rate from the confirmed population to the removed population. All parameters are nonnegative constants. Model (1) always exists the disease-free equilibrium E 0 = ( A µ , 0, 0, 0), whose stability is determined by the basic reproduction number µ(δ+γ+µ+µ 1 ) . If R 0 < 1, then E 0 is globally asymptotically stable; and if R 0 > 1, E 0 becomes unstable and an endemic equilibrium E * (S * , I * , H * , R * ) appears and it is locally asymptotically stable, see [11] for details. However, deterministic model (1) has certain limitation, it cannot describe the effects of random environment. In fact, there are many stochastic factors that can effect the transmission of disease. For example, Jamshidi et al. [14] investigated the impact of mobility, urban density, population, homestay, and mask-wearing separately on COVID-19 by conducting a multiple regression analysis and found that a higher level of population mobility and traveling can increase the transmission rate. Sabbir Hossain et al. [15] studied the impact of weather on COVID-19 in part of South Asian countries through adopting the Autoregressive Integrated Moving Average with Explanatory Variables (ARIMAX) model and revealed that maximum wind speed had a significant negative effect on the transmission rate in India, whlie rainfall, relative humidity, temperature and maximum air pollutants matter PM 2.5 had different influence on COVID-19 in different areas. We also refer the readers to Habeebullah [16], Baniasad [17] and Damette [18] for learning more about the effect of weather, pollution and mobility to the transmission efficiency of COVID-19. Therefore, the effects of these random factors can be translated to the fluctuations in the transmission rate β [19]. It has been well established in literatures that introducing parameter perturbations can affect the dynamic behavior of population. Gray et al. [20] studied the effect of stochastic parameter perturbation on SIS model and fill previous gap. Li et al. [21] applied similar parameter perturbation to investigate the effect of stochastic environmental variability on inter-pandemic transmission dynamics of influenza A. Cai et al. [22] showed that appropriate parameter perturbation to the system is useful in controlling the spread of the disease. Motivated by these, in this paper we suppose that the transmission rate β fluctuates around an average value due to the continuous fluctuations in the environment by the white noise β + σḂ(t), and then obtain the following stochastic model: where B(t) is a standard Brownian motion, which is defined on the complete probability space (Ω, F , {F t } t≥0 , P), and σ represents the noise intensity; the function f (·) : R + → R + satisfies the following assumptions: where f (x) and ( f (x)x) respectively denote the derivative of f (x) and ( f (x)x). The main purpose of this paper is to explore the effect of random variability in the environments on the spread of COVID-19 based on realistic parameters from [11] and [23], and obtain the strict threshold condition of the disease. The main contributions of our study can be summarized as follows.
• It has been shown in [24] that uncertainty is certain in the disease transmission rate of COVID- 19 and there are large variation in its range. Therefore, it seems necessary and important to consider random factors in the context of COVID-19.
• Under the setting of more general population-size dependent contact rate, we obtain the threshold condition of disease extinction and persistence by constructing suitable Lyapunov functions. In this sense, we extend the previous studies such as [25] and [26], where the standard incidence form is used.
• We have proved that the basic reproduction R s 0 for our stochastic model completely determines the extinction or persistence of the disease. This is contrasted with the existing literatures [27][28][29][30][31][32], where except for the conditions required for the basic reproduction number, there still have some additional conditions for noise to ensure the extinction or persistence of disease. Therefore our results can be regard as an significant extension of the previous articles and also can provide effective information to the control of COVID-19.
The organization of this paper is as follows. In Section 2, we prove the existence and uniqueness of the positive solution. In Section 3, using techniques in [33][34][35][36] we derive the threshold R s 0 which completely determines the extinction and strongly permanent of the disease. Finally, in Section 4, numerical simulations are carried out to illustrate our theoretical results by analyzing the effect of key parameters on disease and obtain sensitivity indices of some parameters on R s 0 by sensitivity analysis. We also discuss the impact of some measures (media coverage, government intervention, testing and tracking) on COVID-19 and give a summary.

Existence and uniqueness of the positive solution
In order to investigate the dynamics of stochastic model (2), we first need to show that the model has a unique global positive solution. Denote by Φ(t) = (S (t), I(t), H(t), R(t)) the solution of model (2) and Moreover, for any function V ∈ C 2 (R 4 , ∞); R + ), define the differential operator L associated with model (2) as where V Φ (Φ) and V ΦΦ (Φ) are the gradient and Hessian of V(·);f andg are respectively the drift and diffusion coefficients of model (2). By Itô's formula, We have the following theorem.
Theorem 2.1. For any initial value S (0), I(0), H(0), R(0) ∈ R 4 + , there is a unique solution S (t), I(t), H(t), R(t) of model (2) on t ≥ 0, and the solution will remain in R 4 + with probability 1. Proof. Notice that the coefficients of model (2) satisfy the local Lipschitz condition. Then there is a unique local solution S (t), I(t), H(t), R(t) on [0, τ e ) for any initial value S (0), I(0), H(0), R(0) in R 4 + , where τ e is the explosion time [37]. In order to show the global existence of the positive solution, we need to prove that τ e = ∞ almost surely (a.s.).
It is easy to see from model (4) that ∆ is the positive invariant set of stochastic model (2). Thus, in the sequel of this paper, we only need to consider the dynamics of model (2) constrained in ∆.

Threshold dynamics of the disease
In this section, we perform the persistence and extinction analysis of stochastic model (2). Define We can see from below that R s 0 plays the similar role of basic reproduction number of disease as defined in classical deterministic epidemic models, called stochastic basic reproduction number, which completely determines the dynamics of stochastic model (2). Denote and let λ : It is easy to check that R s 0 = 1 implies λ = 0, and moreover R s 0 < 1 if and only if λ < 0.

Extinction
Consider the following stochastic differential equation of x(t) ∈ R n : where a(·) : R n → R n and b(·) : Assume that x(t) = 0 is the trivial solution of model (11). The following lemma presents a proper adaptation of Theorem 3.1 in Dang and Yin [38], which will be used later to establish the condition for the extinction of disease. The following result is about the extinction of disease.
. Notice the continuity of functions g(·) and f (·). We can take δ 1 and p sufficiently small such that model (16) holds and for any (S , I, H, R) ∈ U δ 1 , the following two inequalities hold: where 0 < ξ 1 , ξ 2 < ξ and ξ 1 + ξ 2 = ξ. Consequently, we have Moreover, we can easily check that Combining models (15), (17) and (18), we know that if we take p and δ 1 both sufficiently small, it then follows from model (14)  According to Lemma 3.1, we know that the disease free equilibrium ( A µ , 0, 0, 0) is asymptotically stable in probability. That is, for any ε > 0, there exists a δ 2 , 0 < δ 2 < δ 1 such that for any (S (0), 3 . Now we are in a position to prove that any solution starting in ∆ will eventually enter U δ 2 .

Persistence
In this section, we prove that the disease will be persistent provided R s 0 > 1. We first present the following useful lemma.
The following theorem is about the persistence of the disease.
4.1. The effect of β, σ and δ on the dynamics of model (2) In the subsection, we numerically simulate the solution of model (2) using MATLAB R2016b to illustrate the theoretical results obtained in Section 3, mainly revealing the effect of σ, β and δ on the dynamics of model (2). The numerical scheme is obtained through Milstein's higher order method [41].  Table 1. In this case one can get that R 0 = 2.3081. It then follows from [11] that deterministic model (1) has a unique endemic equilibrium E * (12824304.8094, 3743.8322, 9661.8313, 16762296.1178), which is stable. It follows from R s 0 = 1.9900 that we have σ = 0.1000. By Theorem 3.4, we can get that the disease is strongly stochastically permanent. The computer simulations shown in Figure 1(a) clearly support the result. To show how the noise affects the dynamics of disease, now we take σ = 0.0100, and other parameters remain unchanged. In this case, we obtain R s 0 = 2.3045 > 1. The similar simulation result is shown in Figure 1(b). We observed that the path of I(t) for model (2) is oscillating around the steady state value I * = 3743.8322. Compared with Figure 1(a), one can get that when R s 0 > 1, the small noise does not change the stability of the equilibrium state of model (2), but with the intensity of white noise increasing, the volatility of I(t) is getting larger. Finally, we consider σ = 0.3000. It is easy to compute that R s 0 = 0.9494 < 1, according to Theorem 3.2, the disease will go extinct almost surely as shown in Figure 1(c). However, deterministic model (1) claims the persistence of the disease. This discrepancy highlights the impact of stochastic environmental to the disease dynamics.  Table 1 (Color figure online).

Case 2. The effect of transmission rate
Here we assume b = 0.3900, σ = 0.1000, and other parameters take values as in Table 1. In this case we consider three different values of β = 0.0500, 0.1860, 0.3500 to see the effect of transmission rate on the spread of infectious disease. The corresponding values of R s 0 are respectively 0.5353, 1.9900 and 3.7473. By Theorems 3.2 and 3.4, we can get that the disease is strongly stochastically permanent when β = 0.1860 and β = 0.3500, while the disease is extinctive when β = 0.0500. The computer simulations are shown in Figure 2(a) which clearly support these results. Figure 2(b) shows the corresponding persistence level of I(t) for various values of β. It is observed that when R s 0 > 1, the persistence level of I(t) is reduced gradually with the decrease of transmission rate. This indicates that decreasing transmission rate is beneficial to the control of the spread of COVID-19. So we can take some measures to reduce the scale of outbreaks by decreasing transmission rate. For example, transmission rate can be reduced through improving the media response rate to reports on the severity of COVID-19 and encouraging citizens to actively prevent disease. Moreover, government can adopt a series of policies including wearing masks, avoiding farm and wild animals, travel restrictions, stay at home, lockdowns, and so on to decline transmission rate. These measures could effectively reduce the number of infected cases and suppress the outbreak of disease.
Case 3. The effect of confirmed rate Here we assume β = 0.1860, b = 0.3900, σ = 0.1000, and other parameters taking values as in Table 1. In this case we choose three different values of δ − 0.1000, 0.1836, 0.6000 to see the effect of confirmed rate on the spread of infectious diseases. The corresponding values of R s 0 are respectively 3.0601, 1.9900, 0.7270. According to Theorems 3.2 and 3.4, the disease persists when δ = 0.1000 and δ = 0.1836, while the disease will be extinct when δ = 0.6000. Figure 3(a) clearly support these results. Figure 3(b) shows that the corresponding persistence level of I(t) for various values of δ. We can see that when R s 0 > 1, the persistence level of I(t) is reduced gradually with the increase of confirmed rate. This indicates that the increase of confirmed rate is beneficial to control the spread of COVID-19, while blind testing is not desirable, it will cause a huge burden on society. Therefore, in order to enhance confirmation rate, a positive tracking and testing strategy should be carried out to control the spread of disease [8]. Frequently testing in smaller scale populations, such as schools, factories, community, etc., where virus is more easier to spread, and test less frequently parts of the population who are not as exposed. Detecting continually close contacts also leads to the increase of confirmed cases.

Parameters Description Sensitivity Index A
The recruitment rate +0.0001477 µ The natural death rate -0.0003257 γ The natural recovery rate -0.0773 β The transmission rate +1.0000 δ The confirmed rate -0.7675 σ The noise intensity -0.2746 The disease-induced death rate of infected individuals -0.0184

Sensitivity analysis
Varying parameter values will have different effects on the output of model (2). In order to qualitatively analyze the influence of some parameters on the output of model (2), the sensitivity analysis method is adopted. The normalized forward sensitivity index of a variable, R, that depends on a parameter, h, is defined as: Using the above formula, we analyze the sensitivity of state variable R s 0 to the following parameters of model (2) Figure 4. The sensitivity indices of state variable R s 0 with respect to some parameters for model (2).
From Table 2, parameters with positive sensitivity index, A and β, indicate that the transmission of COVID-19 increases with the increase of these two parameters. Similarly, parameters with negative sensitivity index, µ, γ, δ, σ and µ 1 , mean that the transmission of COVID-19 decreases with the increase of these parameters. As shown in Figure 4, we observe that β, σ and δ have significant effects on R s 0 . This verifies that our analysis of parameters is meaningful in Subsection 4.1.

Conclusions
In reality, there exist many random environmental factors, like weather, relative humidity, temperature and population mobility, which may have significant effects on the transmission of COVID-19. Therefore, considering stochastic influences into the epidemic model seems necessary and important. In this paper, we propose and investigate a stochastic SIHR epidemic model with the environmental variability in the transmission rate to describe the transmission of COVID-19, and based on which we numerically illustrate the evolution dynamic of COVID-19 using the realistic parameter values from literatures.
The main contribution of our paper can be summarized as the following two aspects. Mathematically, we prove that the stochastic dynamics of stochastic model (2) is completely determined by the reproduction number R s 0 : If R s 0 < 1, the disease will go to extinction ultimately, and if R s 0 > 1, the disease stochastically permanent. Epidemiologically, we can conclude that: (i) The presence of environmental noises can sustain the irregular recurrence of disease and the volatility of infected population increases with the increasing noise intensity if R s 0 > 1. When the noise increases to a certain level such that R s 0 < 1, the disease will go to extinction (See Figure 1). And also it can be seen from Figure 1(c) that white noise may reshape the solution behavior of corresponding deterministic model (1). In other words, noise may change the evolution tendency of disease. (ii) The decrease of transmission rate and the increase of confirmed rate are beneficial to the control of COVID-19 spread (See Figures 2 and 3). (iii) Our sensitivity analysis indicates that the transmission rate β, noise intensity σ and confirmed rate δ are the most sensitive parameters to R s 0 (See Figure 4). More than two years have lasted since the emergence of COVID-19 in the world. It is well known the transmission of disease will necessarily be affected by other factors such as media coverage, seasonal changes and so on [42][43][44][45][46][47]. Considering the seasonal effect or the switching of environments in model (2) will be an interesting research topic. We leave this for our future investigation.