BIOMATHEMATICS

Analysis of Local Stability for Discrete Lotka-Volterra Models

Puma, F. F. C.; Mager, M. M.

DOI 10.5433/1679-0375.2025.v46.53574

Citation Semin., Ciênc. Exatas Tecnol. 2025, v. 46: e53574

Received: August 14, 2025 Received in revised for: October 28, 2025 Accepted: November 20, 2025 Available online: December 16, 2025

Abstract:

This paper applies the stability theory of discrete systems to a predator–prey model with a specific structure, formulated directly in discrete time. This approach offers didactic and computational advantages for modeling ecological systems with non-overlapping generations, contrasting with methods that discretize continuous models. This direct formulation captures the inherently discrete nature of ecological monitoring and non-overlapping generations, while presenting particular analytical challenges. Through linearization and spectral analysis, we obtain explicit stability conditions for the system’s three equilibria: total extinction, which is always unstable; predator exclusion; and coexistence, whose local behaviors depend on conditions among biotic parameters. The results provide practical criteria for predicting population persistence, offering a foundation for applied studies in control and conservation.

Keywords: Lotka-Volterra, mathematical modeling, population dynamics, stability of discrete systems

Introduction

In predator-prey model studies, a pertinent question is: how can we translate the mathematical properties of these equations, often elegant and continuous, into practical answers for ecologists dealing with limited data? This work applies discrete systems stability theory to analyze a predator-prey model formulated directly in discrete time. This approach offers didactic and computational advantages while providing a more natural framework for systems with non-overlapping generations and ecological monitoring at fixed intervals, contrasting with common approaches that discretize continuous models. This sensitivity to discretization is well documented in recent studies, where numerical step size and parameter values may induce qualitative changes in the system dynamics ([axioms - not found]).

This study stems from the interest in deepening theoretical aspects and applications of discrete models in population ecology; this approach involves the qualitative study of difference equations ([diniz2011 - not found]). Although ecological processes occur in continuous time, practical research in population dynamics often deals with discrete data: seasonal censuses, annual monitoring, or non-overlapping generations. This discrete modeling has historical roots: the logistic equation of (Verhulst et al., 1838), proposed as a counterpoint to the unlimited growth of (Malthus et al., 1826), was later popularized by (May et al., 1976) in its discrete version, which revealed its potential to model complex behaviors such as bifurcations and chaos in simple systems.

Contemporary studies have significantly expanded this field: (Braverman et al., 2006) established rigorous stability criteria for perturbed Ricker discrete models, while (Seno et al., 2008) demonstrates how management strategies can lead to counterintuitive effects in discrete population dynamics. These advances have not only consolidated the theory of one-dimensional systems, but have also paved the way for the analysis of more complex interactions. More recent approaches have also investigated stability in Lotka–Volterra systems from global and structural perspectives, including the use of invasion graphs and geometric methods ([Almaraz2024Structural - not found]).

Additionally, discrete predator–prey models incorporating non-linear functional responses and semi-discretization techniques have been studied in connection with stability and bifurcation phenomena ([Lv2024 - not found]). In the national context, (Luiz et al., 2022t) investigate a predator–prey model formulated as a system of partial differential equations of Telegraph type and analyze the stability and convergence of finite-difference discretization schemes. Although our model is formulated directly in discrete time, this work is included to illustrate how discretization of continuous predator–prey systems may affect stability properties and numerical behavior.

In this context, discrete versions of the Lotka-Volterra model emerge as powerful tools for capturing predator-prey dynamics in systems with discrete temporal observations ([Din2013 - not found]; [Khaliq2022 - not found]). While various discrete formulations exist in the literature, it is crucial to distinguish the approach adopted here. Previous works, such as (Din et al., 2013), often investigate models incorporating logistic growth terms for both predator and prey populations.

In contrast, this work adopts a distinct conceptual framework. For the prey, we use the classical form of the discrete logistic equation, as presented by Smith in (Fisher et al., 1979), suitable for populations with non-overlapping generations. For the predators, we assume a simpler growth model, without a logistic term, where the dynamics are governed essentially by the intrinsic mortality rate and predation success. This choice reflects ecological scenarios where predators, being typically less abundant and more mobile, do not experience strong enough intraspecific competition to impose an explicit carrying capacity in the short term.

Thus, the discrete Lotka-Volterra model we study here serves as a tool to analyze observable scenarios in the interaction of two species: \(x_n\) (prey) with logistic growth and \(y_n\) (predator):

\[ \begin{cases} x_{n+1} = x_n \left[1 + r \left(1 - \frac{x_n}{K}\right)\right] - b x_n y_n, \\ y_{n+1} = (1 - c) y_n + d x_n y_n, \end{cases}\tag{1}\]

where \(r\) is the prey growth rate, \(K\) is the carrying capacity, \(b\) and \(d\) are species interaction parameters, and \(c\) is the predator decay rate.

All biotic parameters are strictly positive, and in the absence of predators (\(y_n = 0\)), the system of equations (1) reduces to the discrete logistic model presented by Maynard Smith in 1968, applied to populations with non-overlapping generations (Fisher et al., 1979).

Figure 1 shows the temporal evolution of these populations, obtained by numerical simulation of the discrete-time system (1). The simulation was performed by directly iterating the system’s equations with initial conditions \(x_0 = 100\) and \(y_0 = 6\). The parameter values \(r = 1.4\), \(b = 0.01\), \(c = 0.05\), and \(d = 0.001\) were selected to satisfy the coexistence equilibrium conditions derived in this work, and the simulation was run over 200 iterations.

image

Figure 1 - Oscillations of the populations \(x_n\) and \(y_n\) in the predator-prey system (1).

The notion of equilibrium state is fundamental to the study of discrete systems dynamics. Across various fields such as biology, economics, physics, and engineering, it is often desirable that all states or solutions of a discrete system converge to an equilibrium state. This behavior lies at the heart of stability theory, a field of great importance to scientists and engineers. In this context, a system is considered stable if, when subjected to a small change in its initial conditions, it does not deviate indefinitely but rather maintains or returns to predictable behavior.

The objective is to establish conditions among the biotic parameters that imply stability of the equilibria for system of equation (1), which are identified as:

A systematic theoretical foundation is then established to support future applied investigations. The main contributions are as follows:

The results, obtained via linearization and spectral analysis of the Jacobian matrix, provide a mathematical methodology for interpreting dynamics observed in real systems with discrete sampling. This is the first step for future applications in control and management. In summary, to analyze the stability of the discrete Lotka-Volterra model, equation (1), we adopt a three-step approach: identification of the equilibria, linearization of the system around these points, and application of stability criteria based on the eigenvalues of the Jacobian matrix. This strategy is typical for dynamical systems (Saber et al., 2007) and (Krabs et al., 2010). The results reveal a subtle picture: while the equilibrium \(E_1\) is unconditionally unstable, the equilibria \(E_2\) and \(E_3\) require fine adjustments between growth and interaction to guarantee stability.

Materials and methods

In the study of stability for nonlinear discrete systems, a common approach is the linearization of the system around equilibrium points. This method involves approximating the nonlinear equations by linear equations, which simplifies the analysis of the system’s local behavior. Thus, stability analysis through linearization becomes an essential tool for understanding the behavior of nonlinear systems and evaluating their ability to return to an equilibrium state.

The theoretical foundation that follows, necessary for the local stability analysis, is primarily based on the works of (Krabs et al., 2010); (Saber et al., 2005); (Saber et al., 2007). The definitions and theorems have been adapted in notation and scope to better suit the objectives of this research, particularly for the analysis of system (1). The following definitions and theorems consist of classical results from discrete dynamical systems theory, reformulated in a notation that better suits the context and objectives of the present work.

Local stability and linearization

To analyze the local stability of the system, we first recall key mathematical concepts and formal definitions that will be used throughout work, in which a basic and standard concept in difference equations theory is introduced, for which equivalent formulations can be found in the introductory chapters of (Saber et al., 2005); (Saber et al., 2007), as well as a classical concept adapted from (Saber et al., 2005). Moreover, the linearized system is formalized following the standard procedure exemplified in (Saber et al., 2005); (Saber et al., 2007), and the classical concept of a saddle point in continuous systems is adapted to discrete systems, following the intuition established in (Saber et al., 2007) for hyperbolic points.

Furthermore, we state a theorem essential for the analysis of stability; for the complete proofs of the constituent results, see Theorem 1.5 and its corollary in (Krabs et al., 2010), pp. 47, 49.

Definition 1 (Equilibrium Point). An element \(E\) in the domain of a function \(f\in C^1\) (continuously differentiable on an open set) is called an equilibrium point or state for the system \(x_{n+1}=f(x_n)\) if

\[ \begin{array}{rcl} E &=&f(E). \end{array}\tag{2}\]

Equation (2) states that \(E\) is a fixed point of the function \(f\). Thus, an element \(E=(x^*,y^*)\in \mathbb{R}^2\) is an equilibrium point for system of equations (1) if,

\[ \left\{ \begin{array}{rcl} x^* \left(r\left( 1-\frac{x^*}{K}\right)-by^*\right)&=&0, \\ y^*\left(-c+dx^*\right)&=&0. \end{array}\right.\tag{3}\]

By solving system of equations (3), we identify the equilibria:

Being compatible with biological hypotheses, we consider positive equilibria.

Local Stability means that, if we take an initial value \((x_0,\,y_0)\) sufficiently close to \(E=(x^*,y^*)\), the subsequent iterates \((x_n,\,y_n)\) will remain close to \((x^*,\,y^*)\).

Definition 2. Let \(E\) be an equilibrium point of the system \(x_{n+1} = f(x_n)\).

  1. \(E\) is called locally stable if, given \(\varepsilon > 0\), there exists \(\delta = \delta(\varepsilon) > 0\) such that if \(||x_0 - E|| < \delta\) then \(||x_n - E|| < \varepsilon\) for all \(,\, \forall n \in \mathbb{N}\).

  2. \(E\) is called unstable if it is not stable.

  3. A locally stable equilibrium \(E\) is called asymptotically stable if there exists \(\eta > 0\) such that if \(||x_0 - E||~<~\eta\) then \(\displaystyle\lim_{n\rightarrow\infty}||x_n - E|| = 0\).

To analyze the local stability of nonlinear systems of equations (1), it is appropriate to consider linearization around equilibrium points, which simplifies the analysis by approximating nonlinear equations with linear equations. This technique is essential for understanding the system’s local behavior and its ability to return to a predictable state.

Using vector notation, the system of equations (1) can be written as

\[ (x_{n+1},y_{n+1})=\left(f_1(x_n,y_n),f_2(x_n,y_n)\right),\tag{4}\]

as shown in equation (4), where

\[ f_1(x,y)= x \left[1 + r \left(1 - \frac{x}{K}\right)\right] - b x y \quad \mbox{and} \quad f_2(x,y)=(1 - c) y + d x y.\tag{5}\]

These component functions, defined in (5), capture the prey and predator dynamics, respectively. Thus, the Jacobian matrix of the vector field \(f(x,\,y)\) is defined as

\[ \mathbf{J}_f(x,y)=\left[ \begin{array}{cc} \frac{\partial f_1}{\partial x}(x,y) & \frac{\partial f_1}{\partial y}(x,y) \\ \frac{\partial f_2}{\partial x}(x,y) & \frac{\partial f_2}{\partial y}(x,y) \end{array}\right].\tag{6}\]

The Jacobian matrix (6) is central to our stability analysis. Next, we define the linearized system.

Definition 3 (Linearized System). Let \((x^*,y^*)\) be an equilibrium of system (1). The Linearized System around this equilibrium is defined by the following linear system:

\[ \begin{cases} x_{n+1}&=\ a_{11}x_n+a_{12}y_n\\ y_{n+1}&=\ a_{21}x_n+a_{22}y_n \end{cases}\tag{7}\]
where \([a_{ij}]_{2\times 2}=\mathbf{J}_f(x^*,y^*)\).

If \((x^*,\,y^*)\) is an equilibrium point under the conditions of Definition 3, and denoting by \(\mathbf{X}_n=(x_n,\,y_n)\),   \(\mathbf{J}=\left[a_{ij}\right]_{2\times 2}\), we can rewrite the Linearized System as \[\mathbf{X}_{n+1}\ = \ \mathbf{J}\mathbf{X}_n, \tag{8} \] where

\[ \mathbf{J} = \begin{bmatrix} 1+r(1-\frac{2x^*}{K})-by^* & -bx^*\\ dy^* & 1-c+dx^* \end{bmatrix}. \tag{9}\]

The local stability of an equilibrium can be determined by analyzing the linearized system. The Theorem 1 provides a criterion based on the eigenvalues of the Jacobian matrix (9).

Theorem 1 (Spectral Criterion for Local Stability). Let \(X\) be an open subset of \(\mathbb{R}^k\) and \(x^* \in X\) a fixed point of a function \(f: X \to X\), with \(f\in C^1\), and let \(J_f(x^*)\) be the Jacobian matrix of \(f\) evaluated at \(x^*\). Denote by \(\varrho(J_f(x^*))\) the spectral radius of \(J_f(x^*)\), i.e., the largest modulus of the eigenvalues of \(J_f(x^*)\). Then, the following statements hold:

  1. (Asymptotic Stability) If \(\varrho(J_f(x^*))<1\), and \(x_{n}=f^n(x_0)\), then \(\displaystyle \lim_{n\rightarrow\infty}{||x_n-x^*||}=0\). This establishes an asymptotically stable equilibrium.

  1. (Source-Type Instability) If \(J_f(x^*)\) is invertible, \(x_n=f^n(x_0)\), and every eigenvalue of \(J_f(x^*)\) has modulus greater than \(1\), then there exists \(\delta>0\) such that if \(||f^k(x_0)-x^*||<\delta\) for \(0\leq k\leq n-1\), then \(\displaystyle \lim_{n\rightarrow\infty}{||x_n-x^*||}=+\infty\). This establishes an unstable equilibrium point of source type.

Theorem 1 is a restatement of the classical result from (Krabs et al., 2010). This approach is more natural for discrete systems in \(\mathbb{R}^k\), and is directly applicable to ecological contexts where the behavior of iterates \(x_n=f^n(x_0)\) provides clear dynamic interpretation regarding equilibria.

Definition 4 (Saddle Point). If \(J(E)\) has eigenvalues \(\lambda_1, \lambda_2\) with \(|\lambda_1| < 1 < |\lambda_2|\), then \(E\) is an unstable equilibrium of Saddle Point type.

This definition is consistent with the previous criteria based on spectral radius. Indeed, it suffices for a single eigenvalue with modulus greater than 1 to make the equilibrium point unstable, even if the other eigenvalue lies within the unit disk. The saddle point case precisely illustrates this situation: the presence of a dominant eigenvalue \(|\lambda_2|>1\) guarantees local instability of the equilibrium. A broader theoretical justification can be provided by the Hartman-Grobman Theorem, which in the context of discrete systems ensures the validity of linear analysis around hyperbolic points of nonlinear systems (Saber et al., 2007).

Results and discussion

With the theoretical foundation established, we now turn to the study of the equilibrium points’ behavior in the discrete Lotka-Volterra system. An equilibrium point \(E=(x^*,\,y^*)\) of system (1) is, by definition, a solution of system 3. Next, we will analyze each of these points to better understand their local characteristics.

We emphasize that our analysis will be conducted from the perspective of hyperbolic equilibrium points, that is, those for which the eigenvalues of the associated Jacobian satisfy \(|\lambda_i|\neq 1\). This restriction allows us to employ classical local stability criteria, since non-hyperbolic cases require more in-depth study.

Trivial equilibrium

The equilibrium \(E_1=(0,0)\) corresponds to the scenario of complete extinction, where both populations are absent. The Jacobian matrix, equation (9), at this point takes a simple diagonal form:

\[ \mathbf{J}^*_1= \begin{bmatrix} 1+r\ & 0\\ 0 & 1-c \end{bmatrix} \tag{10} \]

whose eigenvalues are \(\lambda_1=1+r\)   and   \(\lambda_2=1-c\).

Since \(r>0\) by biological hypothesis, it follows that \(|\lambda_1|=1+r>1\). This condition, as previously explained, already allows us to conclude that the equilibrium \(E_1\) is unstable, regardless of the predator mortality rate \(c\). The nature of this instability depends on the behavior of \(\lambda_2\), that is, on the parameter \(c\). Thus, the instability can be classified according to two criteria:

  1. When \(c > 2\): Both eigenvalues of the Jacobian matrix, equation (10), have modulus greater than one, characterizing \(E_1\) as an unstable source according to item b) of Theorem 1.

  2. When \(0 < c < 2\): We have \(|1-c| < 1\), thus \(|\lambda_2| < 1 < |\lambda_1|\) reveals a saddle point, following Definition 4.

 

Predator extinction equilibrium

This equilibrium corresponds to prey survival at their carrying capacity \(K\), with predator extinction. The Jacobian matrix, equation (9), at \(E_2=(K,\, 0)\) is upper triangular:

\[ J^*_2= \begin{bmatrix} 1-r & -bK \\ 0 & 1-(c-dK) \end{bmatrix} \tag{11} \]

whose eigenvalues are \(\lambda_1=1-r\)  and   \(\lambda_2=1-(c-dK)\).

Note that in this case, the eigenvalues of the Jacobian matrix, equation (11) depend on a more complex interaction between the model parameters. Beyond \(r\) and \(c\), they explicitly involve \(b\), \(d\) and the carrying capacity \(K\), making the analysis more sensitive to system variations. This richer parametric structure requires careful examination of possible combinations influencing the local stability of equilibrium \(E_2\).

Depending on these parameter combinations, three types of local behavior can occur: Asymptotic stability, Source-type instability, and Saddle point.

Asymptotic stability

Applying item (a) of Theorem 1 we obtain that \(E_2\) is an asymptotically stable equilibrium if, \[0 < r < 2 \quad \text{and} \quad dK < c < 2+dK.\]

From a biological perspective, this means the populations remain stable at the point \(E_2=(K, 0)\) when the prey growth rate \(r\) is moderate and predators are ineffective at invading the system. That is, the predator mortality rate exceeds the benefits gained from predation at maximum prey density. Note also that the carrying capacity lies in a narrow interval \(K\in\left( \frac{c-2}{d}, \, \frac{c}{d} \right)\).

Figure 2 illustrates the asymptotic stability at the predator-free equilibrium \(E_2 = (K, 0)\). To satisfy the stability conditions for \(E_2\) established in the analytical section (\(0 < r < 2\) and \(dK < c < 2 + dK\)), a different parameter set (\(r = 1.2\), \(b = 0.015\), \(c = 1.45\), \(d = 0.01\), \(K = 120\)) was used. The results were obtained through numerical simulation of system of equations (1), with parameters selected to meet these specific analytical conditions.

image

Figure 2 - Asymptotic stability at \(E_2=(120, \, 0)\).

Figure 2 shows that, as Trajectories 1-4 approach sufficiently close to \(E_2\), they begin to converge toward the equilibrium.

Source-type instability

Considering item (b) of Theorem 1 we obtain that \(E_2\) is a source-type unstable equilibrium if, \[r > 2 \quad \text{and} \quad (c - dK < 0 \quad \text{or} \quad c - dK > 2).\]

Let us consider two relevant scenarios:

In Figure 3, three systems were configured to demonstrate different types of local instabilities at \(E_2~=~(120, 0)\), using the parameter values \(K = 120\), \(b = 0.015\), and \(d = 0.01\), and considering 250 iterations:

image

Figure 3 - Three systems configured to illustrate different local instabilities at \(E_2=(120, \, 0)\).

The source-type instability configuration is represented by the red Trajectory 1 in Figure 3, where we observe the progressive divergence from the equilibrium point \(E_2 = (K, 0)\).

Saddle Point

We verify that this type of instability manifests in two distinct regimes:

Figure 3 illustrates the trajectories of three different systems, obtained through numerical simulation with parameters selected for each instability regime. Each trajectory is initiated near the equilibrium \(E_2 = (120, 0)\) and configured to represent the instability regimes discussed previously. The source-type instability is represented by the red Trajectory 1, which shows progressive divergence from the equilibrium point \(E_2\). In Trajectory 2 (blue), we visualize the predator invasion scenario, where the system begins near \(E_2\) and evolves away from this equilibrium, allowing predator population growth. Finally, Trajectory 3 (green) displays the case where prey exceed \(K\), with predators remaining extinct (or near-zero) while prey density shows abrupt variations above and below \(K\) at each iteration.

Coexistence equilibrium

This equilibrium \(E_3=(x^*, y^*)\) corresponds to the coexistence of both species in the same environment, with positive populations of prey and predators. To determine its coordinates, we start from system of equations 3, assuming \(x^*>0\) and \(y^*>0\), which allows us to rewrite it as:

\[ \begin{cases} 0&= r\left( 1-\dfrac{x^*}{K}\right)- by^* \\ 0&= dx^*\ -\ c \end{cases} \tag{12} \]

The second equation of system, equation (12), directly yields

\[x^*=\frac{c}{d}.\]

Substituting this value into the first equation of system, equation (11), we obtain

\[y^*=\frac{r}{Kb}\left(K-\frac{c}{d}\right).\]

Thus, the coexistence equilibrium point is given by

\[E_3=\left(\frac{c}{d}, \frac{r}{Kb}\left(K - \frac{c}{d}\right)\right) , \ K>\frac{c}{d}.\]

Note that \(E_3\) will have both coordinates positive whenever \(Kd-c>0\), guaranteeing biologically viable populations. Proceeding with the local analysis, the Jacobian matrix, equation (8) evaluated at \(E_3=(x^*, y^*)\) takes the form:

\[ \mathbf{J}^*_3=\begin{bmatrix} 1-\dfrac{rc}{dK} & -\dfrac{bc}{d} \\[0.15cm] \dfrac{rd}{bK}\left(K-\dfrac{c}{d}\right) & 1 \end{bmatrix} \tag{13} \]

Although the eigenvalues of the Jacobian matrix, equation (13), can be obtained analytically, the resulting expression, given by

\[\lambda=\dfrac{2-\dfrac{rc}{dK} \pm \sqrt{\left(2-\dfrac{rc}{dK}\right)^2-4\left(1-\dfrac{rc}{dK}(1+c-dK) \right) }}{2}. \tag{14} \]

is algebraically intricate and not very enlightening. We therefore reserve equation (13) for the analysis of the specific Case II, \(\dfrac{rc}{dK} = 2\).

To further analyze the coexistence equilibrium \(E_3\), we now employ the trace-determinant plane approach, although this classical stability criterion appears in various texts under different presentations, we adopt here the formulation in (Saber et al., 2007) for its particular convenience in our analysis of the coexistence equilibrium \(E_3\). The trace-determinant criterion proves especially suitable for establishing explicit stability conditions in terms of the trace and determinant of the Jacobian matrix.

To support the previous criteria, we now present a classical result for matrices \(A=[a_{ij}]_{2 \times 2}\), which establishes conditions on the trace \(\mathrm{tr} \mathbf{A}=a_{11}+a_{22}\) and determinant \(\det\mathbf{A}=a_{11}a_{22}-a_{12}a_{21}\) to guarantee that \(\varrho(\mathbf{A})<1\). This result, together with Theorem 1, will be instrumental for analyzing the stability of the coexistence equilibrium \(E_3\). For a complete proof and a detailed development of this methodology, see Saber et al. (2007), Chapter 4, Theorem 4.4, pages 200–201.

Theorem 2. (Stability in the Trace-Determinant Plane) Let \(\mathbf{A} = (a_{ij})\) be a \(2\times 2\) matrix. Then the spectral radius \(\varrho(\mathbf{A}) < 1\) if and only if

\[ |\mathrm{tr} \mathbf{A}| - 1 < \det \mathbf{A} < 1. \tag{15} \]

We apply this result to our problem with \(\mathbf{A}=\mathbf{J}^*_3\) to identify relationships between the parameters \(r,\,b,\,c,\,d\) and \(K\) that satisfy inequality (16). For this purpose, we observe that:

\[ | \mathrm{tr} \mathbf{J}^*_3| - 1 < \det \mathbf{J}^*_3 < 1 \tag{16} \]

Inequality (16) expands to the following explicit form:

\[ \left|2- \frac{rc}{dK} \right| \ < 2-\frac{rc}{dK}(1+c-dK)<2. \tag{17} \]
Furthermore, the modulus function in inequality (17) can be expanded as

\[ \left|2- \frac{rc}{dK} \right| =\begin{cases} 2- \dfrac{rc}{dK} &, \text{if} \dfrac{rc}{dK}<2; \\[0.15cm] 0 &,\text{if} \dfrac{rc}{dK}=2 \\[0.15cm] \dfrac{rc}{dK}-2 &, \text{if} \dfrac{rc}{dK}>2. \end{cases} \tag{18} \]

To resolve inequality (17), we use the piecewise expansion of the modulus function (18) to break the problem into cases. This enables a straightforward application of Theorem 2 to each parameter regime.

Applying Theorems 1 and 2, we conclude that the equilibrium \(E_3\) is asymptotically stable when \(K~\in~\left( \dfrac{c}{d}, \dfrac{c+1}{d} \right)\).

The condition \(K \in \left( \dfrac{c}{d}, \dfrac{c+1}{d} \right)\) is equivalent to \(dK-1 < c < dK\), and in biological terms, this means the mortality rate \(c\) cannot be excessively low (which would lead to instability) nor excessively high (which would prevent predator persistence). This result highlights the sensitivity of the system’s stability to small variations in mortality rate, particularly relative to maximum predation (\(dK\)).

To illustrate the local asymptotic stability of \(E_3\), we return to the parameter set used in Figure 1 (\(r = 1.4\), \(b = 0.01\), \(c = 0.05\), \(d = 0.001\)), which satisfies the asymptotic stability conditions established for \(E_3\). Figure 4 shows numerical simulations of system of equation (1), approach to this equilibrium from three different initial conditions: Trajectory 1 (red) starts near the trivial equilibrium \(E_1\); Trajectory 2 (blue) starts near the predator-free equilibrium \(E_2\); and Trajectory 3 (purple) starts near the coexistence equilibrium \(E_3\).

Figure 4 illustrates how the dynamics approach equilibrium, using the parameter values \(r = 1.4\), \(b = 0.01\), \(c = 0.05\), and \(d = 0.001\), showing 250 and 1000 iterations, respectively.

(a)
image
(b)
image
Figure 4 – Asymptotically stable at \(E_3\), with the same parameters as in Figure 1: (a) 250 iterations and (b) 1000 iterations.

Remark: For \(\dfrac{rc}{dK} > 2\), algebraic manipulations do not preserve the equivalence of inequality (16), yielding only partial implications. Investigating this interval requires different methods, whose analysis will be considered in future work.

The local stability analysis revealed explicit conditions for each system. Table 1 summarizes these conditions as functions of biological parameters, highlighting the associated equilibrium regimes.

Table 1 - Classification of the system’s equilibria according to parameter values.
Equilibrium Parameters Classification
\(E_1=(0,0)\) \(c > 2\) Unstable source
\(0 < c < 2\) Unstable saddle point
\(E_2=(K,0)\) \(0 < r < 2\) \(dK < c < 2 + dK\) Asymptotically Stable
\(c < dK \;\; \text{or} \;\; c > 2 + dK\) Unstable saddle point
\(r > 2\) \(dK < c < 2 + dK\)
\(c < dK \;\; \text{or} \;\; c > 2 + dK\) Unstable source
\(E_3=\left(\dfrac{c}{d},\,\dfrac{r}{Kb}\left(K-\dfrac{c}{d}\right)\right)\) \(0 < \dfrac{rc}{dK} < 2\) \(K \in \left(\dfrac{c}{d},\, \dfrac{c+1}{d}\right)\) Asymptotically Stable
\(\dfrac{rc}{dK} = 2\)

Conclusions

This work has presented a local stability analysis for a system of the Lotka-Volterra type formulated directly in discrete time. The main contribution lies in obtaining explicit and interpretable stability conditions for each of the three biologically relevant equilibria, based on linearization and spectral analysis of the Jacobian matrix.

Our stability analysis reveals three distinct dynamical regimes governed by the system’s biological parameters. This classification enables comparative identification of critical intervals for coexistence, exclusion, or population collapse, serving as a practical reference for both theoretical interpretations and future applications in simulations or management. The analysis establishes that:

Although the parameter \(b\), representing the predation rate, does not directly affect local stability criteria, it determines the predator density at the coexistence equilibrium. Thus, its variation strongly influences the configuration of the coexistence equilibrium or ecological viability.

Beyond the theoretical results, this work provides a pedagogical synthesis and application of classical stability methods for discrete systems. This unified framework, which includes linearization, spectral analysis, and the trace-determinant criterion, is computationally efficient as it is developed entirely in discrete time. The numerical simulations, while challenging to calibrate, were designed not only to visualize the dynamic regimes but also to show the practical utility of this methodological integration in an ecological context. This approach, which explicitly bridges analytical results with interpretable dynamics, establishes discrete modeling as a powerful tool for both education and the analysis of ecosystems with sparse monitoring data, providing a foundation for future studies on management strategies.

Finally, although the model depicts natural interactions, its transposition to the human realm reveals a reality difficult to ignore: the same species capable of rapid regeneration and longevity exerts, through its actions, one of the most intense predatory pressures on the very system that sustains it. Recognizing this force is not an invitation to resignation, but to active awareness. Models like this not only outline the fragile contours of stability but also remind us that it endures only when understood and respected.

Acknowledgments

The authors would like to thank Victoria Gattass for her careful review on the English version and translation of this work.

Author Contributions

F. F. C. Puma: conceptualization, project managements, formal analysis, data curation, investigation, methodology, supervision, validation, visualization and writing – preparation of the original draft , revision and editing. M. M. Mager: formal analysis, visualization, supervision, validation and writing – preparation of the original draft , revision and editing.

Conflicts of Interest

The authors declare no conflict of interest.

References

Almaraz, Pablo, Kalita, Piotr, Langa, José A. & Soler-Toscano, Fernando (2024). Structural stability of invasion graphs for Lotka-Volterra systems. Journal of Mathematical Biology 88(64), 1-25. https://doi.org/10.1007/s00285-024-02087-8

Braverman, Elena & Kinzebulatov, Damir (2006). On linear perturbations of the Ricker model. Mathematical Biosciences 202(2), 323-339. https://doi.org/10.1016/j.mbs.2006.04.008

Din, Q. (2013). Dynamics of a discrete Lotka-Volterra model. Advances in Difference Equations 2013(1), 95. https://doi.org/10.1186/1687-1847-2013-95

Diniz, Geraldo Lúcio (2011). Equações de diferenças e sistemas: com aplicações biológicas. SBMAC 54

Fisher, M. E., Goh, B. S. & Vincent, T. L. (1979). Some Stability Conditions for Discrete-Time Single Species Models. Bulletin of Mathematical Biology 41, 861-875. https://doi.org/10.1007/BF02462383

Kekulthotuwage Don, S., Burrage, K., Helmstedt, K. J. & Burrage, P. M. (2023). Stability Switching in Lotka-Volterra and Ricker-Type Predator-Prey Systems with Arbitrary Step Size. Axioms 12(4), 390. https://doi.org/10.3390/axioms12040390

Khaliq, Abdul, Ibrahim, Tarek F., Alotaibi, Abeer M., Shoaib, Muhammad & El-Moneam, Mohammed Abd (2022). Dynamical Analysis of Discrete-Time Two-Predators One-Prey Lotka–Volterra Model. Mathematics 10(21), 4015. https://doi.org/10.3390/math10214015

Krabs, Werner & Pickl, Stefan (2010). Dynamical Systems: Stability, Controllability and Chaotic Behavior. Springer

Luiz, Kariston Stevan, Organista, Juniormar, Cirilo, Eliandro Rodrigues, Romeiro, Neyva Maria Lopes & Natti, Paulo Laerte (2022). Convergência numérica de um sistema Telegraph Predator-Prey. Semina: Ciências Exatas e Tecnológicas 43(1Esp), 51-66. https://doi.org/10.5433/1679-0375.2022v43n1Espp51

Lv, Luyao & Li, Xianyi (2024). Stability and Bifurcation Analysis in a Discrete Predator–Prey System of Leslie Type with Radio-Dependent Simplified Holling Type IV Functional Response. Mathematics 12(12), 1803. https://doi.org/10.3390/math12121803

Malthus, Thomas Robert (1826). An Essay on the Principle of Population. John Murray https://commons.wikimedia.org/w/index.php?title=File:Malthus_(IA_darwin-online_1826_Malthus_A545.2).pdf&page=4

May, Robert M. (1976). Simple mathematical models with very complicated dynamics. Nature 261, 459-464. https://doi.org/10.1038/261459a0

Saber, Elaydi (2005). An Introduction to Difference Equations. Springer

Saber, Elaydi (2007). Discrete Chaos: With Applications in Science and Engineering. CRC Press https://www.taylorfrancis.com/books/9781420011043

Seno, H. (2008). A paradox in discrete single species population dynamics with harvesting/thinning. Mathematical Biosciences 214(1), 63-69. https://doi.org/10.1016/j.mbs.2008.06.004

Verhulst, P. F. (1838). Notice sur la loi que la population poursuit dans son accroissement. In A. Quetelet, Correspondance Mathématique et Physique (Vol. 10, pp. 113-121). https://books.google.com.br/books?hl=fr&id=8GsEAAAAYAAJ&jtp=113&redir_esc=y#v=onepage&q&f=false