Dynamical analysis of a predator-prey system with prey vigilance and hunting cooperation in predators

1.   Introduction

Interactions between individuals are a crucial aspect of life history traits for many species ^[1]. Predator-prey systems have been used to study various ecological population interactions ^[2,3,4,5]. The effect of the presence of predators does not only directly impact prey through killing and consumption, but also induces non-lethal effects such as fear. This fear can drive the prey to use tactics to secure their lives ^[6]. One such tactic is vigilance. Vigilance is the act of keeping a careful watch in an environment for any possible harm. Vigilance involves a concerted effort by the prey population to actively avoid predation by means of avoiding competition, approaching group cohesion, and detecting and avoiding predators, ^[7] thereby reducing predation risks ^[8]. Vigilance behavior has been observed in several species, including elk ^[9], seals ^[10], and mother cheetahs, who protect their young ^[11]. When prey live in a group, vigilance behavior is beneficial to them because, while their proximity to the group increases their conspicuity, their grouping helps reduce predation risks ^[12]. When prey adopt anti-predator behavior such as vigilance, their rate of food intake is reduced due to less time spent foraging ^[13,14]. Interestingly, this also affects the food intake rate of predators as prey decrease their predation risks ^[15]. In intraguild predation, an intermediate predator can exhibit vigilance behavior to avoid predation by a top predator. This behavior can lead to a reduced efficiency in hunting for shared prey species. For example, field experiments conducted by Durant ^[16] revealed that when cheetahs listened to the playbacks of lion vocalizations, they were less likely to make a kill or hunt after hearing the playback. They move just as far from the area of the playback. Mathematical models have been used to gain insights into the effects of prey vigilance. Kimbrell et al. ^[17] studied the influence of vigilance on intraguild populations. Their results showed that when top predators kill intermediate predators without eating them, it can increase the level of vigilance by the intermediate predator or influence the vigilance behavior of the shared prey, which may aid in the stability of the ecological community. Hossain et al. ^[18] also studied vigilance dynamics in a three-species food chain model. Their model produced rich dynamics, including a Hopf bifurcation, shrimp-shaped periodic structures, and multiple coexisting attractors. Their results also suggested that too much vigilance can lead to species extinction. This is because the prey will starve and/or reproduce less, and hence reduce its lifetime reproductive fitness ^[19].

As said earlier, vigilance behavior in prey can impact the food intake rate of predators as a result of a decrease in prey vulnerability. Therefore, many predators enhance their predation efficiency when they engage in cooperative hunting. Hunting cooperation is the combined effort of several individuals to capture and share prey ^[20]. Many predators, especially carnivores, work together to hunt and to forage ^[21]. Carnivores such as lions ^[22], African wild dogs ^[23], chimpanzees, ^[24] and wolves ^[25] have been documented to engage in cooperative hunting. Hunting cooperation comes with its benefits. Included are increased hunting success rates with the number of adults, decreased chasing distance ^[23,26], more effective utilization of food resources ^[27], high likelihood of capturing large prey ^[28], less time finding food ^[29], and also protection of food (carcasses) from being stolen by other predators ^[30]. There are several continuous time models which have studied the impacts of hunting cooperation among predators. Alves and Hilker ^[21] found that hunting cooperation can destabilize the system and lead to a collapse of the predator population. Berec ^[31] studied hunting cooperation effects in relation to population oscillations and concluded that the stability of coexistence states could change due to cooperation. Pal et al. ^[32] studied a modified Leslie-Gower predator-prey model with hunting cooperation among predators and fear effect in prey. Their findings revealed that hunting cooperation can induce both subcritical and supercritical Hopf bifurcations. Spatially explicit models have also been used to explore hunting cooperation effects. A variety of spatio-temporal dynamics such as spots, stripe patterns, and mixed patterns (spots and stripes) were observed for different intensities of the rate of hunting cooperation among predators ^[33]. The spatially explicit model in ^[34] cannot produce Turing patterns when hunting cooperation is absent, whereas the model with hunting cooperation can. Discrete-time models have also been used to study cooperative hunting effects in predator-prey relationships. Pal et al. ^[35] showed that hunting cooperation is able to stabilize a chaotic discrete-time system and can induce strong demographic Allee effects.

Many researchers have studied the impact of hunting cooperation ^{[21,31,32,35,36]} and vigilance ^[17,18,37] in predator-prey systems separately. However, their combined effects in predator-prey dynamics is yet to be studied. The aim of this paper is to explore the dynamics when both vigilance behavior in prey and hunting cooperation in predators are present. We organize our paper as follows: We present our proposed modeling framework with its underlying ecological assumptions in Section 2. Preliminary results on positivity and boundedness of solutions are presented in Section 3. Section 4 is dedicated to finding feasible equilibria and performing stability analysis on our proposed model. We derive local codimension $1$ bifurcation results in Section 5. In Section 6, we provide numerical experiments to validate our theoretical findings. We study the dynamics of our proposed model when predators do not hunt cooperatively in Section 7. We conclude the paper with a discussion of our results and possible future work in Section 8.

2.   Model formulation

Here, we consider an unstructured prey and predator population. We let $x$ and $y$ denote the prey and predator population respectively at any time instant $t$. We take into account the following assumptions in our model formulation:

$({\rm{i}})$ The prey population grows logistically in the absence of predators and vigilance behavior.

$({\rm{ii}})$ We let the parameter $v$ denote the level of prey vigilance where $v\in [0, 1]$. Also, the lethality of predation is $\frac{1}{l}$ when vigilance is absent.

$({\rm{iii}})$ We suppose that predators cooperate when they hunt the prey.

$({\rm{iv}})$ We use the Holling type Ⅱ functional response to describe the relationship between the predator and its prey.

$({\rm{v}})$ We assume natural death rates $\delta$ for the prey and $\delta_1$ for the predator.

The following nonlinear system of ordinary differential equations satisfies our assumptions:

with positive initial conditions $x(0) = x_0$ and $y(0) = y_0$. We assume all parameters used are positive, and their descriptions are provided in Table 1.

3.   Preliminary results

This section provides basic results on the positivity and boundedness of solutions to system (2.1) for biological meaningfulness.

3.1.   Positivity of solutions

We recap the following result which guarantees the positivity of solutions from ^[38,39].

Lemma 3.1. Consider the following system of ODEs:

Non-negativity of solutions is preserved with time, that is

if and only if

and thus we have

3.2.   Boundedness of solutions

The boundedness property of solutions to system (2.1) ensures that populations do not grow unboundedly due to scarce food resources and limited habitat space.

Theorem 3.2. Solutions to system (2.1) are bounded when they initiate from $\mathbb{R}_2^+$.

Proof. By considering Lemma 3.1,

Using the comparison principle and simple calculations,

as $t\rightarrow \infty$.

Let $W(t) = x(t)+\frac{1}{\gamma}y(t)$. Then, for large $t$ we have

Therefore,

As $t\rightarrow \infty$, ${\rm{lim}} \; {\rm{sup}}$ $W(t) \leq \frac{r(1-v)^2K}{\epsilon}$ and hence all solutions starting from $\mathbb{R}_2^+$ are bounded. □

4.   Equilibria and stability analysis

4.1.   Equilibria

To obtain the equilibria for system (2.1), we solve $X(x, y) = 0$ and $Y(x, y) = 0$ simultaneously. The system possesses the following non-negative equilibria:

$({\rm{a)}}$ $E_{0} = (0, 0)$,

$({\rm{b}})$ $E_{1} = (x_{1}^{*}, 0)$ where $x_1^* = K \left(1-v-\frac{\delta}{r} \right)$. We note that $E_{1}$ is feasible when $1-\frac{\delta}{r} > v$.

$({\rm{c}})$ $E_2 = (x^*, y^*)$ where $y^* = \frac{1}{c}\left(\frac{(1+x^*)\delta_1}{Ax^*}-q \right)$, $A = \frac{\gamma}{l+\mu v}$, and $x^*$ is a positive real root of the following third-order equation

This is obtained by substituting $q+cy^* = \frac{\delta_1\left(1+x^* \right)}{Ax^*}$ into the right-hand side of $\frac{dx}{dt}$. $E_{2}$ is feasible when $\frac{(1+x^*)\delta_1}{Ax^*} > q$. Thus, $0 < x^* < \frac{\delta_1}{Aq-\delta_1}$. This implies that $Aq-\delta_1 > 0$.

We consider three cases in determining the number of positive real roots to Eq (4.1) using Descartes' rule of signs. These cases are when

$({\rm{I}})$ $r(1-v)-\delta > 0$

$({\rm{II}})$ $r(1-v)-\delta < 0$

$({\rm{III}})$ $r(1-v)-\delta = 0$.

In each of the cases above, the number of possible non-negative real roots for Eq (4.1) is $2$. Therefore, when $x^*$ is obtained from Eq (4.1) and substituted into $y^*$, we may either have two feasible interior equilibria or one feasible interior equilibrium point.

4.2.   Global stability analysis

Define $\chi = \frac{\gamma r(1-v)^2K}{\epsilon}$. Using results from Theorem 3.2, we have $y(t) \leq \chi.$ We state the following theorem:

Theorem 4.1. The predator-free state $E_{1}$ is globally stable if $1 < \frac{(l+\mu v)\delta_1}{\gamma x_1^{*}(q+c \chi)}$.

Proof. Suppose that $1 < \frac{(l+\mu v)\delta_1}{\gamma x_1^{*}(q+c \chi)}$ and consider the Lyapunov function $V(t) = A_1\left[x-x_1^{*}-x_1^{*} \text{ln}\left(\frac{x}{x_1^{*}}\right)+\frac{y}{\gamma}\right]$ where $A_1$ is a positive constant to be chosen. Clearly, $V = 0$ at $(x, y) = (x_1^{*}, 0)$. Also, $V > 0$ when $(x, y) \neq (x_1^{*}, 0).$ Now, evaluating the derivative of $V$ with respect to $t$ yields

We substitute $r(1-v)-\delta = \frac{rx_1^{*}}{K}$ and $x = x-x_1^{*}+x_1^{*}$. Then,

Since $y(t) \leq \chi$,

We choose $A_1 = \frac{l+\mu v}{x_1^{*}(q+c\chi)}$. Then,

Since our Lyapunov function satisfies the asymptotic stability theorem ^[40,41], then by our theorem, $E_1$ is globally stable. This completes the proof. □

Theorem 4.2. The extinction state $E_{0}$ is globally stable if $v > 1-\frac{\delta}{r}$.

Proof. We provide the proof in the Appendix. □

Remark 1. The conditions stated in Theorems 4.1 and 4.2 are sufficient conditions.

4.3.   Local stability analysis

We compute the Jacobian of system (2.1) to aid in the local stability analysis of the feasible equilibria. The Jacobian is given by

We state the following:

Theorem 4.3. The predator-free state $E_{1}$ is locally stable if $v < 1-\frac{\delta}{r}$ and $\frac{\gamma K q (r (1-v)-\delta)}{(l+\mu v) (K(r(1-v)-\delta)+r)} < \delta_1$.

Proof. Suppose $v < 1-\frac{\delta}{r}$ and $\frac{\gamma K q (r (1-v)-\delta)}{(l+\mu v) (K(r(1-v)-\delta)+r)} < \delta_1.$ Evaluating $J^*$ at $E_1$ gives

Since eigenvalues $\lambda_1 = r(v-1)+\delta < 0$ and $\lambda_2 = \frac{\gamma K q (r (1-v)-\delta)}{(l+\mu v) (K(r(1-v)-\delta)+r)}-\delta_1 < 0,$ $E_1$ is locally stable. □

In ecosystems, it is very common to see coexistence of species. Therefore it is important to study the dynamics pertaining to the stability of the coexistence equilibrium $E_2$ using standard linear stability analysis. The characteristic equation for $J^*_{E_2}$ is given by

where

and

Here, $tr(J^*_{E_2})$ and $det(J^*_{E_2})$ represent the trace and determinant of $J^*$ evaluated at $E_2$. The stability of $E_2$ depends on the signs of $tr(J^*_{E_2})$ and $det(J^*_{E_2})$. Through the Routh-Hurwitz criteria, we state the following theorem in connection to the local stability of $E_2$.

Theorem 4.4. For $0 < x^* < \frac{\delta_1}{Aq-\delta_1}$ with $Aq-\delta_1 > 0$, the coexistence state $E_{2}$ is locally stable if $tr(J^*_{E_2}) < 0$ and $det(J^*_{E_2}) > 0$.

We present numerical results for the findings in Theorem 4.4. We consider the parameters in Figure 1(d). The coexistence equilibrium is $E_2(1.06618, 0.515449)$. Simple calculations show that $1.06618 = x^* < \frac{\delta_1}{Aq-\delta_1} = 2.609$. Evaluating $J^*$ at $E_2$ yields

From Eq (4.4), $tr(J^*_{E_2}) = -0.000857263 < 0$ and $det(J^*_{E_2}) = 0.000578922 > 0$. Therefore, $E_2(1.06618, 0.515449)$ is locally stable.

Remark 2. If $tr(J^*_{E_2})\leq0$ or $> 0$ and $det(J^*_{E_2}) < 0$, then $E_2$ is a saddle.

5.   Bifurcation analysis

Bifurcation analysis plays an important role in providing insights into the qualitative behavior of a system when parameters are varied continuously. We focus on the effects of prey vigilance levels and the rate of hunting cooperation on the dynamics of system (2.1). Therefore, we explore local codimension $1$ bifurcations and find the occurrence of Hopf, saddle-node, and transcritical bifurcations.

Theorem 5.1. Suppose that $E_2$ exists and consider the Jacobian of system (2.1). Then, system (2.1) experiences a Hopf bifurcation at $E_2$ with respect to the bifurcation parameter $c$ if the following hold:

$(i)$ $\text{tr}(J^*_{E_2}) = 0,$

$(ii)$ $\text{det}(J^*_{E_2}) > 0,$

$(iii)$ $\frac{d}{dc} \left(\text{tr}(J^*_{E_2})\right) \neq 0$.

Proof. Simple calculations show that $\text{tr}(J^*_{E_2}) = 0$ when

where $p_1 = \delta _1 K (x^*+1) (l+\mu v)+2 r x^* (x^*+1) (l+\mu v)$. The Jacobian evaluated at $E_2$ with $c = c^*$ is

Now, $det(J^*_{E_2}) = \frac{p_2}{K (x^*+1) (l+\mu v)}$ where

with $p_3 = \delta _1 K (x^*+1) (l+\mu v)-4 r x^* (x^*+1) (l+\mu v)$ and $p_4 = \gamma q x^* (K (\delta +r (v-1))+2 r x^*)$, respectively. We let $p_2 > 0$. We proceed to validate the transversality condition of the Hopf bifurcation theorem ^[42,43] by letting $\omega_1 = -\frac{2r}{K}-\frac{\gamma}{\left(l+\mu v \right)\left(x^*+1 \right)^2}-\frac{\delta_1 q \left(2x^*+1 \right)}{\gamma c^* \left(x^* \left(x^*+1 \right) \right)^2 } + \frac{2\delta_1^2 \left(l+ \mu v \right)}{\gamma^2 c^* x^{*3}}$ and $\omega_2 = \frac{\delta_1^2 \left(l+ \mu v \right)}{\left(\gamma c^* x^* \right)^2}-\frac{\delta_1 q}{\gamma c^{*2} x^* \left(x^*+1\right)}$, and ensuring that

Hence, system (2.1) undergoes a Hopf bifurcation around $E_2$ with respect to the bifurcation parameter $c$. □

For system (2.1) to experience a saddle-node bifurcation around $E_2$ for the parameter $c$, it is required that $det(J^*) = 0$. To obtain this, we give an implicit expression for $c$ as $c = c^* = \frac{\delta _1^2 K (l+\mu v)}{\gamma ^2 x^{*^2} (K (\delta +r (v-1))+2 r x^*)}$ since $E_2(x^*, y^*)$ depends on $c$. We state the following theorem accordingly:

Theorem 5.2. Suppose that $E_2$ exists. Then, system (2.1) experiences a saddle-node bifurcation around the coexistence equilibrium $E_2$ at $c = c^*$ when $tr(J^*) < 0$ and $det(J^*) = 0$ are satisfied by system parameters.

Proof. Let $c = c^* = \frac{\delta _1^2 K (l+\mu v)}{\gamma ^2 x^{*^2} (K (\delta +r (v-1))+2 r x^*)}$ and $\frac{(K (\delta +r (v-1))+2 r x^*) \left(\gamma q x^*-2 \delta _1 (x^*+1) (l+\mu v)\right)}{\delta _1 K (x^*+1) (l+\mu v)} < \frac{\gamma q x^*}{(x^*+1) (l+\mu v)}-\delta_1$. We use Sotomayor's theorem ^[43] to show that system (2.1) experiences a saddle-node bifurcation at $c = c^*$. At $c = c^*$, we can have $det(J^*) = 0$ and $tr(J^*) < 0$ when $\frac{(K (\delta +r (v-1))+2 r x^*) \left(\gamma q x^*-2 \delta _1 (x^*+1) (l+\mu v)\right)}{\delta _1 K (x^*+1) (l+\mu v)} < \frac{\gamma q x^*}{(x^*+1) (l+\mu v)}-\delta_1$. This shows that $J^*$ admits a zero eigenvalue. Define $G = (g_1, g_2)^T$ and $H = (h_1, h_2)^T$ to be the nonzero eigenvectors of $J^*$ and $J^{*T}$ corresponding to the zero eigenvalue, respectively. Then,

$G = \left(-\frac{\delta _1 K}{\gamma (K (\delta +r (v-1))+2 r x^*)}, 1 \right)^T$ and $H = \left(\frac{\gamma \left(\gamma q x^*-\delta _1 (x^*+1) (l+\mu v)\right)}{\gamma q x^*-2 \delta _1 (x^*+1) (l+\mu v)}, 1 \right)^T$.

Furthermore, let $Z = (z_1, z_2)^T$ where

Now,

provided $1-\gamma q x^* +\delta_1 (x^*+1)(l+\mu v)\neq 0$ and $\gamma q x^* - 2 \delta_1 (x^*+1)(l+\mu v) \neq 0$. Furthermore,

Therefore, by Sotomayor's theorem, system (2.1) experiences a saddle-node bifurcation at $c = c^*$ around $E_2$, which concludes the proof. □

Similarly we can give an implicit expression for $v = v^* = \frac{\frac{\gamma q x^*}{\delta _1 (x+1)}-l}{\mu }$ or $v = v^* = \frac{c \gamma ^2 x^{*2} (K (r-\delta)-2 r x^*)+\delta _1^2 K l}{K \left(c \gamma ^2 r x^{*2}-\delta _1^2 \mu \right)}$ to ensure that $det(J^*) = 0$ since $E_2(x^*, y^*)$ depends on $v$. The conditions under which $tr(J^*) < 0$ can easily be found. Next, we state the following theorem:

Theorem 5.3. Suppose that $E_2$ exists. Then, system (2.1) experiences a saddle-node bifurcation around the coexistence equilibrium $E_2$ at $v = v^*$ when $tr(J^*) < 0$ and $det(J^*) = 0$ are satisfied by system parameters.

Proof. The proof is similar to Theorem 5.2 and is therefore omitted. □

Theorem 5.4. Suppose that $E_1$ exists. Then, system (2.1) experiences a transcritical bifurcation around the predator-free state $E_1$ when the level of vigilance is $v = v^* = 1-\frac{\delta}{r}$.

Proof. Suppose $v = v^* = 1-\frac{\delta}{r}$. An evaluation of the Jacobian matrix for system (2.1) at $E_1$ with $v^*$ is

The eigenvalues of the Jacobian matrix in Eq (5.5) are $\lambda_1 = 0$ and $\lambda_2 = -\delta_1$. Next, we represent the eigenvectors corresponding to the zero eigenvalue of the matrices $J_{E_1}^{*}$ and $J_{E_1}^{*T}$ respectively by $L = (l_1, l_2)^T$ and $M = (m_1, m_2)^T$. Simple calculations show that $L = \left(1, 0 \right)^T$ and $M = \left(1, 0 \right)^T$. Now, let $Z = (z_1, z_2)^T$ as defined in Eq (5.4). We proceed to validate the transversality conditions using Sotomayor's theorem ^[43]. Now,

Also,

and

Therefore, by the Sotomayor's theorem, system (2.1) experiences a transcritical bifurcation at some $v = v^*$ around $E_1$. □

6.   Numerical experiments

In this section, we provide numerical simulations which support our theoretical results using the Python programming language, Wolfram Mathematica 13.0, MATLAB version R2019a, and MATCONT ^[44]. We show the existence of biologically feasible equilibria when the nullclines for the prey and predator populations intersect for certain parameter choices for system (2.1). See Figure 1. We provide experiments in Figures 2 and 3 respectively to validate sufficient conditions for global stability results for the predator-free and extinction states. Applications of these results in this section are discussed in Section 8. We show local codimension one bifurcations for the level of vigilance parameter $v$ and hunting cooperation parameter $c$.

Figure 4(a) shows the existence of Hopf, saddle-node, and transcritical bifurcations when $v$ is varied for certain parameter choices. When the vigilance level $v$ is increased, the coexistence state gains stability at the critical threshold $v^* = 0.386322$ around $E_2 = (14.448794, 4.769456)$, and the disappearance of oscillatory dynamics is observed. We used MATCONT to compute the Lyapunov coefficient. This value is given by $\sigma_1 = -8.447014e^{-4}$, and thus the bifurcation is supercritical. A slight increase in $v$ causes the system to experience a saddle-node bifurcation at $v^* = 0.405605$. At this level, two coexistence equilibria (a saddle and a node) collide and disappear. This bifurcation occurs around $E_2 = (23.146242, 2.795123)$. A transcritical bifurcation occurs when the stable coexistence equilibrium $E_2$ collides and interchanges its stability property with the unstable predator free state $E_1$. Hence $E_2$ becomes unstable and $E_1$ gains stability. Here, this bifurcation is observed at vigilance level $v^* = 0.399556$ around $E_1 = (29.022207, 0)$. Similar bifurcations are seen in Figure 4(b).

System (2.1) experiences saddle-node and Hopf bifurcations for the hunting cooperation parameter $c.$ See Figure 4(b). The saddle-node bifurcation occurs at $c^* = 0.032175$ around $E_2 = (22.139420, 2.488556)$. The Hopf bifurcation is seen to occur at $c^* = 0.033090$ around $E_2 = (19.226665, 2.534502)$, and the calculated Lyapunov coefficient is $\sigma_2 = -1.758377e^{-3}$. The Hopf bifurcation is supercritical.

7.   Case where there is no hunting cooperation ($c=0$)

We study a special case where predators do not cooperate when hunting. Thus, system (2.1) reduces to

7.1.   Equilibria

The feasible equilibria for system (7.1) are

${\rm{(a}})$ $E_{0}' = (0, 0)$,

$({\rm{b}})$ $E_{1}' = \left(K \left(1-v-\frac{\delta}{r}\right), 0\right)$. $E_{1}'$ is feasible when $1-\frac{\delta}{r} > v$.

$({\rm{c}})$ $E_{2}' = \left(\frac{\delta _1 (l+\mu v)}{\gamma q-\delta _1 (l+\mu v)}, \frac{\gamma (l+\mu v) \left(\delta _1 (l+\mu v) (\delta K+r (K (v-1)-1))-\gamma K q (\delta +r (v-1))\right)}{K \left(\gamma q-\delta _1 (l+\mu v)\right)^2}\right)$.

$E_{2}'$ exists when $\gamma q > \delta _1 (l+\mu v)$ and $\delta _1 (l+\mu v) (\delta K+r (K (v-1)-1)) > \gamma K q (\delta +r (v-1))$.

We state the following theorem pertaining to the global stability of the unique coexistence equilibrium $E_{2}'$. From $E_{2}'$, we let $y^* = \frac{\gamma (l+\mu v) \left(\delta _1 (l+\mu v) (\delta K+r (K (v-1)-1))-\gamma K q (\delta +r (v-1))\right)}{K \left(\gamma q-\delta _1 (l+\mu v)\right)^2}$.

Theorem 7.1. The coexistence state $E_{2}'$ is globally stable if $y^* < \frac{r(l+ \mu v)}{Kq}$.

Proof. Suppose that $y^* < \frac{r(l+ \mu v)}{Kq},$ and consider the Lyapunov function $V(t) = A\left[x-x^*-x^* \text{ln} \left(\frac{x}{x^*}\right)\right] + \frac{B}{\gamma} \left[y-y^*-y^* \text{ln} \left(\frac{y}{y^*}\right)\right]$ where $A, B$ are positive constants to be determined. Clearly, $V = 0$ at $(x, y) = (x^{*}, y^{*})$. Also, $V > 0$ when $(x, y) \neq (x^{*}, y^{*}).$ Now, evaluating the derivative of $V$ with respect to $t$ yields

Using the results $r(1-v) - \delta = \frac{rx^{*}}{K} + \frac{qy^{*}}{(1+x)(l+ \mu v)}$ and $\frac{\gamma q x^{*}}{(1+x^*)(l+ \mu v)} = \delta_1$, we have

Here, we choose $A = \frac{1}{1+x^{*}}$ and $B = 1$. Thus,

Here also, the Lyapunov function satisfies the asymptotic stability theorem ^[40,41], and by our theorem, $E_2'$ is globally stable. This completes the proof. □

We omit the local stability analysis of all the equilibria for system (7.1) as well as global stability results for the extinction state and the predator-free state for brevity.

8.   Discussion and conclusions

In ecosystems, many species exhibit anti-predator behaviors such as vigilance to mitigate threats and predation risks. When prey populations are vigilant, it makes predators spend more time and energy in capturing them. In order to capture prey efficiently, predators cooperate during their hunt. In this work, we explore the impacts of prey vigilance and hunting cooperation in a predator-prey system. Our results show that, for certain parameter choices, an increase in the level of vigilance can stabilize the system via a Hopf bifurcation for a fixed hunting cooperation rate. However, for a fixed prey vigilance level, an increase in the rate of hunting cooperation can cause the system to destabilize. See Figures 1 and 4. We also observed from Figure 4(a) that too much vigilance by prey can have a negative effect, causing the extinction of the population due to a continuous decrease in population density. This is because they trade-off between foraging and staying alert. This will cause starvation and reduce lifetime reproductive fitness ^[19]. For example, the Nubian Ibex is now known to be vulnerable to extinction ^[45] and is very vigilant when obstructed between their safety region and food patch ^[46]. Other bifurcation results such as saddle-node and transcritical bifurcations were observed. We obtained sufficient conditions for the global stability of the predator-free state and the extinction state with rigorous proofs. Refer to Theorems 4.2 and 4.1. The transcritical and global stability results will provide ecosystem managers with information on how best to provide structures and develop strategies in conserving endangered species and thus promote their persistence. Furthermore, our results show that hunting cooperation can change the stability of a coexistence state. Refer to Figure 1(d) and (e). This supports the results obtained by Berec in ^[31]. Our proposed system exhibited rich dynamical behavior including bi-stability between a stable limit cycle and the predator-free equilibrium. See Figure 1(b). Prey and predator populations will go between oscillatory populations and stable levels. In this case, prey vigilance levels and cooperative hunting play a role in maintaining ecosystem stability. Therefore, the sensitivity to initial conditions will play a significant role in determining whether the two species will continue to coexist or the predator population will die out. We also found that hunting cooperation when intensified can cause a decrease in the population densities of both prey and predators when vigilance levels are fixed. See Figure 1(d)–(f). When prey are at low densities and predators hunt cooperatively, it can lead to a reduction in the growth rate of the predator population and hence induce an Allee effect. It will be interesting to study an extension of our temporal model by incorporating Allee effects into both prey and predator populations. A study of such a mechanism will be useful in biocontrol and species conservation programs. We will extend our temporal model to include spatial effects to explore the possible occurrence of Turing patterns which provide insights on how hunting cooperation and prey vigilance contribute to the patchy spread of species in space.

Use of AI tools declaration

The authors declare they have not used Artificial Intelligence (AI) tools in the creation of this article.

Acknowledgments

EMT, CO, MC and NK would like to acknowledge valuable support from the National Science Foundation via grant number 1851948.

Conflict of interest

The authors declare there is no conflict of interest.

Appendix

We provide the proof for Theorem 4.2.

Proof. Suppose that $v > 1-\frac{\delta}{r},$ and consider the Lyapunov function $V(t) = \gamma x(t)+y(t)$. Clearly, $V = 0$ at $(x, y) = (0, 0)$. Also, $V > 0$ when $(x, y) \neq (0, 0).$ Now, evaluating the derivative of $V$ with respect to $t$ yields

□