Time­delayed and stochastic effects in a predator­prey model with ratio dependence and Holling type III functional response

In this paper we derive and analyse a novel predator-prey model with account for maturation delay in predators, ratio dependence, and Holling type III functional response. Analysis of the system’s steady states reveals conditions on predation rate, predator growth rate and maturation time that can result in a prey-only equilibrium, or facilitate simultaneous survival of prey and predators in the form of a stable coexistence steady state, or sustained periodic oscillations around this state. Demographic stochasticity in the model is explored by means of deriving a delayed chemical master equation. Using a system size expansion, we study the structure of stochastic oscillations around the deterministically stable coexistence state by analysing the dependence of variance and coherence of stochastic oscillations on system parameters. Numerical simulations of stochastic model are performed to illustrate stochastic ampliﬁcation, where individual stochastic realisations can exhibit sustained oscillations in the case, where deterministically the system approaches a stable steady state. These results provide a framework for studying realistic predator-prey systems with Holling type III functional response in the presence of stochasticity, where an important role is played by non-negligible predator maturation delay.

In this paper we derive and analyse a novel predator-prey model with account for maturation delay in predators, ratio dependence, and Holling type III functional response. Analysis of the system's steady states reveals conditions on predation rate, predator growth rate and maturation time that can result in a prey-only equilibrium, or facilitate simultaneous survival of prey and predators in the form of a stable coexistence steady state, or sustained periodic oscillations around this state. Demographic stochasticity in the model is explored by means of deriving a delayed chemical master equation. Using a system size expansion, we study the structure of stochastic oscillations around the deterministically stable coexistence state by analysing the dependence of variance and coherence of stochastic oscillations on system parameters. Numerical simulations of stochastic model are performed to illustrate stochastic amplification, where individual stochastic realisations can exhibit sustained oscillations in the case, where deterministically the system approaches a stable steady state. These results provide a framework for studying realistic predator-prey systems with Holling type III functional response in the presence of stochasticity, where an important role is played by non-negligible predator maturation delay.
Recent years have witnessed an explosion of interest in different aspects of modelling biological interactions, driven both by new ecological evidence, and by theoretical advances. One particular class of models that have attracted particular attention are models with ratio dependence, where the per-capita rate of predation depends on the ratio between predator and prey populations. Motivated by recent work on plant diseases based on plant-insect interactions, in which insects are the predators feeding on plants playing the role of food source, we have proposed a new predator-prey model with ratio dependence, and a Holling type III (sigmoidal) functional response. Predator population in the model is assumed to be maturing with a non-negligible maturation time, which is explicitly included in the model in the form of a time delay. We explore the role of different parameters and the time delay, and show that the model always supports a prey-only equilibrium, which is stable if the predation rate is sufficiently small, or when the predators take too long to mature. In the opposite case, the model can exhibit coexistence, where both prey and predators are present, and the conditions for biological feasibility and stability of this steady state are established depending on different parameters. The model is also reformulated as a stochastic delayed system to explore the role of demographic stochasticity in the dynamics of predator-prey interactions. We derive equations for stochastic fluctuations around the deterministically stable coexistence steady state and use these to quantify the variance and coherence of these fluctuations. Numerical simulations illustrate how even when the coexistence steady state is deterministically stable, in individual stochastic realisations, the model still exhibits stochastic oscillations around it, which can have major implications a) Electronic mail: k.blyuss@sussex.ac.uk for understanding the dynamics of interactions between real biological species.

I. INTRODUCTION
Ever since the pioneering work of Lotka 1 and Volterra 2 , mathematical models of predator-prey type have provided tremendous insights into dynamics of interactions between different species, or more widely, between interacting agents, that have found application not only in biology, but also in a diversity of other areas, from immunology to economics 3-6 . The starting point for many of many of these models in ecological context is a general predator-prey model of Gause-Kolmogorov type 7-9u where u(t) and v(t) are abundances or population densities of prey and predator, respectively, f (u) describes the intrinsic per-capita growth rate of prey in the absence of predator, and d is predator's natural death rate. In order to guarantee boundedness of prey population, the function f (u) is often chosen in the form of a monotonically decreasing linear function f (u) = r(1 − u/K), which biologically describes an intra-specific prey competition and corresponds to the logistic growth of prey with a linear growth rate r and the carrying capacity K. Function g(u, v) is known as the trophic function or the so-called functional response of the predator 10 , it quantifies how efficiently predators are consuming prey, and how this enhances their own reproduction. To account for the details of the process of predators searching for prey, and more specifically, for the idea that at higher predator densities, predators would have to share some of the prey, Arditi and Ginzburg 11,12 proposed that at the time scales of population dynamics, the overall rate of predation is better represented by a function that depends on the ratio of prey to predators. To this end, they suggested choosing the trophic function g(u, v) in the form g(u/v):u which is known as ratio dependence. Subsequently, besides a large number of theoretical models with ratio dependence [13][14][15][16][17][18][19] , several field and experimental studies 12,20,21 have also provided support for using ratio-dependence predation as a more realistic representation of predator-prey dynamics. While in the classical Lotka-Volterra model, the trophic function g is simply proportional to the number of prey g(u) = au, a more realistic representation that has become one of the most commonly used in ecology was suggested by Holling 22,23 . It accounts for two distinct aspects of interaction between prey and predators, namely, predators searching for prey, and predators handling (i.e. chasing, killing and digesting) the prey. Holling proposed three types of functional response g(u), which all satisfy g(0) = 0, they all approach some constant for large values of u, and the difference is in their behavior for smaller prey numbers/densities. Type I response is linearly increasing for small prey densities, whereas for large prey numbers it saturates at some constant value; type II and type III are also functions that are saturating at high prey numbers, and are, respectively, concave and sigmoidal. In this paper we are interested in Holling type III functional response, which is a sigmoidal function, and is often represented in the form 24,25 g(z) = az n 1 + ahz n , n > 1.
A large number of authors have analysed predator-prey models with this type of functional response and ratio dependence, both without [26][27][28] and with time delays [29][30][31] , as well as with spatial dependence [32][33][34] . Holling type III response with prey dependence and ratio dependence has been observed in a number of experimental settings [35][36][37][38] . In this paper, we will instead consider the following form of the Holling type III functional response g(z) = ae −α/z , which satisfies the conditions of g(0) = 0, is monotonically increasing and settles at a constant value as z → ∞. This functional response is reminiscent of the Ivlev (Holling II) trophic function g(z) = a (1 − e −αz ) 39,40 , as well as of the Ricker model 41 for single-species populations. It is also not entirely dissimilar from a logistic-type term of the form g(u, v) = a(1 − v/u) used in models of plant-insect interactions, where plants serve as hosts for insect predators, and current plant population represents the carrying capacity for insect population 42,43 . With account for ratio dependence, we, therefore, have the trophic function g(u, v) in the form g(u, v) = ae −αv/u . (1) A recent work has used a delayed version of this trophic function to analyse the dynamics of vector-plant interactions in the context of modelling plant mosaic disease 44 .
The last aspect we want to include in our model is the idea that predators normally take some time to mature, and only mature predators can reproduce. Mathematically, this is often represented by stage-structured models that separate populations into immature and mature individuals and include a maturation delay, which describes the period of time it takes for immature individuals to reach maturity, when they start to reproduce [45][46][47][48] . Including maturation delay in the above model, we obtain the following systeṁ which represents a predator-prey model with ratio dependence and Holling type III function response, and which will be studied in this paper. Here, a is the predation or consumption rate, b is the conversion rate from prey into predator biomass, τ is the maturation delay, v(t) represents the population of mature predators, and the factor e −dτ represents the fraction of newly born predators that have survived to maturation. It is easy to show that this model is well-posed in that for any non-negative initial condition, its solutions will remain nonnegative and bounded.
Rescaling the variables and parameters the model (2) can be rewritten as followṡ where we have dropped hats for notational convenience. The remainder of this paper is organized as follows. In Sec. II we find steady states of model (3), derive conditions for their biological feasibility, and study their stability, both analytically and using numerical computations. To understand the role of intrinsic stochasticity of finite population sizes, in Sec. III we develop and analyse a stochastic version of the model, with an emphasis on deriving analytical expressions for spectra of fluctuations around the deterministically stable coexistence steady state. This will be used to numerically compute variance and coherence of stochastic oscillations depending on system parameters. In that section we will also solve the stochastic model numerically and compare its solutions to that of the deterministic model to illustrate the phenomenon of stochastic amplification. The paper concludes in Sec. IV with a discussion of results.

II. STEADY STATES AND THEIR STABILITY
For any values of parameters, the model (3) has a prey-only steady state E p = (1, 0). In the neighbourhood of the point (u, v) = (0, 0) in the first quadrant, the term ve −αv/u is welldefined and positive, and in the limit we have which shows that E 0 = (0, 0) is another steady state of the system that exists for any parameter values, and biologically represents the extinction of both species.
Besides the above two steady states, the model can also have a coexistence steady state whereb = be −τ . In order for this steady state to be biologically feasible, i.e. to have both of its components positive, one has to requireb > 1 (which immediately implies b > 1) and This last conditions can be written as h(z) > 0, where Function h(z) is convex, with two branches going up from a global minimum, which is located at z 0 = ln(b) − 1. Depending on whether z 0 and h(z 0 ) are positive or negative, the function h(z) can have zero, one or two real positive roots. This gives the following restrictions on maturation delay τ, for which the coexistence steady state is biologically feasible where and W 0 (·) and W −1 (·) are Lambert functions with k = 0 (principal branch) and k = −1.
Characteristic equation for linearisation of the system (3) near the steady state E p has the form which shows that one of the eigenvalues is λ = −r < 0, and other eigenvalues satisfy a transcendental equation For τ = 0, the condition for stability, i.e. Re(λ ) < 0, becomes b < 1, which simplifies to b < 1. Let us assume the condition also holds for τ > 0. To check, whether or not the steady state E p is stable, we look for roots of equation (6) in the form λ = ρ + iκ. Substituting this into (6) and separating real and imaginary parts yields Sinceb < 1, the first of these equations has no roots with ρ > 0, which implies that the steady state E p is stable. Ifb > 1, due to the definition ofb, this means that b > 1, and hence, the steady state E p is unstable already at τ = 0, though the analysis we have just performed shows that this steady state can get stabilised for sufficiently large τ, provided the condition (7) holds. For a fixed value of b > 1, the loss of stability of the steady state E p occurs at with E p being stable for τ > τ c , and unstable for τ < τ c . We note that the conditionb > 1 for instability of E p is also one of the conditions for biological feasibility of the coexistence steady state E * . Since parameters in the model (3) have been rescaled using d, whose inverse is the average life expectancy of predators, this suggests an additional biological constraint of τ ≤ 1 to avoid a possibility of maturation time exceeding the overall life expectancy. This restriction implies that for b > e, the prey-only equilibrium E p is unstable within the entire range of possible maturation delays 0 ≤ τ ≤ 1. For the coexistence steady state E * , the characteristic equation has the form (9) For τ = 0, this simplifies to a quadratic equation and with the last term being positive to ensure feasibility of the steady state E * , the condition for stability of E * at τ = 0 can be readily found as Since the coefficients of the characteristic equation themselves depend on τ throughb = be −τ , for τ > 0 it does not prove possible to obtain closed form explicit conditions for stability of the steady state E * . Hence, we compute characteristic eigenvalues numerically using traceDDE 52 .
For the extinction steady state E 0 , we note that the system (3) is not differentiable at this steady state, hence stability analysis of this steady state cannot be performed in a standard fashion. A standard approach to analysis of stability of extinction steady state in systems with ratio dependence is the so-called Briot-Bouquet transformation that replaces the variable v by a new variable z defined as v = zu [49][50][51] . In the particular case of system (3), due to the fact that the ratio-dependent term is instantaneous in one equation and time-delayed in another equation, such a transformation would not remove the singularity in the Jacobian at the point (0, 0). However, since solutions of the system (3) are always non-negative and bounded, and the axes u = 0 and v = 0 are invariant manifolds of the system, one can use Poincaré-Bendixson theorem for time-delayed systems 53,54 to conclude that in the parameter region, where the prey-only equilibrium E p is unstable, and the coexistence steady state E * is infeasible, the system approaches a stable extinction steady state E 0 . Figure 1 shows bifurcation diagram of the system (3) depending on parameters a, b, and the maturation delay τ. In plot 1(a) we observe that since b < e, there is a critical value of maturation delay τ c , as determined by Eq. (8), such that the prey-only equilibrium E p is stable for τ > τ c (region D), unstable for τ < τ c , and at τ = τ c it loses stability through a steady-state bifurcation. For τ < τ c , which biologically means that predators are maturing quite fast, and sufficiently small predation rate a, we have a stable coexistence equilibrium E * (region C), where both prey and predator populations are maintained at some steady levels given by Eq. (4). For higher values of a, as the maturation delay τ is decreased, the steady state E * loses its stability through a supercritical Hopf bifur-cation, resulting in the appearance a periodic solutions (region B). For even smaller values of the time delay (or higher values of a), once the feasibility boundary of the coexistence steady state is crossed, which is also stability boundary for this steady state, the system approaches a stable extinction steady state E 0 (region A). Numerical solutions of the system (2) in each of these parameter regions are illustrated in Fig. 2.
Similar behavior is observed in Fig. 1(b), where due to the fact that now b > e, the prey-only equilibrium E p is unstable in the entire admissible range of τ values. For smaller values of predation rate a and small maturation time τ, we have a stable coexistence steady state E * . If the predation rate is sufficiently small, the coexistence steady state remains stable in the entire range of τ values. However, starting with some value of a, as the maturation time increases, this again results in the loss of stability of E * through a supercritical Hopf bifurcation. In this case the boundary of feasibility of the steady state E * at τ = τ c is outside the admissible range of 0 ≤ τ ≤ 1, once the periodic solution emerges, it remains present for all values of τ up to τ = 1. Increasing the predation rate eventually makes the coexistence steady state biologically infeasible (and unstable), which signifies that the rate of predation is so high that prey is not replenished fast enough to be maintained at some steady level, which, in turn, also drives the predator population down, and the system approaches a stable extinction steady state E 0 .
We also note that as b is reduced to 1, the stability boundary τ = τ c of the steady state E p comes closer to zero, and for b < 1, this steady state is stable for all values of τ, which biologically corresponds to a case, where predator fecundity is too small to maintain predator population, regardless of how quickly it matures. As a result, the coexistence steady state does not exist, and the system always approaches the preyonly equilibrium E p .
Maturation delay τ plays a dual role in the dynamics of the system (3) -it reduces the fecundity of predators by means of replacing b withb, and it also enters the model through delayed values of prey and predator populations in the predator growth term. To explore separate roles of these two effects, in Fig. 1(c) we plot bifurcation diagram of the model in terms of maturation delay τ and a rescaled predator fecundityb. As discussed above, the prey-only steady state E p loses stability once the value ofb exceedsb = 1, at which point the coexistence steady state E * becomes biologically feasible. For lower values ofb, the steady state E * is stable for all admissible values of maturation delay. In contrast, for much higher values ofb, increasing maturation delay τ results in a destabilisation of E * through a supercritical Hopf bifurcation, and the emergence of a stable periodic solution around this steady state. Biologically, this could be interpreted as high fecundity of predators, which results in the higher number of new adult predators, is counterbalanced by the longer time it actually takes for those predators to reach maturity, and it is the balance of those two factors that determines whether the system will exhibit a stable coexistence steady state, observed for smallerb and τ, or a periodic solution around this steady state for higherb and τ. We also note that increasing the predation rate a has the effect of increasing the size of the parameter region, where the coexistence steady state is stable, with instability occurring at higher values ofb and τ. This can be attributed to the fact that increasing a results in a decreased growth of prey population, which means that for the same values of other parameters, lower prey availability facilitates the maintenance of stable coexistence, as is observed in Figs. 1(a) and (b).

III. STOCHASTIC MODEL
In order to explore the role of stochasticity in the dynamics, we begin by constructing a continuous-time Markov chain corresponding to the deterministic system (2). To this end, we introduce X u , X v ∈ {0, 1, 2, . . .} as discrete random variables representing, respectively, the numbers of prey and predators, with the initial condition X(s) = φ (s), s ∈ [−τ, 0]. Interpreting prey and predator populations as chemical reactants, we can derive a delayed chemical master equation (DCME) 55,56 corresponding to system (2), in which we separately account for non-delayed reactions, where there is a single time point, where the update happens both for reactants and the products terms, and the consuming delay terms, where there are two update points: original reactants are updated at the start of reaction, while the products are updated at the end of the reaction [56][57][58] . In the particular case of system (2), logistic growth term of prey, and natural death of predators are nondelayed reactions, whereas predation is a delayed reaction.
If we denote by P(n,t) the probability of finding the system in the state n = (n u , n v ), n u , n v ∈ {0, 1, 2, . . .} at time t, i.e. we have then it satisfies the following DCME where I(n) is the set of all possible past states of the system, from which the state n can be reached via a chain of transitions, P(n,t; m,t − τ) is the joint probability of finding the system in state n at time t and in state m at time t − τ, the shift operators ε ± u,v are defined as follows ε ± u f (n u , n v ) = f (n u ± 1, n v ), ε ± v f (n u , n v ) = f (n u , n v ± 1). To make further progress in the analysis, and in particular, to get a handle on stochastic fluctuations around deterministically stable steady state, we perform the system size expansion of the DCME 59 , which will allow us to separate deterministic and stochastic components. To apply system size expansion to equation (10), we anticipate n u and n v to be of order Ω, with fluctuations of order Ω 1/2 : n u (t) = Ωu(t)+Ω 1/2 ξ 1 (t), n v (t) = Ωv(t)+Ω 1/2 ξ 2 (t), (11) where (u(t), v(t)) are determined by the original deterministic model (2), and (ξ 1 (t), ξ 2 (t)) describe random fluctuations around the deterministic solution. Similarly, we can express delayed variables 60 One can then express probability distributions P(n,t) and P(n,t; m,t − τ) in terms of ξ = (ξ 1 , ξ 2 ) T and η = (η 1 , η 2 ) T as P(n,t) = P(Ω(u, v) T + Ω 1/2 ξ T ,t) = Π(ξ,t), P(n,t; m,t − τ) = Π(ξ,t; η,t − τ), which then implies Expanding this equation in powers of Ω, at the order Ω 1/2 one recovers the system (2) describing deterministic dynamics, and at the next order, i.e. at order Ω 0 one obtains a delayed Fokker-Planck equation that describes stochastic oscillations around the deterministic trajectory (see Appendix A). Following the methodology of Galla 60 (see also [61][62][63][64] ), we can use this delayed Fokker-Planck equation to obtain the following system of Langevin equations describing the dynamics of fluctuations around a deterministic steady state (u * , v * ) of the model (2) where and ζ(t) = (ζ 1 (t), ζ 2 (t)) T is a vector of two independent Gaussian white noise variables with zero mean, and the noise correlators given by Deterministic part of these equations is exactly the same as what can be obtained directly from the system (2) after linearising it near a steady state (u * , v * ), while the noise covariance matrix describing large-scale oscillations around this steady state can only be derived from the system size expansion of the DCME 68,69 . Fourier transform of the system (14) can be found as Introducing the matrix of spectra S(ω) as S i j (ω) = ξ i (ω) ξ j (ω ) , we then have where † denotes Hermitian conjugate, and This then gives the power spectra of fluctuations in prey and predators around the coexistence steady state as The power spectrum of fluctuations in prey population is shown in Fig. 3. One can observe that increasing the predation rate a makes the dominant frequency of oscillations more pronounced, while increasing the maturation delay makes the dominant frequency smaller, though this effect is less pronounced for higher values of a, where the system is closer to stability boundary of the coexistence steady state. When looking at the dynamics of fluctuations around a steady state, the covariance matrix Ξ defined as and it can be found from the matrix of spectra S as follows 65 while for stochastic systems without time delay one could use a Lyapunov equation to obtain the covariance matrix 61,66,67 . The level of fluctuations around the dominant spectral frequency for each of the two populations X u (t) and X v (t) around their steady-state values X * u and X * v can be quantified using their respective mean-square variances Focusing on the particular frequency interval [ω 1 , ω 2 ] around the peak frequency in the distribution of P u,v (ω), one can then compute the quantity and use this quantity to define the coherence of stochastic fluctuations as c u,v = A p u,v /A u,v 61,62,68 . Introducing Y(t) = (Y 1 (t),Y 2 (t)) T as a vector of two continuous random variables, we can derive the following Itô stochastic delay differential equation model where the drift vector µ and the matrix Q are given by and W(t) = [W 1 (t),W 2 (t)] T is a vector of two independent Wiener processes. Equation (17) was solved using the strong predictor-corrector method with the implicitness in the drift coefficient being chosen to be equal to 1/7, which ensures the method has the largest stability region. Figure 4 illustrates the result of a comparison between 20,000 of stochastic realisations obtained by solving the SDDE (17). For chosen parameter values, the coexistence steady state is deterministically stable, with the largest eigenvalues being a pair of complex conjugate eigenvalues with a negative real part. As a result of this, the deterministic trajectory, which coincides with the average of those stochastic simulations, shows the behaviour of decaying oscillations towards a stable steady state E * . In contrast, an individual stochastic realisation can still exhibit stochastic oscillations around this deterministically stable steady state. This phenomenon is known as endogenous 'stochastic resonance', or stochas-tic amplification 68-70 of demographic stochasticity. The fundamental point is that sustained oscillations around the deterministically stable steady state are driven not by external noise, but rather by intrinsic finite size effects. Ecological importance of this result lies in the observation that even in the in the absence of other external drivers, such as seasonality or environmental variability, the system can exhibit sustained stochastic oscillations. We also note that similarly to the deterministic case, predator population in the stochastic simulation also lags behind the prey population.
In Fig. 5 we explore how the variance and coherence of stochastic fluctuations in prey population around the coexistence steady state change with parameters. For b < e, there appears to be little variation in the level of variance of stochastic fluctuations with the value of maturation delay, while increasing predation rate is associated with higher variance. In contrast, for b > e we observe a spilt of the stable parameter region of the coexistence steady state in two parts: for smaller values of maturation delay variance is decreasing with a, while for higher values of maturation it increases with a, particularly in the neighbourhood of stability boundary. In terms of coherence, we observe the split of stability region of coexistence steady state into two parts. The first part is represented by the right region in plot (a) or left region in plot (b), where the largest eigenvalue of linearisation around E * is real and negative, hence deterministically there are now decaying oscillations around the steady state E * , and we observe that there is very little difference in the coherence of stochastic os-cillations with maturation delay, though it does slowly grow with a. The second part of stability region of E * , which corresponds to the left region in plot (a) or right region in plot (b), is characterised by the steady state E * being stable but with the leading eigenvalues being complex. In this case, the system deterministically exhibits decaying oscillations towards this steady state, and there are significantly larger values of coherence, which also demonstrate larger variation with the time delay. As the values of a are increased, and we approach the deterministic boundary of Hopf bifurcation of E * , coherence increases and reaches the value of 1 at the boundary, where stable deterministic oscillations around E * emerge. When plotted in τ-b parameter plane, the part of stability region of E * , where stability eigenvalues for this steady state are complex, is confined to intermediate value ofb and larger values of τ. In this part of the parameter plane, one observes the lowest values of variance and coherence of stochastic oscillations around the deterministically stable steady state. Similarly to what was discussed above, approaching the boundary of Hopf bifurcation, the coherence of stochastic oscillations increases.

IV. DISCUSSION
In this paper we have proposed and analysed a predatorprey model with ratio dependence, Holling type III functional response of Ricker type, and maturation time delay in predators. The model always admits a prey-only equilibrium, and depending on the values of system parameters, it can also have a coexistence steady state with positive values of prey and predator populations. If the predator growth rate is rather small, and the predation rate is large, then there is some minimum maturation delay, for which the coexistence steady state is feasible, while for small predation rates it is always feasible. In contrast, for higher rates of predator growth and intermediate values of predation, the coexistence steady state is only feasible for small and for large values of maturation delay, and we have identified these critical values of time delays required for feasibility in terms of Lambert W function. Stability conditions have been derived for both steady states, and they show that the coexistence steady state is only feasible, when the prey-only equilibrium is unstable. The coexistence steady state can undergo a supercritical Hopf bifurcation, giving rise to stable periodic solutions, and numerical simulations have been performed both to identify stability regions of the coexistence steady state, and to illustrate the dynamics of the model around it.
To explore the role of demographic stochasticity, we have derived a stochastic counterpart of the model, and then applied a van Kampen system size expansion to the delayed chemical master equation to obtain a system of equations describing stochastic fluctuations around the deterministically stable coexistence steady state. Due to linearity of these Langevin equations, it proved possible to obtain spectra of stochastic oscillations in a closed form, and these were then used to study the dependence of variance and coherence of fluctuations on system parameters and the maturation delay. Comparison of numerical solution of Itô SDDE model with the deterministic analogue illustrates that, while in the limit of very large system size the system may approach a coexistence steady state that is deterministically stable, individual stochastic realisations still exhibit stochastic oscillations around this state. This result can be important in the context of applying models of this kind to analysing data emerging from ecological observations of real predator-prey interactions.
Substituting these expressions into the DCME (10), rewriting the left-hand side in terms of Π(ξ,t), and using the relation (13), we collect terms of the same order of Ω. At the highest order, i.e. Ω 1/2 , we havė which is nothing else but the original deterministic model (2) describing macroscopic dynamics.