Abstract
In this study, a deterministic SEQIR model with standard incidence and the corresponding stochastic epidemic model are explored. In the deterministic model, the reproduction number is given, and the local asymptotic stability of the equilibria is proved. When the reproduction number is less than unity, the disease-free equilibrium is locally asymptotically stable, whereas the endemic equilibrium is locally asymptotically stable in the case of a reproduction number greater than unity. A stochastic expansion based on a deterministic model is studied to explore the uncertainty of the spread of infectious diseases. Using the Lyapunov function method, the existence and uniqueness of a global positive solution are considered. Then, the extinction conditions of the epidemic and its asymptotic property around the endemic equilibrium are obtained. To demonstrate the application of this model, a case study based on COVID-19 epidemic data from France, Italy, and the UK is presented, together with numerical simulations using given parameters.
1. Introduction
At the end of 2019, COVID-19 was reported in Wuhan, China [1], setting the beginning of a global epidemic. Now, there are more than 200 countries trapped in the epidemic disaster, which has seriously affected the development of the society. It was found that COVID-19 was caused by the SARS-CoV-2 coronavirus [2, 3]. This coronavirus mainly spreads in three ways: (i) direct transmission, which refers to the infection caused by patient’s sneezing, coughing, and talking droplets, and the direct inhalation of exhaled gas at a short distance; (ii) aerosol transmission, which refers to droplet mixing in the air to form an aerosol that can lead to infection after inhalation; (iii) contact transmission, which refers to the droplet deposited on the surface of an object, contact with contaminated hands, and then contact with the oral cavity, nasal cavity, eyes, and other mucous membranes, resulting in infection.
In the absence of a specific therapeutic method, governments can only do their best to hinder the spread of the epidemic through measures such as population lockdown, reasonable allocation of medical resources, and calls for public health (mask wearing, avoidance of places with many people, meetings with reduced attendance, etc.) [4, 5]. These measures can effectively prevent epidemics and delay the arrival of epidemic peaks [6]. For now, through the aforementioned measures and vaccine production, the epidemic situation has been controlled in many countries; however, there are still some countries with severe epidemic situations.
Mathematical modeling is a suitable approach to describe the spread of diseases and predict their epidemic trend [7]. Based on the spread characteristics of diseases, many epidemic models, e.g., SIR, SEIR, and SIQR [8–12], were proposed. Using mathematical models to describe infectious diseases can help to better understand the spread trends of diseases and propose reasonable control methods and suggestions. Novel coronavirus pneumonia models were reported by scholars who have been involved in this task since the outbreak of the epidemic. Tang et al. constructed an SEIR model in [13]; their results indicate that the number of infected people and the propagation of infectious diseases can be reduced when governments take certain measures. In [14], the author predicted that the spread of infectious disease would be costly. After that, the governments of the UK and USA established strict movement restrictions. Many models, related suggestions, and opinions on COVID-19 have been provided [15–21].
Although there are many mathematical models, few considered the inherent uncertainties in these models. Manski and Molinari [22] showed that there is an uncertainty in the epidemic transmission process. Therefore, assuming that stochastic disturbance applied to a model is suitable, some stochastic epidemic models were proposed [18, 23–26]. In [18], the authors studied an SLIRD with the deterministic and stochastic models to analyze actual COVID-19-infected cases in Spain. In the numerical simulations, they fitted data from Spain and obtained the daily basic reproduction number. Thus, they predicted the epidemic situation for the next few months. In [23], by extending the SIR epidemiological model, the authors established a random propagation model for additional modeling of whether an individual is far from a crowded area. The results showed that the propagation of the epidemic in Japan would gradually slow down by reducing the time spent in crowded zones. Anwarud et al. analyzed a stochastic SIQ model and obtained sufficient conditions for disease extinction, disease existence, and stationary distribution in [24]. The numerical simulations were divided into two parts: fitting data from Khyber Pakhtunkhwa, Pakistan, and verification of the theoretical results for the given parameters. In [25], the authors studied an SEQIR model with Markovian switching, established a random threshold of disease extinction and persistence, and used the data from Indian states to confirm their conclusions. In [26], an epidemic system was considered through an additive fractional white noise. It was indicated that the epidemic under fractional random settings may be more suitable for modeling than deterministic modeling. Based on the above studies, an SEQIR model with standard incidence is analyzed in this study, together with deterministic and stochastic models. A case study of COVID-19 epidemic data from France, Italy, and the UK is presented, and numerical simulations using given parameters are reported.
The remainder of this paper is organized as follows. In the next section, deterministic and stochastic SEQIR models with standard incidence are described. The dynamic analysis of the deterministic model is presented in Section 3, and an analysis of the stochastic system is presented in Section 4. Numerical simulations and application to the COVID-19 epidemic are described in Section 5. The paper ends with a discussion in Section 6.
2. Model Formulation
The entire population is divided into five classes according to its states, namely, susceptible , exposed , quarantined , hospitalized infected , and recovered . Then, . In [19], Mandal et al. constructed an SEQIR model with optimized control and provided advice to Indian cities. Based on their work, Brahim et al. considered a model with Markovian switching and generalized incidence in [25]. Interestingly, they considered that people who are exposed and asymptomatic are infectious. The standard incidence rate is more suitable for infectious disease models with a large population. Based on these results, we consider an SEQIR model with standard incidence as follows:
The constant input of the susceptible population is ; represents the transmission rate of the disease; is the rate of the exposed class removed from the infected class; is the portion of the exposed class moving to the quarantined class; and indicate the rate at which the quarantined people become susceptible and infected types, respectively; the recovery rate of hospitalized infected people is ; is the natural mortality rate; and is the diseased death rate.
A popular technique to incorporate stochasticity into a deterministic model is parametric perturbation. Following this method, the transmission rate, , is replaced by . Therefore, we have a stochastic model defined as follows:
Note that is a standard independent Brownian motion, and the intensity of the white noise is . An initial condition should be satisfied when analyzing systems (1) and (2):
In the following analysis, we declare that the standard Brownian motion in the system is defined on a complete probability space with a filtration satisfying the usual conditions (i.e., it is right continuous and increasing while contains all P-null sets).
3. Analysis of the Deterministic Model
In this section, we analyze system (1) without stochastic perturbation. In this part, some of the basic dynamics are analyzed, including the basic reproduction number and stability of the equilibria.
3.1. The Basic Reproduction Number
An important parameter in infectious diseases is the basic reproduction number which can determine the outbreak or extinction of the disease. Therefore, we give it by the next-generation matrix method [27]. Because the infection classes involve , , and classes, we consider the three classes:
Infection subsystem (4) has transmission matrix and transition matrix , where
So, according to this method, we can get that the spectral radius of is .
3.2. The Equilibria
System (1) has two equilibria, the disease-free equilibrium and the endemic equilibrium , respectively, where
Here, . It is worth noting that the endemic equilibrium exists when .
3.3. The Stability of Two Equilibria
In this part, we study the local asymptotic stability of two equilibria.
Theorem 1. The disease-free equilibrium is locally asymptotically stable if .
Proof. The Jacobian matrix at the disease-free equilibrium of system (1) is given byThe characteristic equation of system (1) at its disease-free equilibrium is given byIt is easy to see that the characteristic values of the Jacobian matrix are negative if and only if . Hence, the disease-free equilibrium of system (1) is locally asymptotically stable if .
Theorem 2. The endemic equilibrium is locally asymptotically stable for .
Proof. The Jacobian matrix at the endemic equilibrium of system (1) is given byThe characteristic equation of system (1) at its endemic equilibrium is given bywhere , , and . It is clear that the characteristic polynomial has two negative eigenvalues, and we study the last cubic polynomial. By calculation, we can verify and . Hence, according to the Routh–Hurwitz criterion, we can know that the endemic equilibrium of system (1) is locally asymptotically stable for .
4. Analysis of the Stochastic Model
Theorem 3. For any initial values , there exists a unique global solution of system (2) for all , and the solution will remain in with probability 1, i.e., for all almost surely.
Proof. Since the coefficients of system (2) are locally Lipschitz continuous on , for any initial value , there exists a unique local solution on , where represents the explosion time [28]. If we can prove a.s, then the solution is global. To this end, let be sufficiently large such that , and all lie in the interval . For each integer , we can define the stopping time [28] through , where, throughout this paper, we set (as usual represents the empty set). Obviously, is increasing as . Let , whence a.s. If we can show , then and a.s for all . That is to say, to complete the proof, all we need is to show that a.s. If this assertion is false, then there exists a pair of constants and such that . There is an integer such that for all . Define a Lyapunov function by .
Since for any , the function is nonnegative. Making use of It’s formula (see [28]) to , we getwhere is a positive constant. Integrating from 0 to and then taking the expectation on both sides, we obtainThus, we haveSet for , and we get . Note that, for every , there exists that equals either . Therefore, is no less than either . Therefore, we can seeThen, we getwhere is the indicator function of . Letting ,is a contradiction, and thus, we can get ., which implies that for all . This completes the proof.
Remark 1. Define . By system (2), we can get . Next, we have . Therefore, we can obtain the solution of system (2) in the region . Then, we can give some general conclusions and results.
Theorem 4. Let be the positive solution of system (2) with the initial value . If the assumption holds, then the disease will go to extinction exponentially with probability one, i.e.,
Proof. Applying It’s formula, we haveIntegrating (18) from 0 to , we getwhere is a continuous local martingale [28] whose quadratic variation isApplying Theorem 7.4 in [28], p44, we can getwhere and is a random integer. Using the Borel–Cantelli lemma [28], we obtain that, for almost all , there exists a random integer such that, for all ,That is to say,Taking it to (19), we obtainObserving the integrand function, we can do the following:Consequently, we getTherefore, for , one can see thatLetting , , , and . We can obtainLetting giveswhich shows that .
On the other hand, there exists a constant such that . Applying the third equation of system (2), we getUsing the comparison theorem, we haveLetting , we get . Similarly, in the case that . and ., we can getThis completes the proof.
Theorem 5. Assume that , , , , , , , and are all positive constants. Then, for any initial value , the solution of system (2) has the property
Proof. System (1) has the endemic equilibrium, and we can getWe construct a functionwhere is a constant. By It’s formula, we obtainAnd we getTaking (37)–(41) into (36), one can obtain thatwhere , , and are defined in Theorem 5.Integrating (39) from 0 to and taking the expectation on both sides we haveThen, letting , dividing on both sides of (44), and taking the superior limit, we get
5. Application of the Model
5.1. Numerical Simulations
The numerical simulations are divided into two aspects. We first set values of the parameters to verify the theorem results obtained in this study. Specifically, we set , , , , , , , , and . Using different initial values for such as , , and , we numerically solve system (1). Figure 1 shows that when , the epidemic will go to extinction, and the solutions with different initial values will converge to the disease-free equilibrium expressed as , which coincides with the conclusion of Theorem 1. For and the other parameters set as before, it is verified that when , the epidemic will continue to spread, and the solutions will converge to the endemic equilibrium expressed as in Figure 1.

If we set and , the other parameters remain unchanged, and the initial value is given by , while these parameters satisfy , the condition of Theorem 4 is satisfied as shown in Figure 2, and in this case, the epidemic will die out. When , , and , the other parameters and initial value are the same as before, and we obtain the asymptotic property of system (2) and the solutions with stochastic oscillation near the endemic equilibrium of system (1), as shown in Figure 3.


5.2. Case Study
In this section, we present some cases in which the proposed research model can be applied. After the WHO reported the first case on December 31, 2019, the COVID-19 pandemic has spread rapidly worldwide. To date, 220 countries have suffered from this pandemic. By May 30, 2021, there were 171,474,925 cases of global confirmed cases, and the cumulative mortality reached 3,565,330 cases. The infectious, noxious, and transmission speeds of the epidemic are shocking. Some countries are still trapped in the epidemic; here, we analyze data from France, Italy, and the UK. The data used in our study were extracted from Worldometers [29], an online freely available repository. In this study, we consider 95 days of time-series data of currently infected patients from three countries from February 24 to May 30, 2021. The source of the parameters is as follows. The population of Italy [30] is approximately 60,500,000, which we set it as the initial value of the total population. The value of can be set on the basis of the number of births in Italy, 1427. The average life expectancy in Italy is 84.01 [31]. Given that is the average life expectancy, the value of can be calculated as . Through the same calculations, we derive some parameters for France and the UK. We also set the range of , which represents the rate of exposure to infection according to [32]. Because the average recovery period in France is approximately 13.6 days [33], the parameter is set to . According to [34], we take the average recovery period in Italy as 13.15 days; thus, in this case, we set the parameter . The recovery rate of the UK, , and are estimated by the least square method.
We fit the parameters of the system by using the daily active COVID-19 cases in France, Italy, and the UK. The parameters resulting from fitting the actual active cases are shown in Table 1; the results are shown in Figure 4. We consider hospitalized infected patients as active cases in this part. Note from Figure 4 that the simulations are correct, and the epidemic development of these three countries shows a decreasing trend. It is evident that the people infection will not stop shortly. The government should take measures to reduce the number of newly infected people and maintain a decreasing trend.

(a)

(b)

(c)
According to the previous analysis, stochastic perturbation can change the solution of the system and the stage of disease development. The epidemic becomes extinct when the noise intensity is sufficiently large. However, if the noise intensity is small, the overall trend of disease development will not change. This situation can be seen in Figure 5, which shows the solutions of the deterministic and stochastic models with , illustrating that the stochastic solution fluctuates around the deterministic solution. This leads to conclude that the government should take some measures to prevent new infections. Figure 6 shows the increase in , which represents the strengthening isolation, and , which means that raising nucleic acid detection will effectively reduce the risk of disease transmission.


(a)

(b)
6. Conclusions
In 2020, new coronavirus pneumonia seriously affected people’s daily life and hindered social development. More than 200 countries worldwide suffered from this epidemic. Prevention of this coronavirus pneumonia has become the most important issue worldwide. Infectious disease models constitute an effective tool for analyzing and understanding the trend of disease spread; therefore, we built a mathematical model with standard incidence based on the spread of COVID-19. We obtained the basic reproduction number, which affects the spread of the disease, and proved that when , the disease-free equilibrium is locally asymptotically stable, and when , the endemic equilibrium is locally asymptotically stable. Then, we incorporated stochastic perturbation in the deterministic model because there exists stochastic influence in the process of disease spread. Specifically, we obtained the conditions under which the epidemic goes to extinction and the asymptotic property near the endemic equilibrium of the deterministic model.
Finally, we verified the results in two aspects. First, we set some parameters to verify the proposed theorem. Second, we fitted the active data cases using the least square method. In Figure 4, we effectively fitted the active cases in France, Italy, and the UK. The active cases of these three countries all decreased. In the first part of the numerical simulations, we took into account that adding a stochastic disturbance will prevent the spread of the disease. Therefore, we considered stochastic disturbance in the model for comparison with the deterministic model. The parameters in Table 1 do not satisfy the condition of disease extinction, so the epidemic will not go to extinction. It is well known that white noise with sufficiently large intensity can lead to the extinction of the disease. Hence, we can control the development of the disease by increasing the noise intensity.
Random perturbance can be considered as a human-issue disturbance. Therefore, governments must take measures to decrease the number of infections, e.g., population lockdown, increase of media publicity, rational allocation of medical resources, and appeal to the public to pay attention to personal hygiene. Some suggestions, such as increasing the isolation and raising nucleic acid detection, should be made. With these suggestions, the number of infections will be reduced. Currently, novel coronavirus pneumonia vaccines have been developed. It is believed that a reasonable combination of government control and vaccination strategies will eliminate the epidemic. Next, we will also consider the model with immunity to analyze the impact of immune intensity on the spread of infectious diseases.
Data Availability
The data used to support the findings of this study are available from the corresponding author upon request.
Conflicts of Interest
The authors declare that they have no conflicts of interest.
Acknowledgments
The authors would like to express their heartfelt appreciation for the financial support provided by the Research Project supported by the National Nature Science Foundation of China (nos. 12071445 and 12001501), Fund for Shanxi (1331KIRT), and the Outstanding Youth Fund of North University of China.