Spatial Predator-Prey Model

• The paper presents a spatial predator-prey model that extends classical systems by incorporating diffusion, nonlinear prey movement, and habitat refuges to explain coexistence via bifurcation analysis.
• It employs numerical techniques like finite element discretization and pseudo-arclength continuation to capture pattern formation and track multiple steady states.
• The study demonstrates that spatial heterogeneity and nonlocal interactions fundamentally alter predator-prey dynamics, offering insights for ecosystem management and conservation.

A spatial predator-prey model is a mathematical framework describing the dynamics of interacting consumer (predator) and resource (prey) populations that are distributed over space and allowed to move, reproduce, interact, and die. These models generalize classical non-spatial systems (e.g., Lotka–Volterra, Rosenzweig–MacArthur) to capture the effects of diffusion, directed or density-dependent movement, habitat heterogeneity, and spatial refuges. They provide essential insights into pattern formation, species persistence, coexistence mechanisms, and effects of spatial heterogeneity, which are inaccessible to models lacking spatial dimension.

1. Fundamental Model Structure and Biological Motivation

Spatial predator–prey models are typically formulated as systems of partial differential equations (PDEs) for prey density $u(x, t)$ and predator density $v(x, t)$ over a spatial domain $\Omega \subset \mathbb{R}^d$:

$$\begin{cases} \partial_t u = \nabla \cdot \bigl( D_u(u, x) \nabla u \bigr) + \text{Prey growth} - \text{Predation loss} \ \partial_t v = \nabla \cdot (D_v(x) \nabla v) + \text{Predator gain from predation} - \text{Natural mortality} \end{cases}$$

This structure allows for the representation of:

• Diffusive dispersal (linear when $D_u, D_v$ constant; nonlinear if movement rates depend on density).
• Local demographic processes: logistic prey growth, predator conversion gains, Holling type II or III functional response.
• Habitat structure: spatially variable parameters, e.g., predation rates $b(x)$, impenetrable refuges.
• Variants can include nonlocal interactions, delay (e.g., spatial memory), or stochasticity.

This formulation captures a diverse range of ecological phenomena including spatial aggregation, the formation of prey refuges, density-dependent movement, and spatial patterning, all driven by ecological mechanisms such as crowding (encouraging faster prey dispersal) or habitat fragmentation (Quinones et al., 2019, Quinones et al., 2020, Täuber, 2024, Lv, 2021).

2. Nonlinear Diffusion, Refuge Zones, and Steady-State Bifurcation

A representative model incorporating density-dependent prey diffusion and prey refuge is the spatial Rosenzweig–MacArthur system (Quinones et al., 2019, Quinones et al., 2020). On a domain $\Omega$ containing a refuge $\Omega_0$, the dimensional equations are:

$$\begin{cases} \partial_t u = D_u \nabla \cdot (u \nabla u) + r u(1 - u/\Lambda) - \dfrac{b(x) u v}{1 + m u} & x \in \Omega \ \partial_t v = D_v \Delta v - \mu v + \dfrac{c u v}{1 + m u} & x \in \Omega_1 := \Omega \setminus \overline{\Omega}_0 \ v \equiv 0 & x \in \Omega_0 \end{cases}$$

with Neumann (zero-flux) boundary conditions and nonnegative initial data. Key components:

• Prey diffusion $D_u(u) = D_u u$: higher mobility at greater density (crowding-avoidance).
• Predator diffusion $D_v$ (constant): simple Brownian movement.
• Refuge: $b(x)$ is zero in $\Omega_0$ (no predation), positive elsewhere.
• Holling II functional response: captures predator handling time and saturation.
• Bifurcation analysis: With predator mortality rate $\mu$ as bifurcation parameter, a positive-coexistence steady state emerges via a transcritical bifurcation at:

$\mu_\lambda = \frac{c \lambda}{1 + m \lambda}$

where $\lambda = r\Lambda / D_u$ is the rescaled prey carrying capacity.

At steady state, the system admits either extinction, prey-only occupancy, or positive coexistence, with the structure and number of solutions determined by the aforementioned bifurcation parameter.

3. Effects of Diffusion, Nonlocality, and Hopf/Turing Instabilities

Spatial predator–prey models predict diverse pattern-forming instabilities:

• Turing instability is possible when differences in species diffusion (e.g., $D_v \ll D_u$) destabilize a spatially homogeneous equilibrium, resulting in stationary patterns.
• Nonlocal prey competition: Integral terms representing competition at a distance can induce stable, spatially inhomogeneous periodic oscillations via Hopf bifurcations, even when local terms do not permit such patterns (Chen et al., 2018). The critical parameter regimes for stability and pattern formation are shifted due to the additional feedback from nonlocality.
• Delay and memory/taxis: The inclusion of delays in predation or predator-taxis terms introduces spatial memory, generating oscillatory inhomogeneous modes via Hopf bifurcations absent in local or non-delayed models (Lv, 2021).
• Stochastic and individual-based extensions: Fluctuations and demographic noise lead to erratic oscillations, extinction thresholds, and non-equilibrium phase transitions (e.g., to directed percolation universality) (Täuber, 2024, Tauber, 2011).

A crucial finding is that nonlinear prey diffusion strongly amplifies predator density for small values of $\mu$ (predator-loss rate), even though it does not shift the primary bifurcation threshold for predator persistence when compared to linear diffusion (Quinones et al., 2019).

4. Numerical Methods and Spatial Simulation Insights

Discretization and computational approaches for these models are essential due to the complexity of PDEs with spatial heterogeneity and nonlinear terms:

• Finite difference or finite element spatial discretization: Standard for domains $\Omega \subset \mathbb{R}^2$ with Neumann or mixed boundary conditions.
• Newton–Raphson iteration and continuation: Employed to track steady states and follow branches as parameters (e.g., $\mu$) are varied. Pseudo-arclength continuation handles fold and turning points robustly.
• Stochastic or individual-based simulations: Used when internal reaction noise or spatial discreteness (finite population numbers) are important (Täuber, 2024). These reveal persistent, noisy spatio-temporal structures (activity fronts, clustered patches) not captured by mean-field PDEs.

Time-dependent simulations highlight the stabilizing effects of prey refuges and nonlinear diffusion against both oscillatory and stationary pattern-forming instabilities. In models with density-dependent prey dispersal, prey density can concentrate sharply at the refuge boundary when $\mu$ is low, while linear diffusion leads to more homogeneous spatial profiles (Quinones et al., 2019).

5. Extensions: Multi-Patch, Stochasticity, and Eco-Evolutionary Feedback

Spatial predator–prey frameworks extend beyond homogeneous continuous domains:

• Multi-patch or metapopulation models: Abstract spatial structure as coupled ODEs for discrete habitat patches with migration or dispersal dependent on predation strength or local population densities (Kang et al., 2015, Messan et al., 2015). Such models capture bistability, synchrony/asynchrony of oscillations, and dispersal-induced stability/instability transitions.
• Stochastic partial differential equations (SPDEs): Incorporate environmental noise, spatially heterogeneous coefficients, and diffusion, yielding permanence/extinction criteria dependent on both reaction and noise parameters (Nhu et al., 2018).
• Eco-evolutionary extensions: When predator perceptual range or attack rate is a heritable trait, spatial mixing generates feedback loops mediating species coexistence times and extinction probabilities (Colombo et al., 2019). Such models require tracking spatial and trait distribution dynamics, often via individual-based implementations.

6. Biological and Ecological Interpretations

Key biological conclusions drawn from spatial predator–prey models include:

• Refuges and density-dependent dispersal stabilize coexistence by dampening oscillations and reducing predator-driven prey extinction risk.
• Nonlocal and memory effects can both generate and stabilize complex spatial oscillatory regimes, highlighting the need to incorporate long-range feedback and delay when modeling real ecological systems.
• Spatial heterogeneity and stochasticity are not merely refinements but fundamentally alter dynamical regimes, enabling coexistence at parameter values where deterministic or mean-field models predict extinction.
• Management and conservation strategies must account for the spatial scale of refuges, the strength of density-dependent movement, the possibility of pattern-forming instabilities, and the influence of demographic/environmental noise on critical thresholds for persistence.

In summary, spatialization and biologically motivated movement rules profoundly affect the dynamics, resilience, and diversity-maintenance capacity of predator–prey systems (Quinones et al., 2019, Chen et al., 2018, Täuber, 2024, Lv, 2021, Tauber, 2011, Colombo et al., 2019).