Dynamics of a modified Leslie-Gower predator-prey model with double Allee effects

Mengyun Xing; Mengxin He; Zhong Li; Mengyun Xing; Mengxin He; Zhong Li

doi:10.3934/mbe.2024034

Mathematical Biosciences and Engineering

2024, Volume 21, Issue 1: 792-831. doi: 10.3934/mbe.2024034

Previous Article Next Article

Research article

Dynamics of a modified Leslie-Gower predator-prey model with double Allee effects

1.
School of Mathematics and Statistics, Fuzhou University, Fuzhou 350108, China
2.
College of Mathematics and Data Science, Minjiang University, Fuzhou 350108, China

Academic Editor: Xiangsheng Wang

Received: 14 October 2023 Revised: 04 December 2023 Accepted: 11 December 2023 Published: 20 December 2023

In this paper, we investigate the dynamic behavior of a modified Leslie-Gower predator-prey model with the Allee effect on both prey and predator. It is shown that the model has at most two positive equilibria, where one is always a hyperbolic saddle and the other is a weak focus with multiplicity of at least three by concrete example. In addition, we analyze the bifurcations of the system, including saddle-node bifurcation, Hopf bifurcation and Bogdanov-Takens bifurcation. The results show that the model has a cusp of codimension three and undergoes a Bogdanov-Takens bifurcation of codimension two. The system undergoes a degenerate Hopf bifurcation and has two limit cycles (the inner one is stable and the outer one is unstable). These enrich the dynamics of the modified Leslie-Gower predator-prey model with the double Allee effects.

Keywords:

Citation: Mengyun Xing, Mengxin He, Zhong Li. Dynamics of a modified Leslie-Gower predator-prey model with double Allee effects[J]. Mathematical Biosciences and Engineering, 2024, 21(1): 792-831. doi: 10.3934/mbe.2024034

Related Papers:

[1]	Ryo Kinoshita, Sung-mok Jung, Tetsuro Kobayashi, Andrei R. Akhmetzhanov, Hiroshi Nishiura . Epidemiology of coronavirus disease 2019 (COVID-19) in Japan during the first and second waves. Mathematical Biosciences and Engineering, 2022, 19(6): 6088-6101. doi: 10.3934/mbe.2022284
[2]	Avinash Shankaranarayanan, Hsiu-Chuan Wei . Mathematical modeling of SARS-nCoV-2 virus in Tamil Nadu, South India. Mathematical Biosciences and Engineering, 2022, 19(11): 11324-11344. doi: 10.3934/mbe.2022527
[3]	Yukun Tan, Durward Cator III, Martial Ndeffo-Mbah, Ulisses Braga-Neto . A stochastic metapopulation state-space approach to modeling and estimating COVID-19 spread. Mathematical Biosciences and Engineering, 2021, 18(6): 7685-7710. doi: 10.3934/mbe.2021381
[4]	Christian Costris-Vas, Elissa J. Schwartz, Robert Smith? . Predicting COVID-19 using past pandemics as a guide: how reliable were mathematical models then, and how reliable will they be now?. Mathematical Biosciences and Engineering, 2020, 17(6): 7502-7518. doi: 10.3934/mbe.2020383
[5]	Qian Shen . Research of mortality risk prediction based on hospital admission data for COVID-19 patients. Mathematical Biosciences and Engineering, 2023, 20(3): 5333-5351. doi: 10.3934/mbe.2023247
[6]	Antonios Armaou, Bryce Katch, Lucia Russo, Constantinos Siettos . Designing social distancing policies for the COVID-19 pandemic: A probabilistic model predictive control approach. Mathematical Biosciences and Engineering, 2022, 19(9): 8804-8832. doi: 10.3934/mbe.2022409
[7]	Damilola Olabode, Jordan Culp, Allison Fisher, Angela Tower, Dylan Hull-Nye, Xueying Wang . Deterministic and stochastic models for the epidemic dynamics of COVID-19 in Wuhan, China. Mathematical Biosciences and Engineering, 2021, 18(1): 950-967. doi: 10.3934/mbe.2021050
[8]	Enahoro A. Iboi, Oluwaseun Sharomi, Calistus N. Ngonghala, Abba B. Gumel . Mathematical modeling and analysis of COVID-19 pandemic in Nigeria. Mathematical Biosciences and Engineering, 2020, 17(6): 7192-7220. doi: 10.3934/mbe.2020369
[9]	Xinyu Bai, Shaojuan Ma . Stochastic dynamical behavior of COVID-19 model based on secondary vaccination. Mathematical Biosciences and Engineering, 2023, 20(2): 2980-2997. doi: 10.3934/mbe.2023141
[10]	Hamdy M. Youssef, Najat A. Alghamdi, Magdy A. Ezzat, Alaa A. El-Bary, Ahmed M. Shawky . A new dynamical modeling SEIR with global analysis applied to the real data of spreading COVID-19 in Saudi Arabia. Mathematical Biosciences and Engineering, 2020, 17(6): 7018-7044. doi: 10.3934/mbe.2020362

Abstract

1. Introduction

The novel coronavirus SARS-CoV-2 has caused a global pandemic of unprecedented viral pneumonia ^[1,2]. This infection is known as coronavirus disease 2019 (COVID-19) ^[3]. Because of the high human-to-human transmissibility, SARS-CoV-2 has spread rapidly around the world ^[4,5,6,7]. In mainland China, the outbreak started in December 2019, reached the peak in February and then the number of new confirmed cases decreased. On March 18, 2020, there were no new cases of infection for the first time and the economy and daily life gradually returned to normal ^[8,9,10,11]. In Hong Kong Special Administrative Region (SAR), the first case was reported on January 18, 2020. The outbreak peaked in late March 2020 with no new confirmed cases on April 23, 2020. However, spread of COVID-19 in the world is continuing and the outbreak is ongoing globally ^[12,13]. The infection has been confirmed in about 190 countries up to now. As of February 21, 2021, there have been over 110.6 million cases and 2.45 million deaths reported globally since the beginning of the pandemic. The European Region has the largest new cases and new deaths. The United States accounts for the greatest proportion of cumulative cases and deaths ^[14]. In China, although the epidemic has been under control, confirmed cases have been found occasionally in different places, which raised significant concerns on the resurgence of future waves of COVID-19.

If there are no confirmed cases in a region for a long time, then the risk of disease re-emergence might be mainly from imported cases or viruses. The major COVID-19 transmission pathway is from human to human through respiratory droplets ^[15,16]. In particular, asymptomatic individuals who do not have COVID-19 symptoms can still spread the virus. Transmission from asymptomatic individuals poses a significant public health challenge ^{[17,18,19,20,21,22,23]}. The cases imported from other high-risk places are another path of viral spread. To reduce the potential of imported cases, many countries have issued travel restrictions, for example, reducing the frequency of flights from abroad ^[24,25]. However, as the infection is still prevailing in many places, imported cases still represent a tremendous risk, which may lead to new local outbreaks ^{[26,27,28,29]}. Another possible path of SARS-CoV-2 transmission might be through the food supply chain, surfaces and environment. In China, the coronavirus was detected on frozen foods, including their packaging materials and storage environment in July 2020. There seemed to have two outbreaks related to the transmission via frozen food ^[30]. In view of this, interventions that reduce foodborne transmission of pathogens need to be considered ^[32].

Non-pharmaceutical control measures implemented so far are mainly wearing mask, hand washing, social distancing, quarantine and city/region lockdown ^[31]. These interventions were gradually lifted in consideration of the trade-off between economic sustainability and public health. An agent-based model was developed to evaluate the possibility of a second-wave emergence under different extents and timing of intervention relaxation ^[32]. More work assessed the risk of secondary waves since the control measures like lockdowns were relaxed ^{[33,34,35,36,37,38,39,40]}. The study ^[41] compared the data of the epidemiological pattern of COVID-19 in 53 countries or regions where the pandemic experienced two waves, and analyzed the differences between the two outbreaks. Their results suggested that there was a shift of infection to younger age groups, which may make it more difficult to control the pandemic.

In this work, we focus on the COVID-19 spread in serval places where the epidemic has been under control but new cases have been reported occasionally. To study the impact of imported cases on the dynamics of COVID-19 in China under different scenarios of prevention and control measures, Jia et al. developed an impulsive epidemic model to describe imported cases from abroad ^[42]. In their model, the time when the exposed cases were imported was fixed. However, the exposed cases who carry virus without symptoms are usually unknown. When an infected case is identified, the virus has probably been spreading for a period of time. In the beginning of a new wave of epidemic, the infection might be induced by a small number of infected cases. The disease transmission in this stage can be affected by many random factors. In addition, the data of new/accumulated cases were reported on every day. All of these motivate us to develop a stochastic discrete-time compartmental model that considers randomness, epidemic data, as well as the impact of input virus/cases and the initial entry time. By fitting the model to the two waves of outbreaks in two places in mainland China (Beijing and Xinjiang) and Hong Kong SAR, we evaluate the risk factors that can affect the second or future wave of COVID-19.

2. Methods

2.1. Model

We develop a stochastic discrete-time model based on the classic compartmental model. Individuals who have no clinical manifestations such as fever, cough, sore throat and other symptoms that can be self-perceived or clinically recognized, but test positive in serological or blood test are referred to as asymptomatic infection. This population includes two types of individuals. One is asymptomatic infection in the incubation period. They will later develop clinical symptoms or become a confirmed case by screening test or CT (Computed Tomography) examination. The other has no symptoms until the nucleic acid test turns negative. The total population is divided into five epidemiological classes, including susceptible ( $S$ ), exposed ( $E$ ), asymptomatically infected ( $A$ ), symptomatically infected ( $I$ ), and recovered ( $R$ ). Due to quarantine, the susceptible and exposed states are further divided into $S_q$ and $E_q$ . With hospitalization, the infected class (both asymptomatic and symptomatic) can be further divided into $H_A$ and $H_I$ . Because the infection and disease progression can be affected by random factors, we assume that the flow between any two compartments is a stochastic process ^[43,44,45]. For example, $D_{11}(t)$ is the number of susceptible individuals who become newly infected and this process obeys a binomial distribution. The diagram of the model is shown in Figure 1 and the corresponding stochastic discrete-time model is given by the following system:

$\begin{equation} \begin{array}{l} S_{t+1} = S_t-D_{11}(t)-D_{12}(t)+D_{51}(t), \\ E_{t+1} = E_t+(1-q)D_{11}(t)-D_{21}(t), \\ A_{t+1} = A_t+(1-\rho)D_{21}(t)-D_{31}(t)-D_{32}(t)-D_{33}(t), \\ I_{t+1} = I_t+\rho D_{21}(t)-D_{41}(t)-D_{42}(t)-D_{43}(t)+D_{33}(t), \\ S_{qt+1} = S_{qt}+D_{12}(t)-D_{51}(t), \\ E_{qt+1} = E_{qt}+qD_{11}(t)-D_{61}(t), \\ H_{At+1} = H_{At}+D_{31}(t)+(1-\rho)D_{61}(t)-D_{71}(t)-D_{72}(t), \\ H_{It+1} = H_{It}+D_{41}(t)+\rho D_{61}(t)+D_{71}(t)-D_{81}(t)-D_{82}(t), \\ R_{t+1} = R_t+D_{32}(t)+D_{42}(t)+D_{72}(t)+D_{82}(t), \\ \end{array} \end{equation}$

(2.1)

Figure 1. Flow diagram of the COVID-19 infection model.

DownLoad: Full-Size Img PowerPoint

where $D_{ij}(t)$ obeys a binomial distribution $Bin(n, p)$ with the parameters $(n, p)$ , and the specific form is as follows

$\begin{array}{l} D_{11}(t)\thicksim Bin(S_t, P_{11}(t)), \ \ D_{12}(t)\thicksim Bin(S_t, P_{12}(t)), \ \ D_{21}(t)\thicksim Bin(E_t, P_{21}), \\ D_{31}(t)\thicksim Bin(A_t, P_{31}(t)), \ \ D_{32}(t)\thicksim Bin(A_t, P_{32}), \ \ D_{33}(t)\thicksim Bin(A_t, P_{33}), \\ D_{41}(t)\thicksim Bin(I_t, P_{41}), \ \ D_{42}(t)\thicksim Bin(I_t, P_{42}), \ \ D_{43}(t)\thicksim Bin(I_t, P_{43}), \\ D_{51}(t)\thicksim Bin(S_{qt}, P_{51}), \ \ D_{61}(t)\thicksim Bin(E_{qt}, P_{61}), \ \ D_{71}(t)\thicksim Bin(H_{At}, P_{33}), \\ D_{72}(t)\thicksim Bin(H_{At}, P_{32}), \ \ D_{81}(t)\thicksim Bin(H_{It}, P_{43})\ \ D_{82}(t)\thicksim Bin(H_{It}, P_{42}). \end{array}$

$\begin{array}{l} P_{11}(t) = 1-exp[-\beta c(t)\frac{(I+\theta A)}{N}], \ \ P_{12}(t) = 1-exp[-(1-\beta)qc(t)\frac{(I+\theta A)}{N}], \\ P_{21} = 1-exp(-\sigma), \ \ P_{31} = 1-exp(-\delta_A), \\ P_{32}(t) = 1-exp(-\gamma_A), \ \ P_{33} = 1-exp(-k), \\ P_{41} = 1-exp(-\delta_I), \ \ P_{42} = 1-exp(-\gamma_I), \\ P_{51} = 1-exp(-\lambda), \ \ P_{61} = 1-exp(-\delta_q).\\ \end{array}$

Here $exp[-\beta c(t)\frac{(I+\theta A)}{N}h]$ is the probability of staying in the compartment $S$ . The time period $h$ is chosen to be one, so it is omitted in the expression. Thus, $P_{11}$ is the probability of individuals leaving the susceptible compartment. The other $P$ functions can be explained in a similar way. The meaning of each parameter in the model is summarized in Table 1.

Table 1. Parameter estimates for the COVID-19 epidemics in Beijing, Xinjiang and Hong Kong SAR.

Parameters		Definition	Values			Source
Parameters		Definition	Beijing	Xinjiang	Hong Kong SAR	Source
$c(t)$	$c_0$	Contact rate at the initial time	17.0142	14.0025	11.9175	Estimated
	$c_b$	Minimum contact rate under the most strict control strategies	1.0869	2.0468	3.0460	Estimated
	$r_1$	Exponential decreasing rate of contact rate in the first period	0.2196	0.2049	0.0703	Estimated
	$r_2$	Exponential increasing rate of contact rate in the second period	0.0506	0.0200	0.0324	Estimated
	$r_3$	Exponential decreasing rate of contact rate in the third period	0.2241	0.3096	0.1445	Estimated
$\beta$		Probability of transmission per contact	0.2801	0.2722	0.1504	Estimated
$q(t)$	$q_{0}$	Quarantined rate at the initial time	–	0.2820	–	Estimated
	$q_{m}$	Maximum quarantined rate with control strategies	–	0.7083	–	Estimated
	$r_4$	Exponential increasing rate of quarantined rate in the first period	–	0.2027	–	Estimated
	$r_5$	Exponential decreasing rate of quarantined rate in the second period	–	0.1010	–	Estimated
	$r_6$	Exponential increasing rate of quarantined rate in the third period	–	0.2028	–	Estimated
$q$		Quarantined rate	0.2935	–	0.4866	Estimated
$\rho$		Ratio of symptomatic infection	0.5142	0.5520	0.5529	Estimated
$\sigma$		Transition rate of exposed individuals to the infected class	1/5	1/5	1/5	^[47]
$\lambda$		Rate at which the quarantined uninfected contacts were released into the wider community	1/14	1/14	1/14	^[47]
$\delta_I$		Constant transition rate of symptomatic infected	0.3474	0.2008	0.0999	Estimated
$\delta_A$		Constant transition rate of asymptomatic infected	0.2860	0.2006	0.3027	Estimated
$\delta_q$		Constant transition rate of quarantined exposed	0.3599	0.2816	0.2571	Estimated
$\theta$		Correction factor for transmission probability of asymptomatic infectious	0.5919	0.5031	0.5041	Estimated
$k$		Conversion rate from asymptomatic infected to symptomatic infected	0.6124	0.6221	0.5026	Estimated
$\gamma_I$		Recovery rate of infected individuals	0.0701	0.1632	0.0799	Estimated
$\gamma_A$		Recovery rate of asymptotic infected individuals	0.0906	0.1629	0.2393	Estimated
$\tau$		The time of importation of the first case in the second wave	17	16	–	Estimated
$p_E(T_2-\tau)$		The number of exposed cases entered at the time $T_2-\tau$ in the second wave	6	7.02	–	Estimated
$\alpha$		Disease-induced death rate	0	0	0	Assumed
– means the parameter is not included in that place.

| Show Table

DownLoad: CSV

Due to limited pharmaceutical interventions, wearing mask and social distancing play a critical role in the control of the COVID-19 pandemic. As the epidemic is gradually controlled, people's vigilance will decrease. Strict intervention measures may have to be lifted because of economic consideration. We use a time-varying function for the contact rate to describe this change. When the pandemic began and spread rapidly, control measures such as city lockdown, wearing masks and social activity reduction greatly reduced the contact between people. We denote the time of strict control implementation by $T_0$ . When the number of infected cases gradually decrease after the peak, the control measures are relaxed and people's lives gradually return to normal. We denote this time by $T_1$ . When new cases are reported again, people's vigilance increases, and prevention and control measures are implemented again. We denote this time by $T_2$ . The following time-varying function for the contact rate $c(t)$ is used to describe the change of human behavior and effect of control measures during the epidemic.

$\begin{equation} c(t) = \left\{\begin{array}{l} c_{0}, \quad t < T_0, \\ (c_{0}-c_{b})e^{-r_1(t-T_0)}+c_b, T_0\leq t < T_1, \\ (c_{b1}-c_{0})e^{-r_2(t-T_1)}+c_{0}, T_1\leq t < T_2, \\ (c_{01}-c_{b})e^{-r_3(t-T_2)}+c_{b}, t\geq T_2. \end{array}\right. \end{equation}$

(2.2)

Here $c_{b1} = (c_{0}-c_{b})e^{-r_1(T_1-T_0)}+c_{b}$ and $c_{01} = (c_{b}-c_{0})e^{-r_2(T_2-T_1)}+c_{0}$ .

We define the quarantine rate $q(t)$ in a similar way. The quarantine rate increases as the epidemic gets worse and decreases as it improves. Thus, we assume that the quarantine rate is a time-dependent piecewise function, given by

$\begin{equation} q(t) = \left\{\begin{array}{l} q_{0}, \quad t < T_0, \\ (q_{0}-q_{m})e^{-r_4(t-T_0)}+q_m, T_0\leq t < T_1, \\ (q_{m1}-q_{0})e^{-r_5(t-T_1)}+q_{0}, T_1\leq t < T_2, \\ (q_{01}-q_{m})e^{-r_6(t-T_2)}+q_{m}, t\geq T_2, \end{array}\right. \end{equation}$

(2.3)

where $q_{m1} = (q_{0}-q_{m})e^{-r_4(T_1-T_0)}+q_{m}$ and $q_{01} = (q_{m}-q_{0})e^{-r_5(T_2-T_1)}+q_{0}$ . The functions $c(t)$ and $q(t)$ are shown in Figure 2(a–b), Figure 3(e–f) and Figure 4(a–b) for three different places.

Figure 2. The fitting of model (2.1) to the data of COVID-19 in Beijing from January 26 to July 31 in 2020. (a-c): the change of contact rate

$c(t)$ , quarantine rate

$q(t)$ , and the total population

$N(t)$ , respectively. The daily confirmed cases are shown in (d), daily asymptomatic cases are shown in (e), recovered cases are shown in (f). Stochastic fit was performed 50 times with light blue lines. The data and the fitting are represented by black points and the deep blue lines, respectively.

DownLoad: Full-Size Img PowerPoint

Figure 3. The fitting of model (2.1) to the data of COVID-19 in Xinjiang from January 23 to August 18, 2020. The daily confirmed cases are shown in (a), daily asymptomatic new cases are shown in (b), daily asymptomatic cases are shown in (c), cumulative confirmed cases are shown in (d). Stochastic fit was performed 50 times with light blue lines. The data and the fitting are represented by black points and deep blue lines, respectively. (e-g): the change of contact rate

$c(t)$ , quarantine rate

$q(t)$ , and the total population

$N(t)$ , respectively.

DownLoad: Full-Size Img PowerPoint

Figure 4. The fitting of model (2.1) to the data of COVID-19 in Hong Kong SAR from February 23 to July 31, 2020. (a-c): the change of contact rate

$c(t)$ , quarantine rate

$q(t)$ , and the total population

$N(t)$ , respectively. The daily confirmed cases are shown in (d), recovered cases are shown in (e). Stochastic fit was performed 50 times with light blue lines. The data and the fitting are represented by black points and deep blue lines, respectively.

DownLoad: Full-Size Img PowerPoint

2.2. Data

We collected the data of Beijing and Xinjiang from the local health commissions in mainland China, and the data of Hong Kong SAR from the Centre for Health Protection. It includes the time series data of confirmed COVID-19 cases, recovered cases, and asymptomatic coronavirus carriers. On December 26, 2019, a respiratory and critical care physician in Wuhan reported the pneumonia of unknown cause for the first time. The epidemic then spread rapidly in mainland China, and the number of newly confirmed cases reached the peak on February 4, 2020. As of March 18, the number of newly confirmed cases in mainland China became 0 and the number of confirmed cases fell below 20,000. After that, the reported cases in mainland China were mainly imported cases. A few months later, infected cases began to rise again in some places. On June 11, 2020, a confirmed case was reported in Beijing, without history of traveling outside Beijing and without close contact with suspected infection in the past two weeks. This ended a 56-day streak of no local infection in Beijing. On July 15, 2020, i.e., 149 days since the previous confirmed cases, one confirmed case and three asymptomatic cases were found in Xinjiang. In Hong Kong SAR, there were sporadic confirmed cases after April 20. On July 5, a second-wave outbreak emerged. This paper will focus on the data from these three places to study the risk of the emergence of a future wave of COVID-19. The switching time $T_0$ , $T_1$ and $T_2$ in the piecewise function are determined by the responding time in each place.

2.3. Cause of resurgence

If there are no cases for a long period of time, e.g., several months, after a wave of COVID-19 outbreak, then the new infection is likely to be caused by imported cases or exposure to the virus. The virus that caused a second wave can be summarized by the following three sources: (1) imported cases from abroad. Despite strict regulations on international travel and border inspections, there are still some reported cases imported from abroad. There is no guarantee that 100% of the infected or exposed cases entering the country will be isolated. The incubation period of the infection is not well known and may not be the same for all infected people. With fixed-duration quarantine implemented, the infected individual may become a confirmed case after the quarantine is over. This may be a risk for a second wave in mainland China. (2) Asymptomatic cases. These people carry the virus but cannot be identified if they do not have the nucleic acid test. However, they can infect other people. Therefore, asymptomatic carriers represent another risk for the occurrence of the second wave. (3) Virus from the environment. Some studies have shown that low temperature can greatly promote the persistence of coronaviruses. Frozen foods are potential carriers. Transmission occurs via touching contaminated objects that mediate the infection through mouth, nose, or eyes. This seems to be another risk of transmission that have been ignored.

The potential causes summarized above can be described by new exposed individuals added to our model at a certain time. The time point when the new confirmed case was reported is $T_2$ but when the exposed individual was introduced remains unknown. Here we assume that the number of input exposed individuals is $p_E(T_2-\tau)$ where $\tau$ represents the time lag from the entry of the exposed individual to the later confirmation of infection. Thus, $T_2-\tau$ is the time point when the exposed individuals entered. The equation of $E(t)$ in model (2.1) can be replaced by the following equation

$E_{t+1} = E_t+(1-q)D_{11}(t)-D_{21}(t)+p_E(T_2-\tau).$

It is noted that the reported case and the imported case may not be the same person.

The increase in the susceptible population due to lifted interventions may also contribute to the second wave. The first wave of COVID-19 emerged in Wuhan in early January of 2020. The time happened to be about ten days before the Lunar New Year. This made most people stay at home and take the longest vacation, which greatly reduced the probability of contact. In addition, public transportation was terminated and schools and restaurants were all closed. This series of strict measures reduced the number of susceptible people to a very small level. In our model, we assume that the number of susceptible people in the first wave of outbreak is $S_{01}$ . After the first wave, social activities gradually returned to normal and the size of susceptible population increases to $S_{02}$ when the second wave emerges. The time of the susceptible population change, denoted by $T_3$ , depends on the region. For Beijing and Xinjiang, we let it be the same as $T_1$ . For Hong Kong SAR, it is the time when the local restriction policy is released. Thus, the number of susceptible is given by the following piecewise function

$\begin{equation} S_0 = \left\{\begin{array}{l} S_{01}, \quad t < T_3, \\ S_{02}, \quad t\geq T_3. \end{array}\right. \end{equation}$

(2.4)

3. Results

We use the discrete stochastic model (2.1) with the input parameter $p_E(T_2-\tau)$ to fit the data of the two waves of outbreaks in Beijing and Xinjiang using the least square method. The data fitted include the number of reported confirmed cases, asymptomatic cases and recovered cases. For the epidemic in Hong Kong SAR, there were still sporadic reports of confirmed cases after the first wave. The reason for the second wave in Hong Kong SAR is likely the increase in the number of susceptible population due to lifted restriction of interventions. We use the model (2.1) without the input parameter $p_E(T_2-\tau)$ to fit the data in Hong Kong SAR. Parameter values obtained from the fitting are listed in Tables 1 and 2. The population size of susceptible in the three places is less than the entire population of those places. Here the susceptible population refers to those who may contact with the infected cases. The stochastic simulations provide good fits to the data in these three palaces, see Figure 2(d–f), and . The corresponding contact rate $c(t)$ , quarantine rate $q(t)$ and the susceptible population change $S_0(t)$ are shown in Figure 2(a–c), Figure 3(e–g) and Figure 4(a–c), respectively.

Table 2. Initial values for the COVID-19 epidemics in Beijing, Xinjiang and Hong Kong SAR.

Initial values	Definition	Values			Source
Initial values	Definition	Beijing	Xinjiang	Hong Kong SAR	Source
$S_{01}$	The value of susceptible population in the first wave	$5.0014\times 10^3$	$9.4767\times 10^3$	$7.3327\times 10^3$	Estimated
$S_{02}$	The value of the susceptible in the second wave	$6.0119\times 10^3$	$4.9377\times 10^4$	$2\times 10^4$	Estimated
$E(0)$	The initial value of exposed population	8.0747	8.0209	6.0401	Estimated
$I(0)$	The initial value of infected symptomatic population	4.0902	3.0368	6.0316	Estimated
$A(0)$	The initial value of infected asymptomatic population	5.0848	4.0312	2.0328	Estimated
$S_q(0)$	The initial value of quarantined susceptible population	49.7473	49.9653	45.7412	Estimated
$E_q(0)$	The initial value of quarantined exposed population	20.0155	5.0235	13.1394	Estimated
$H_I(0)$	The initial value of confirmed and hospitalized symptomatic population	1	3	5	Data
$H_I(0)$		1	3	5	Data
$H_A(0)$	The initial value of confirmed and hospitalized asymptomatic population	0	1	0	Data
$H_A(0)$		0	1	0	Data
$R(0)$	The initial value of recovered population	0	0	0	Data

| Show Table

DownLoad: CSV

The emergence of the second wave is influenced by the number of input exposed individuals and how long the infection has been spreading before the report of confirmed cases. We conduct numerical simulations to study the risk of having a second wave. The occurrence of a second wave is evaluated by the maximum number of confirmed cases in 500 stochastic simulations. We denote the average number by $MH$ and choose a threshold value 30. If the $MH$ value exceeds 30, it will be regarded as a second wave. The result shows that not all the scenarios result in a second wave. From 500 stochastic simulations, we calculate the probability of the occurrence of a second wave, which is denoted by $Prop$ .

In and , we explore the effect of varying the input parameter on the risk of second wave in Beijing and Xinjiang, respectively. The range of the parameter $p_E$ is set to $[0, 30]$ at time $T_2-\tau$ , and the time delay parameter $\tau$ is within the range $[0, 20]$ . From , we find that both the number of input exposed individuals and the time between initial entry and subsequent confirmation affect the severity of the second wave. The average maximum value of the second wave peak can reach 1600 cases in Beijing. Increasing the number of input exposed individuals can expand the scale of the disease spread. A larger time delay $\tau$ implies that the disease had spread for a longer time without any interventions before its detection. This poses a substantial challenge for the subsequent control of the disease.

Figure 5. Maximum cases and the probability of occurrence of the second wave with different imported cases and time to detection in Beijing. (a, c): the results when imported cases are 1 to 30 from top to bottom. The time delay before detection is 1 to 20 from left to right. (b): The red point means that the maximum cases

$MH$ excess the threshold 30 and the blue point means

$MH$ is below the threshold.

DownLoad: Full-Size Img PowerPoint

Figure 6. Maximum cases and the probability of occurrence of the second wave with different imported cases and time to detection in Xinjiang. The other is the same as that in Figure 5.

DownLoad: Full-Size Img PowerPoint

We provide the parameter region of a second wave occurrence in Figure 5(b). The deep blue points represent the parameter range of the occurrence of a second wave, while the deep red points represent the parameter range of no second wave. The simulation shows that a second outbreak would not take place when less than three exposed cases were imported. If the infection induced by the imported cases can be quickly identified, then the chance of having a second wave decreases. further shows the probability of the occurrence of the second wave under the same parameter range in . Large values of $p_E$ and $\tau$ will make a second wave inevitable. We have the similar conclusion from the simulation for Xinjiang (see Figure 6). The scale of the second wave is larger than Beijing with the same parameter range because the average maximum value of possible second wave peak can reach 5000 cases in the worst scenario.

Figure 9(a) shows the average result of 500 stochastic simulations of the model (2.1) with six different susceptible populations in Hong Kong SAR. As the susceptible population increases, the average maximum value of the second wave peak also increases. Interestingly, the probability of the occurrence of the second wave remains almost the same for different susceptible populations (see Figure 9(b)). Numerical results on the effect of varying the susceptible population size in Beijing and Xinjiang are shown in Figures 7 and 8, respectively. Based on the simulations in these two places, we have a conclusion similar to Hong Kong SAR. This analysis suggests that the susceptible population size plays a minor role in leading to the second wave when the other parameters are fixed.

Figure 7. (a) Simulation of confirmed cases under six scenarios in Beijing when the susceptible population is chosen to be

$0.5S_{01}$ ,

$0.8S_{01}$ ,

$S_{01}$ ,

$S_{02}$ ,

$1.2S_{02}$ and

$1.5S_{02}$ . (b) The probability of occurrence of a second wave under the six scenarios.

DownLoad: Full-Size Img PowerPoint

Figure 8. Simulation of confirmed cases under five scenarios in Xinjiang when the susceptible population is chosen to be

$0.5S_{01}$ ,

$S_{01}$ , middle value of

$S_{01}$ and

$S_{02}$ ,

$1.5S_{02}$ . (b) The probability of occurrence of a second wave under the five scenarios.

DownLoad: Full-Size Img PowerPoint

Figure 9. Simulation of confirmed cases under five scenarios in Hong Kong SAR when the susceptible population is chosen to be

$0.5S_{01}$ ,

$S_{01}$ , middle value of

$S_{01}$ and

$S_{02}$ ,

$1.5S_{02}$ . (b) The probability of occurrence of a second wave under the five scenarios.

DownLoad: Full-Size Img PowerPoint

4. Conclusions

COVID-19, a highly contagious disease first reported in December 2019, has been spreading globally for more than one year. Some countries/regions have mitigated the outbreak by various measures but are still at risk of recurrence. In this paper, we constructed a stochastic discrete-time compartmental epidemic model to analyze the risk of the occurrence of a second or future wave of outbreak. Compared with a deterministic system, a stochastic model is able to include the random factors in the spread of an infectious disease, particularly when the number of initial infected individuals is small. This is the case when a new wave of outbreak occurs. This discrete model can more intuitively describe the flow between any two compartments. The transition between two compartments is not deterministic and assumed to obey binomial distributions in our model. The change between two compartments in one time step corresponds to the daily data. Thus, using the discrete stochastic model facilitates full use of the data from multiple sources, thereby improving the reliability of the parameter estimation results.

To describe the change in the intensity of control measures in response to the COVID-19 pandemic, we adopt time-varying contact rate and quarantine rate in the model. There are a few possible factors that may lead to a second wave, including import exposed cases, asymptomatic cases, and the presence of viruses in the environment such as the frozen food chain. The common characteristic of these factors is that the transmission is silent and difficult to be identified. We find that the time between the exposed case entry and the confirmation of subsequent infection plays a critical role in the occurrence of the second wave.

The cause of the second-wave outbreak in Beijing and Xinjiang is mainly the imported cases and an increase in the susceptible population due to relaxed interventions. The model provided a good fit to the data of the second wave in Beijing in June 2020. Based on the fitting, the value of input exposed cases is estimated to be 6 and the time from exposed individual entry to the detection of infection is 17 days. The size of susceptible population increases from $5.001\times 10^3$ in the first wave to $6.012\times 10^3$ in the second wave. For Xinjiang where the second wave of the epidemic occurred in July 2020, the value of input exposed cases is estimated to be 7 and the time from entry to detection is 16 days. The change in the number of susceptible people is greater than in Beijing.

Hong Kong SAR also experienced a second wave in July 2020. Unlike Beijing and Xinjiang, there were occasional reports of infected cases all the time in Hong Kong SAR after the first wave and the main cause of the second wave is likely to be the increase in the number of susceptible people. Our modeling result suggests that in a region where the infection is not cleared (e.g., in Hong Kong SAR) susceptible people will increase as the control measures are lifted and this may lead to a second wave. If there is no case for a long time (e.g., in Beijing and Xinjiang), it is necessary to screen imported cases and viruses (e.g., via the food chain), which may be the major cause of the second wave.

On the basis of the fitting to the data in Beijing and Xinjiang, we further evaluated the possibility of having a future outbreak and its severity. Because there were no confirmed cases for a long time after the first wave in Beijing and Xinjiang, the contact rate returned to the normal level, as shown in Figures 2 and 3. If there are imported exposed cases, the time to detect the infection is shown to be critical in leading to the second wave. The simulation shown in Figure 5 and Figure 6 indicates that the second wave is determined by the number of imported exposed individuals and the time needed to detect them. The results suggest that if the imported exposed cases are less than three, then the number of confirmed cases will be below the threshold 30 we set, which would not be considered as a second wave. If the values of imported exposed individuals and the time lag in detection are larger (e.g., in the red region in Figures 5 and 6), a second wave will emerge. The more imported exposed cases and the longer for the infection to be detected, the more likely a second wave will occur. Once a confirmed case is found, it is imperative to track the trajectory of that case and the contacted persons. The information obtained from this study can be used to evaluate the possibility (i.e., the possibility of infected cases above a threshold level) and scale of a future wave of outbreak.

By investigating the effect of the susceptible population on the second wave in Beijing, Xinjiang and Hong Kong SAR, we found that the larger the susceptible population size, the more infections if the second wave occurs. However, the susceptible population size itself does not affect the probability of the occurrence of a second wave. This result suggests that imported cases might be an important factor leading to the occurrence of a second wave in a place where the epidemic has been well controlled. Once a case is found, reducing the number of susceptible people will help control the disease spread in the second wave.

Our study cannot predict when a second or future wave of COVID-19 would take place. When a new wave occurs, the model can be used to predict the scale or severity of the outbreak. This is based on the fitting of the model to existing data. If the data are not sufficient for fitting, then the power of the model prediction would be limited. Lastly, the model does not include the influence of vaccination on the disease spread. How the vaccine rollout influences the emergence of future waves remains to be further investigated.

In summary, we established a stochastic modeling framework that incorporates control measures at different stages of the epidemic and potential causes for the second wave emerged in Beijing, Xinjiang, and Hong Kong SAR. Because infected people without symptoms are contagious and the virus attached to goods is difficult to be detected, comprehensive measures are still imperative to curb the COVID-19 pandemic. It is necessary to screen the imported cases in flights and to detect the virus that may be transported by various routes. If a confirmed case is found, the contact of the case should be thoroughly tracked and the close contacts should be quarantined. Finally, it is important to continue protective measures such as wearing masks and avoiding large-scale gathering to reduce the number of susceptible people. This will make the future wave of outbreak less severe if it takes place.

Acknowledgments

This work was finished when the first author visited the University of Florida in 2020. This research was partially supported by the National Natural Science Foundation of China (grant numbers: 12031010(ST), 11631012(ST)) and by the Fundamental Research Funds for the Central Universities (grant numbers: 2018CBLZ001(SH), GK201901008(ST)). L. Rong is supported by the National Science Foundation (grant number: DMS-1950254).

Conflict of interest

No conflict of interest.

References

[1]	P. H. Leslie, Some further notes on the use of matrices in population mathematics, Biometrika, 35 (1948), 213–245. https://doi.org/10.2307/2332342 doi: 10.2307/2332342
[2]	P. H. Leslie, J. C. Gower, The properties of a stochastic model for the predator-prey type of interaction between two species, Biometrika, 47 (1960), 219–234. https://doi.org/10.2307/2333294 doi: 10.2307/2333294
[3]	A. Korobeinikov, A Lyapunov function for Leslie-Gower predator-prey models, Appl. Math. Lett., 14 (2001), 697–699. https://doi.org/10.1016/s0893-9659(01)80029-x doi: 10.1016/s0893-9659(01)80029-x
[4]	S. B. Hsu, T. W. Huang, Global stability for a class of predator-prey systems, SIAM J. Appl. Math., 55 (1995), 763–783. https://doi.org/10.1137/s0036139993253201 doi: 10.1137/s0036139993253201
[5]	T. Lindstrom, Qualitative analysis of a predator-prey system with limit cycles, J. Math. Biol., 31 (1993), 541–561. https://doi.org/10.1007/bf00161198 doi: 10.1007/bf00161198
[6]	I. Hanski, L. Hansson, H. Henttonen, Specialist predators, generalist predators and the microtine rodent cycle, J. Anim. Ecol., 60 (1991), 353–367. https://doi.org/10.2307/5465 doi: 10.2307/5465
[7]	M. A. Aziz-Alaoui, M. D. Okiye, Boundedness and global stability for a predator-prey model with modified Leslie-Gower and Holling-type II schemes, Appl. Math. Lett., 16 (2003), 1069–1075. https://doi.org/10.1016/s0893-9659(03)90096-6 doi: 10.1016/s0893-9659(03)90096-6
[8]	A. F. Nindjin, M. A. Aziz-Alaoui, M. Cadivel, Analysis of a predator-prey model with modified Leslie-Gower and Holling-type II schemes with time delay, Nonlinear Anal. Real World Appl., 7 (2006), 1104–1118. https://doi.org/10.1016/j.nonrwa.2005.10.003 doi: 10.1016/j.nonrwa.2005.10.003
[9]	C. Xiang, J. Huang, H. Wang, Linking bifurcation analysis of Holling-Tanner model with generalist predator to a changing environment, Stud. Appl. Math., 149 (2022), 124–163. https://doi.org/10.1111/sapm.12492 doi: 10.1111/sapm.12492
[10]	Y. L. Zhu, K. Wang, Existence and global attractivity of positive periodic solutions for a predator-prey model with modified Leslie-Gower Holling-type II schemes, J. Math. Anal. Appl., 384 (2011), 400–408. https://doi.org/10.1016/j.jmaa.2011.05.081 doi: 10.1016/j.jmaa.2011.05.081
[11]	M. K. Singh, B. K. Singh, Poonam, C. Cattani, Under nonlinear prey-harvesting, effect of strong Allee effect on the dynamics of a modified Leslie-Gower predator-prey model, Math. Biosci. Eng., 20 (2023), 9625–9644. https://doi.org/10.3934/mbe.2023422 doi: 10.3934/mbe.2023422
[12]	J. Gine, C. Valls, Nonlinear oscillations in the modified Leslie-Gower model, Nonlinear Anal. Real World Appl., 51 (2020), 103010. https://doi.org/10.1016/j.nonrwa.2019.103010 doi: 10.1016/j.nonrwa.2019.103010
[13]	M. M. Chen, Y Takeuchi, J. F. Zhang, Dynamic complexity of a modified Leslie-Gower predator-prey system with fear effect, Commun. Nonlinear Sci. Numer. Simulat., 119 (2023), 107109. https://doi.org/10.1016/j.cnsns.2023.107109 doi: 10.1016/j.cnsns.2023.107109
[14]	W. C. Allee, Animal Aggregations: A Study in General Sociology, University of Chicago Press, Chicago, 1931. https://doi.org/10.5962/bhl.title.7313
[15]	C. Arancibia-Ibarra, The basins of attraction in a modified May-Holling-Tanner predator-prey model with Allee affect, Nonlinear Anal., 185 (2018), 15–28. https://doi.org/10.1016/j.na.2019.03.004 doi: 10.1016/j.na.2019.03.004
[16]	C. Arancibia-Ibarra, J. Flores, Dynamics of a Leslie-Gower predator-prey model with Holling type II functional response, Allee effect and a generalist predator, Math. Comput. Simulation, 188 (2021), 1–22. https://doi.org/10.1016/j.matcom.2021.03.035 doi: 10.1016/j.matcom.2021.03.035
[17]	W. Q. Yin, Z. Li, F. D. Chen, M. X. He, Modeling Allee effect in the Leslie-Gower predator-prey system incorporating a prey pefuge, Int. J Bifurcat. Chaos, 6 (2022), 2250086. https://doi.org/10.1142/s0218127422500869 doi: 10.1142/s0218127422500869
[18]	Y. Z. Liu, Z. Li, M. X. He, Bifurcation analysis in a Holling-Tanner predator-prey model with strong Allee effect, Math. Biosci. Eng., 20 (2023), 8632–8665. https://doi.org/10.3934/mbe.2023379 doi: 10.3934/mbe.2023379
[19]	N. Martinez-Jeraldo, P. Aguirre, Allee effect acting on the prey species in a Leslie-Gower predation model, Nonlinear Anal. Real World Appl., 45 (2019), 895–917. https://doi.org/10.1016/j.nonrwa.2018.08.009 doi: 10.1016/j.nonrwa.2018.08.009
[20]	P. J. Pal, P. K. Mandal, Bifurcation analysis of a modified Leslie-Gower predator-prey model with Beddington-DeAngelis functional response and strong Allee effect, Math. Comput. Simulation, 97 (2014), 123–146. https://doi.org/10.1016/j.matcom.2013.08.007 doi: 10.1016/j.matcom.2013.08.007
[21]	Y. J. Li, M. X. He, Z. Li, Dynamics of a ratio-dependent Leslie-Gower predator-prey model with Allee effect and fear effect, Math. Comput. Simulation, 201 (2022), 417–439. https://doi.org/10.1016/j.matcom.2022.05.017 doi: 10.1016/j.matcom.2022.05.017
[22]	Z. C. Shang, Y. H. Qiao, Multiple bifurcations in a predator-prey system of modified Holling and Leslie type with double Allee effect and nonlinear harvesting, Math. Comput. Simulation, 205 (2023), 745–764. https://doi.org/10.1016/j.matcom.2022.10.028 doi: 10.1016/j.matcom.2022.10.028
[23]	P. Feng, Y. Kang, Dynamics of a modified Leslie-Gower model with double Allee effects, Nonlinear Dyn., 80 (2015), 1051–1062. https://doi.org/10.1007/s11071-015-1927-2 doi: 10.1007/s11071-015-1927-2
[24]	M. T. Alves, F. M. Hilker, Hunting cooperation and Allee effects in predators, J. Theor. Biol., 419 (2017), 13–22. https://doi.org/10.1016/j.jtbi.2017.02.002 doi: 10.1016/j.jtbi.2017.02.002
[25]	S. Zhou, Y. Liu, G. Wang, The stability of predator-prey systems subject to the Allee effects, Theor. Popul. Biol., 67 (2005), 23–31. https://doi.org/10.1016/j.tpb.2004.06.007 doi: 10.1016/j.tpb.2004.06.007
[26]	Z. Zhang, T. Ding, W. Huang, Z. Dong, Qualitative Theory of Differential Equations, Translations of Mathematical Monographs, New Jersey, 1992. https://doi.org/10.1090/mmono/101
[27]	J. C. Huang, Y. J. Gong, J. Chen, Multiple bifurcations in a predator-prey system of Holling and Leslie type with constant-yield prey harvesting, Int. J Bifurcat. Chaos, 23 (2013), 1350164. https://doi.org/10.1142/s0218127413501642 doi: 10.1142/s0218127413501642
[28]	L. Perko, Differential Equations and Dynamical Systems, Springer, New York, 1996. https://doi.org/10.1007/978-1-4684-0392-3
[29]	M. Lu, J. C. Huang, S. G. Ruan, P. Yu, Bifurcation analysis of an SIRS epidemic model with a generalized nonmonotone and saturated incidence rate, J. Differ. Equation, 267 (2019), 1859–1898. https://doi.org/10.1016/j.jde.2019.03.005 doi: 10.1016/j.jde.2019.03.005
[30]	Y. F. Dai, Y. L. Zhao, B. Sang, Four limit cycles in a predator-prey system of Leslie type with generalized Holling type III functional response, Nonlinear Anal. Real World Appl., 50 (2019), 218–239. https://doi.org/10.1016/j.nonrwa.2019.04.003 doi: 10.1016/j.nonrwa.2019.04.003

This article has been cited by:

1.	Timotei István Erdei, Rudolf Krakó, Géza Husi, Design of a Digital Twin Training Centre for an Industrial Robot Arm, 2022, 12, 2076-3417, 8862, 10.3390/app12178862
2.	Cheng-Long Wang, Shasha Gao, Xue-Zhi Li, Maia Martcheva, Modeling Syphilis and HIV Coinfection: A Case Study in the USA, 2023, 85, 0092-8240, 10.1007/s11538-023-01123-w
3.	Yujie Sheng, Jing-An Cui, Songbai Guo, The modeling and analysis of the COVID-19 pandemic with vaccination and isolation: a case study of Italy, 2023, 20, 1551-0018, 5966, 10.3934/mbe.2023258
4.	Hao Wang, Di Zhu, Shiqi Li, Robert A. Cheke, Sanyi Tang, Weike Zhou, Home quarantine or centralized quarantine? A mathematical modelling study on the COVID-19 epidemic in Guangzhou in 2021, 2022, 19, 1551-0018, 9060, 10.3934/mbe.2022421
5.	Israel Edem Agbehadji, Bankole Osita Awuzie, Alfred Beati Ngowi, COVID-19 Pandemic Waves: 4IR Technology Utilisation in Multi-Sector Economy, 2021, 13, 2071-1050, 10168, 10.3390/su131810168
6.	Tianheng Xie, Jianfang Zhang, Xiangtao Li, Data-Driven Intelligent Risk System in the Process of Financial Audit, 2022, 2022, 1563-5147, 1, 10.1155/2022/9054209
7.	Katica Tomic, 2023, Chapter 11, 978-3-031-13752-5, 301, 10.1007/978-3-031-13753-2_11
8.	Shasha Gao, Mingwang Shen, Xueying Wang, Jin Wang, Maia Martcheva, Libin Rong, A multi-strain model with asymptomatic transmission: Application to COVID-19 in the US, 2023, 565, 00225193, 111468, 10.1016/j.jtbi.2023.111468
9.	德玉郭, A Dynamics Model of COVID-19 with Transition from the Asymptomatic to the Symptomatic Infected Individuals, 2023, 12, 2324-7991, 2457, 10.12677/AAM.2023.125248
10.	Kaijing Chen, Fengying Wei, Xinyan Zhang, Hao Jin, Zuwen Wang, Yue Zuo, Kai Fan, Epidemiological feature analysis of SVEIR model with control strategy and variant evolution, 2024, 9, 24680427, 689, 10.1016/j.idm.2024.03.005
11.	Shasha Gao, Pant Binod, Chidozie Williams Chukwu, Theophilus Kwofie, Salman Safdar, Lora Newman, Seoyun Choe, Bimal Kumar Datta, Wisdom Kwame Attipoe, Wenjing Zhang, P. van den Driessche, A mathematical model to assess the impact of testing and isolation compliance on the transmission of COVID-19, 2023, 8, 24680427, 427, 10.1016/j.idm.2023.04.005

Reader Comments

Your name:*

Email:*
© 2024 the Author(s), licensee AIMS Press. This is an open access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0)