The spatiotemporal dynamics of a diffusive predator-prey model with double Allee effect

1.   Introduction

The dynamical relationship between prey and predator stands as a pivotal issue in ecology. To unravel the intricacies of biological systems, establishing mathematical models has proven to be an invaluable tool. Inspired by the seminal work of Lotka ^[1] and Volterra ^[2], biologists have ventured into describing predator-prey interactions through the lens of diverse mathematical models. When delving into a modified predator-prey model, the selection of an apt functional response emerges as a crucial aspect, offering a more comprehensive understanding of the intricate interplay between prey and predator. After an exhaustive review of the literature, we chose to focus our research on the Allee effect.

The Allee effect, initially introduced to a predator-prey model by W.C. Allee in 1930 ^[3], offers a realistic portrayal of species interactions. This effect is widely recognized as a pivotal manifestation of beneficial density dependence in sparse populations, which significantly heightens the risk of extinction in low-density populations ^[4]. Notably, two distinct types of Allee effects exist: The strong and weak variants. The strong Allee effect establishes a critical threshold, where the per capita growth rate plummets into negative territory, thereby augmenting the likelihood of extinction. In contrast, the weak Allee effect results in a diminished, yet still positive, growth rate as population size or density dwindles ^[5]. A plethora of mechanisms contributes to the emergence of the Allee effect, encompassing invasion, sex ratio imbalance, environmental alterations, and intricate social interactions, among others ^[6,7]. Numerous studies have underscored the pivotal role played by the Allee effect in shaping population dynamics and exacerbating the extinction risk for low-density populations. Consequently, the implications of the Allee effect within predator-prey systems have garnered extensive research attention ^{[8,9,10,11,12]}. Particularly noteworthy are the findings presented in ^[8], which reveal that models incorporating the Allee effect exhibit a heightened complexity in dynamical behaviors compared to their counterparts devoid of this effect. Furthermore, the author demonstrates that the Allee effect has the capacity to alter the direction of the Hopf bifurcation, thereby bolstering the system's stability. In another study ^[9], the author showcases how the Allee effect can elicit spatially homogeneous bifurcating periodic solutions, with the wavelength escalating in tandem with an increase in Allee effect parameters. Additionally, the predator-prey model has been subjected to various considerations, including the incorporation of noise ^[13], delay ^[14], and diffusion ^[15], among others, to further elucidate its intricate dynamics and implications.

Diffusion is ubiquitous in biological systems, and its incorporation into models often leads to the emergence of diverse pattern formations. These patterns, in turn, underpin complex dynamical behaviors that have been extensively studied in the literature ^{[16,17,18,19]}. To delve into the mechanisms of pattern formation, researchers have developed a range of mathematical tools, including Turing instability analysis, amplitude equation formulations, and multiple scale approaches, among others ^[20,21]. In a notable study ^[22], a mathematical model incorporating the Allee effect was proposed to investigate the Turing instability and pattern formation phenomena in predator-prey systems. This work contributes to a growing body of research that seeks to unravel the intricacies of spatial dynamics in ecological systems. A mini-review on pattern transitions in spatial epidemics has been presented, offering valuable insights into disease prevention and control strategies ^[23,24]. This work highlights the importance of understanding pattern formation not only in ecological systems but also in the context of public health. Moreover, an interdisciplinary approach, combining biological experiments and mathematical modeling, has been demonstrated to provide insights into the mechanisms underlying single-cell wound patterning ^[25]. This approach underscores the value of multidisciplinary collaboration in advancing our understanding of complex systems. In conclusion, pattern formation has emerged as a crucial tool in understanding specific systems, particularly in the context of predator-prey dynamics. By visually representing the distribution of predators and prey, pattern formation provides a powerful means of elucidating the underlying ecological processes and their implications for population dynamics.

To investigate the impact of the Allee effect and diffusion on the dynamic behaviors of predator-prey systems, we first present a predator-prey model with Beddington-DeAngelis functional response as follows.

where $u$ and $v$ describe the population density of prey and predator, respectively, $b > 0$ is the birth rate of the prey, $e$ is the death rate of the prey, $s$ is the death rate of the predator, $h$ is the intensity of competition within the prey population. $\frac{luv}{1+mu+pv}$ is the predatory behaviors of the predator, referred to as the Beddington-DeAngelis functional response and $0 < q < 1$ is the energy proportion obtained by the predator from the prey.

In the realm of predator-prey systems, despite the abundance of research focusing solely on either the Allee effect or diffusion, there remains studies that concurrently consider both of these ecological phenomena within the same framework. This gap in the literature underscores the pressing need for further exploration into the intricate and multifaceted interactions that exist between the Allee effect and diffusion in such systems.

Drawing upon the framework of system (1), Alan J. Terry ^[11] proposed an insightful mechanism for prey resurgence, emphasizing the pivotal role of a strong Allee effect based on specific ecological attributes. Subsequently, the author has delved into the existence and stability of predator-prey model solutions that incorporate a double Allee effect ^[12]. The findings of these investigations underscore the crucial significance of the Allee effect in ensuring the survival of both predator and prey species. The corresponding system, as formulated, is presented below.

where $u$ and $v$ describe the population density of prey and predator, respectively. The biological interpretation of $b, e, s, h, q$ and $\frac{luv}{1+mu+pv}$ is the same as the above model. $\frac{u}{c+u}$ and $\frac{v}{r+v}$ are known as Allee effects, $c$ and $r$ are Allee parameters, in ^[11], the author proved that the predator subject to a strong Allee effect.

However, a crucial aspect overlooked in the aforementioned system is the influence of diffusion. Given the inherent tendency of systems to experience temporal and spatial instabilities due to factors such as migration, immigration, and natural disasters, it is imperative to consider diffusion in our analysis. Based on the above analysis, this article primarily focuses on the following two issues:

(1) In the absence of the Allee effect, how do the dynamical behaviors of system (1) evolve when diffusion is incorporated?

(2) How do the dynamical behaviors of system (1) alter when both the Allee effect and diffusion are simultaneously taken into account?

To tackle these questions, we incorporate diffusion into system (2) and meticulously explore the respective roles played by the Allee effect and diffusion in shaping the dynamical behaviors of prey and predator populations. The organization of this paper is as follows: In the second section, we present a modified model that incorporates diffusion and conduct a rigorous linear stability analysis to derive the sufficient conditions for Turing bifurcation. In the third section, we derive the amplitude equation and conduct a thorough stability analysis of Turing patterns, shedding light on their underlying mechanisms. In the fourth section, we validate our theoretical findings through numerical simulations, ensuring the accuracy and robustness of our conclusions. Last, in the fifth section, we summarize our key discoveries and draw conclusions, highlighting the implications of our research for understanding the complex interplay between Allee effect, diffusion, and population dynamics in predator-prey systems.

2.   The modified model and linear stability analysis

To address the two issues raised in the previous text, we further incorporate the influence of diffusion into system (2) and propose the following model.

where

and $u$ and $v$ describe the population density of prey and predator, respectively, $d_1, d_2$ are diffusion coefficients, $\Omega$ is a bounded connected region with smooth boundaries. The biological interpretation of the remaining parameters in the system is the same as system (1) and (2) described in the introduction.

The above reaction-diffusion equation is subject to both non-zero initial conditions and Neumann boundary conditions, namely

where $\theta$ represents the outward normal vector on the boundary.

Now, we conduct a linear stability analysis at the equilibrium point $(u_{0}, v_{0})$ for Eq (3), which satisfies $f(u_0, v_0) = g(u_0, v_0) = 0$.

Considering a perturbation at $(u_{0}, v_{0})$, we denote the perturbed solution as $(u(t), v(t))$

The corresponding linearized system is as follows

where

where

Assuming $U$ can be written in the following form

where $\lambda_k$ is the growth rate of perturbation, $k$ is the wave number, $r$ is the spatial vector.

The characteristic equation of system (4) is

where

The eigenvalues are

For system (3), the necessary conditions for Turing instability to occur are $Tr_{0} < 0, \delta_{0} > 0$. At this point, the only way for Turing instability to occur is when the system has $\delta_k < 0$ for certain values of $k$. By treating $\delta_k$ as a quadratic function of $k^2$ and finding the minimum value of $\delta_k$ with respect to $k^2$, we can obtain the most dangerous mode amplitude of the system in response to perturbations, namely

By substituting the above equation back into $\delta_k$, the necessary conditions for Turing instability to occur become

Therefore, $(\delta_{k})_{min} = 0$ if $a_{11} = a = \frac{d_1d_2a_{22}+2\sqrt{-d_1d_2a_{12}a_{21}}}{d_{2}^2}$, where $a$ represents the Turing threshold.

Based on the above analysis, we obtain the following theorem:

Theorem. The necessary conditions for Turing bifurcation to occur at the equilibrium point $(u_{0}, v_{0})$ in system (3) are

(ⅰ) $Tr_{0} = a_{11}+a_{22} < 0$,

(ⅱ) $\delta_{0} = a_{11}a_{22}-a_{12}a_{21} > 0$,

(ⅲ)$a_{11}d_{2}+a_{22}d_{1} > 0$,

(ⅳ)$a_{11}d_{2}+a_{22}d_{1} > 2\sqrt{d_{1}d_{2}(a_{11}a_{22}-a_{12}a_{21})}$.

3.   The amplitude equation and the stability of the pattern

Then, just as we know that the system with diffusion is not inherently unstable, the dynamic behaviors should be investigated after the onset of instability. In this section, the properties of secondary order dynamics will be studied using multiple scale analysis.

First, we perform a third-order Taylor expansion of Eq (3) around the equilibrium point $(u_{0}, v_{0})$

where $a_{13} = f_{uu}(u_0, v_0), a_{14} = f_{vv}(u_0, v_0)$, $a_{15} = f_{uv}(u_0, v_0), a_{16} = f_{uuu}(u_0, v_0)$, $a_{17} = f_{vvv}(u_0, v_0), a_{18} = f_{uuv}(u_0, v_0)$, $a_{19} = f_{uvv}(u_0, v_0), a_{23} = g_{uu}(u_0, v_0)$, $a_{24} = g_{vv}(u_0, v_0), a_{25} = g_{uv}(u_0, v_0)$, $a_{26} = g_{uuu}(u_0, v_0), a_{27} = g_{vvv}(u_0, v_0)$, $a_{28} = g_{uuv}(u_0, v_0), a_{29} = g_{uvv}(u_0, v_0)$.

When $a_{11} = a$, Eq.(6) can be abbreviated as

where

where $a_{12}^{*}, a_{21}^{*}, a_{22}^{*}$ correspond to the values of $a_{12}, a_{21}, a_{22}$ when $a_{11} = a$, respectively.

Second, we investigate the dynamical behaviors of Eq.(6) when $u, v$ approach $0$ infinitely, which is also the dynamical behaviors of system (3) after the onset of instability. In this case, $u, v$ could be expanded as

For $\varepsilon$, the system can be written as

For $\varepsilon^2$, the system can be written as

Now, we consider the two systems together

and $|\lambda_k-B+Ck^2| = |\lambda_k-A+Dk^2|^2 = 0$.

Although they share the same dynamics as (5), their pattern formations differ, which means that secondary order dynamics play an important role in the system, and further investigation is necessary.

As we all know, amplitude is a slow variable, and the system in different order can be found in the following.

For $\varepsilon$

For $\varepsilon^2$

For system (9), the solution can be expressed as

where $W_j (j = 1, 2, 3)$ is the amplitude corresponding to the mode of $e^{ik_{j}r}$, $c.c.$ denotes complex conjugate, and we can get $B_2 = \frac{a_{12}^{*}}{a+d_1k^2}$.

The adjoint operator $L^+$ of $L$ is defined as ^[21]

The zero eigenvectors of $L^+$ are

where $A_1 = \frac{a+d_1k^2}{a_{21}^{*}}$.

The amplitude equation can be obtained by the Fredholm solvability condition ^[26]

For simplicity, it can be written as

where $\tau_0 = \frac{B_2+A_1}{B_2}$, $\mu = a-a_{12}^{*}$, $g = \frac{a_{13}B^2_2+a_{14}+2a_{15}+a_{23}A_1B^2_2+a_{24}A_1+2a_{25}B_2A_1}{B_2}$.

Assuming that the amplitude in Eq.(12) has the following form

where $\rho_{j} = |W_{j}|$ and $\phi_{j}$ is the corresponding phase angle.

Substituting $Z_{j} = \rho_{j}e^{i\phi_{j}}$ into (12), we get

where $\phi = \phi_{1}+\phi_{2}+\phi_{3}$.

If we consider only the stable solution of Eq.(14), we have

It is easy to find that the stationary solutions of Eq.(15) are $\rho_1 = \rho_2 = \rho_3 = 0$ and $\rho_1 = \rho_2 = \rho_3 = -\frac{\mu}{|g|}$.

In general, the characteristic matrix is

and the characteristic equation is $|\lambda I-M| = 0$.

If $\rho_1 = \rho_2 = \rho_3 = 0$, the characteristic value is $\lambda = \mu$, when $\mu < 0$, the solution is stable, otherwise the solution is unstable.

If $\rho_1 = \rho_2 = \rho_3 = -\frac{\mu}{|g|}$, the characteristic value can be obtained from $\left(\lambda+\mu \right) \left(\lambda-2\mu \right) ^{2} = 0$, which means it is always unstable.

4.   Numerical simulation

Now, we present the results of computational examples for system (3), utilizing a time step of 0.02 and a space step of 1.

First, we consider the system in the absence of the Allee effect, with the parameter set $s = 0.2, b = 2, c = 0, e = 0.3, h = 0.5, l = 0.6, m = 0.2, p = 0.1, q = 0.6, r = 0, d_1 = 2, d_2 = 0.2$. It's easy to find that the system without the Allee effect exhibits stability [Figure 1(a)], indicating a homogeneous distribution of prey and predator, which subsequently leads to pattern formation Figure 1(b). Subsequently, we alter the coefficients to $s = 1, b = 5, c = 0, e = 1, h = 3, l = 5, m = 1, p = 0, q = 0.5, r = 0, d_1 = 2, d_2 = 0.2$. Under these conditions, the system without the Allee effect becomes unstable Figure 2(a), accompanied by the emergence of spot pattern formation Figure 2(b). By comparing the pattern formations under stable Figure 1(b) and unstable Figure 2(b) states, we observe that the degree of aggregation is higher in the stable state, aligning with natural laws.

Second, we delve into the system that incorporates multiple Allee effects and diffusion, utilizing the parameter set $s = 0.2, b = 3, c = 0.5, e = 0.2, h = 0.4, l = 0.5, m = 1, p = 0, q = 0.56, r = 1$. When $d_1 = 0.01, d_2 = 10$, the system is unstable Figure 3(a), yielding the pattern formation depicted in Figure 3(b). Conversely, with $d_1 = 0.1, d_2 = 1$, the system attains stability Figure 4(a), accompanied by a distinct pattern formation Figure 4(b). Notably, the Allee effect enhances the degree of aggregation in the pattern formations, emphasizing its crucial role in shaping the distribution of both species and boosting the system's survival rate.

Furthermore, we investigate the secondary order dynamics systems (7) and (8) under conditions $s = 0.2, b = 3, c = 0.5, e = 0.2, h = 0.4, l = 0.5, m = 1, p = 0, q = 0.56, r = 1$. The resulting pattern formations exhibit homogeneity Figure 5(a) and 5(b), suggesting that the long-term development of prey and predator populations converges to a steady state.

Last, the pattern formations Figure 6(a) and 6(b) was obtained based on the stability of the amplitude equation. Additionally, our simulations reveal the absence of coexistence solutions when $b, q$ are reduced, implying that slower prey growth and decreased conversion efficiency contribute to instability in prey and predator populations, thereby jeopardizing their survival.

5.   Conclusions

In the seminal work of Alan J. Terry, referenced in ^[11], the Allee effect was introduced into a predator-prey system with a Beddington-DeAngelis functional response, aiming to elucidate its influence on the system's dynamical behaviors. However, the author's analysis notably excluded the impact of diffusion, a crucial aspect in biological systems known for its significant role. Consequently, we extend the foundational framework of Alan J. Terry by integrating the effect of diffusion and meticulously dissecting the distinctive impacts of both the Allee effect and diffusion on the system's dynamical behaviors.

In our theoretical analysis of the system, we introduce a novel methodology for deriving the amplitude equation and secondary-order dynamics. These analytical tools serve as powerful vehicles, enabling us to delve deeper into the subsequent dynamics beyond linear stability analysis.

By bridging this gap in the literature and providing novel insights into the interplay between the Allee effect, diffusion, and predator-prey dynamics, we arrive at the following key observations:

(1) The impact of diffusion on the predator-prey system

The system, devoid of the Allee effect, remains stable under specific parameter sets, yielding a homogeneous distribution of prey and predator. Nevertheless, when diffusion is incorporated, pattern formation emerges, and the system may become unstable under different parameter configurations, resulting in spot pattern formation. Notably, the degree of aggregation is higher in a stable system compared to an unstable one. From a biological standpoint, predators may adjust their foraging behaviors in response to the resource distribution characteristics in spotty patterns. For example, they tend to congregate in 'hotspots' rich in resources to enhance foraging efficiency. Conversely, prey species within spotty patterns may seek refuge in relatively safe areas to evade predation by natural enemies. This spatial selection strategy effectively reduces the risk of predation. Furthermore, as the degree of spatial patterns diminishes, the species distribution within ecosystems becomes more uniform and stable, potentially mitigating ecosystem fluctuations and reducing the risk of collapse.

(2) The impact of Allee effect on the predator-prey system

In the diffusive system devoid of the Allee effect, spot pattern formation occurs, signifying a relatively scattered distribution of species. In contrast, the diffusive system incorporating the Allee effect exhibits more intensive pattern formation, resulting in a relatively higher degree of aggregation compared to the former. From a biological perspective, as patterns become more intensive, species tend to congregate, enhancing reproduction rates and thereby improving resilience against environmental fluctuations, as well as mitigating the risk of species extinction. Furthermore, these intensive pattern formations facilitate interactions among species, strengthening the stability and persistence of the ecosystem.

Additionally, our simulations have revealed that under conditions of slower prey growth and smaller conversion efficiency, no coexistence solution exists, leaving both prey and predator unable to maintain stability and at risk of extinction.

In contrast to most previous studies, which either solely investigated the influence of the Allee effect or diffusion on the dynamic behaviors of predator-prey systems, and often lacked relevant biological explanations, this article comprehensively considers both the Allee effect and diffusion in the predator-prey system. It elucidates their impact on the system from two pivotal aspects: System stability and pattern distribution. These findings underscore the significant influence of diffusion and the Allee effect on species distribution and survival rates for both species, emphasizing the crucial importance of incorporating both factors in prey-predator models for a more realistic analysis and to mitigate extinction risks.

Author contributions

Lingling Li: Formal analysis, Methodology, Writing-review & editing, Project administration. Xuechen Li: Software, Writing-review & editing, Project administration.

Conflict of interest

All authors declare no conflicts of interest in this paper.