#### Abstract

In this paper, we present and analyze predator-prey system where prey population is linearly harvested and affected by fear and the prey population has grown logistically in the absence of predators. The predator population follows hunting cooperation, and it predates the prey population in the Holling type II functional responses. Based on those assumptions, a two-dimensional mathematical model is derived. The positivity, boundedness, and extinction of both prey and predator populations of the solution of the system are discussed. The existence, stability (local and global), and the Hopf bifurcation analysis of the biologically feasible equilibrium points are investigated. The aim of this research is to explore the effect of fear on the prey population and hunting cooperation on the predator population, and both prey and predator populations are harvested linearly and taken as control parameters of the model. If the values of , then both prey and predator populations are extinct and also fear parameter has a stabilizing effect on system 4. From the numerical simulation, it was found that the fear effect, hunting cooperation, prey harvesting, and predator harvesting change the dynamics of system 4.

#### 1. Introduction

One of the most challenging tasks in ecological model is understanding ecosystem dynamics. Population ecology is concerned with the developing theory and insight about the persistence, structure, and dynamics of biological communities. The predator-prey dynamics are one of the most important topics in theoretical ecology due to their universal existence and importance [1]. One must carefully construct a biologically meaningful and mathematically tractable population model [2]. The most important element in a population model is a predator-prey model which describes the number (density) of prey consumed per predator per unit time for given quantities (density) of prey and predator communities [3].

The harvesting of the prey species or the predator species or both is another important issue from an ecological and economic point of view. It generally has a strong impact on the dynamics of biological resources. Sometimes, harvesting of predators is used to control the predator population so as to prevent the extinction of the prey species. Concerning the conservation for the long-term benefits of humanity, there is a wide range of interest in the use of bioeconomic modeling to gain insight in the scientific management of renewable resources like fisheries and forestries [4].

The emerging view is that indirect effect (fear effect) on prey population is even more powerful than direct effect [5–7]. Due to fear of predation, prey animals can change its grazing zone to a safer place and sacrifice their highest intake rate areas [8–11]. At the last time, many researchers have proposed and studied fear effect of predator-prey model. Here, due to fear of predation, growth rate of prey reduces, i.e., increase fear effect implies decrease growth rate of prey [7, 11]. Recently, however, [10] conducted a manipulation on song sparrows during an entire breeding season to determine whether perceived predation risk could affect reproduction even in the absence of direct killing. For example, [7] have investigated the modeling of the fear effect in predator-prey model; it is considered as mortality rate of predator and the Holling type I functional responses and showed the local and global stability of the feasible equilibria and the existence of the Hopf bifurcation.

The role of mutual beneficial interactions (e.g., cooperation) among social animals is rapidly growing in the research field in population dynamics and conversion biology cooperative hunting among vertebrates (lions [12], African wild dogs [13], chimpanzees [14], wolves [15], birds [16], primates [17], etc.). Different animals adopt different strategies for making a group; for example, wolves have the strategies to get close to the prey and catch it for the acquisition of large and faster prey [18]. In general, animals of the same species compute with each other for resources, but because some specific mechanism is at work, they help each other for mutual benefit and show cooperative behavior. The hunting cooperation in predator-prey model has been investigated by many researchers [10, 11, 19–24].

In recent years, more attention has been paid to the study of the dynamics of hunting cooperation and fear in a predator-prey interaction [24–26]. For example, [26] design the mathematical model of predator-prey model with fear and hunting cooperation with the Lessie-Gower functional responses, and they investigated the local stability and bifurcation analysis. However, most of the works focus on the nonharvesting of hunting cooperation and fear effect of a predator-prey system. Thus, the main aim of this paper is to study the impact of fear effect in prey, hunting cooperation among predator, and prey and predator harvesting of a predator-prey system.

The organization of this paper is as follows: Section 2 is where the mathematical model formulation is discussed; Section 3 deals with the analysis, existence, and boundedness of the model, extinction criteria, the local and global stability analysis of the biologically feasible equilibrium points, and bifurcation analysis of the system; Section 4 deals with the numerical simulation results that are discussed; Section 5, lastly, deals with the conclusion part.

The novelty of this paper is the consideration of fear effect in prey, hunting cooperation among predator, and both prey and predator harvesting with nonlinear rate with the Holling type II functional responses.

#### 2. Mathematical Model Formulation

In this paper, the following assumptions are taken and a prey population is denoted by and predator population is denoted by . (1)The prey species grows logistically in the absence of predator with carrying capacity and internsic growth rate (2)The predator species feed their prey in the Holling type II functional response(3)The birth rate of prey population is reduced due to fear induced by predators(4)The predator population uses hunting cooperation among themselves for capturing the prey population(5)The predator natural death rate is assumed to be constant (6)The prey and predator are assumed to be harvested by an external force according to the linear type of functional harvesting (proportionate harvesting)

Based on the above assumptions and parameter definitions in Table 1, the model equation is given by the following system of nonlinear differential equations. and . All these parameters in equation (1) are assumed to be positive. To reduce the number of parameters, we nondimensionalized the model system using the following scaling:

The nondimensionalized form of equation (1) can be written as follows: where with .

#### 3. The Model Analysis

##### 3.1. Existence and Positive in Variance

For , let and , where . Then, system (4) can be written as , where with . Here, for . Thus, the function is locally Lipschitzian and completely continuous on . Therefore, the solution of system (4) with nonnegative initial condition exists and is unique. Moreover, it can be shown that these solutions exist for all and stay nonnegative. Hence, the region is an invariant domain of system (4).

##### 3.2. Boundedness

In theoretical ecology model, boundedness of a system implies that the system is biologically well-behaved. The following theorem ensures the boundedness of system (4).

Theorem 1. *The solution of system (4) which initiate in interior of is uniformly bounded.**Proof. The first equation of system (4) that represents the prey equation is bounded through, So by comparison theorem, we obtain for all , where . Thus, , which implies that for sufficiently large . Since is bounded, without any loss of generality it follows that the solutions of is bounded. To show this reality, we assure this mathematically as follows: we define a positive definite function as follows:
**From definition of , is differential in some maximal interval . Thus, differentiating with respect to time and substituting the reaction terms of system (4) yields
**Here, we have
**Applying the theory of differential inequality, we obtain
**For , we have . Hence, all the solution of system (4) starting in for any are contained in the region *

A system is said to be dissipative if all population initiating in are uniformly limited by their environment. Thus, the following theorem establishes the dissipativeness of system (4).

Theorem 2. *System (4) is dissipative.**Proof. System (4) is uniformly bounded. Therefore, system (4) is dissipative.*

##### 3.3. Extinction Criteria

Theorem 3. *If , then .**Proof. We have
**Integrating both sides and it becomes
**Therefore, , if .*

Theorem 4. *If , then .**Proof. We have
**Since , we have
**Integrating both sides and it becomes
**Therefore, we have , if .*

##### 3.4. Model Analysis

###### 3.4.1. Equilibrium Points

System (4) has the following biologically feasible equilibrium points. (1)Trivial equilibrium point (2)Axial equilibrium point which exists for (3)The positive interior equilibrium point exists for where

If exists, , and is the unique positive root of the polynomial equation with

One can see that exists for using the Discartes rule of sign; change equation (17) that has a unique positive solution; the following conditions are satisfied.

##### 3.5. Stability Analysis

In order to investigate the local stability property of system (4), we shall calculate the Jacobian matrix at each equilibrium. The Jacobian matrix of system (4) at any arbitrary equilibrium is given as where

Theorem 5. *The trivial equilibrium point is stable if .**Proof. The eigenvalues of the Jacobian matrix of system (4) at are . All the eigenvalues are negative if and only if . Therefore, system (4) at the equilibrium point is locally asymptotically stable.*

Theorem 6. *The equilibrium point of system (4) is locally asymptotically stable if and .**Proof. The eigenvalues of the Jacobian matrix of system (4) at the axial equilibrium point are . Therefore, the equilibrium point is locally asymptotically stable if and .*

Theorem 7. *The coexistence equilibrium point is locally asymptotically stable if
where
**Proof. The Jacobian matrix of system (4) around the interior equilibrium point is given by where are given by (23). Let be the roots of the characteristics equation of which is given by
where
**By the Routh-Hurwitz criterion, is locally asymptotically stable if and only if and 0. Then, we have , for . Therefore, the coexistence equilibrium point is locally asymptotically stable if the conditions in (22) hold.*

##### 3.6. Global Stability Analysis of the Coexistence Equilibrium Point

Theorem 8. *System (4) does not admit any periodic solution around the positive equilibrium point provided
**Proof. Let be solution of system (4). Define the Dulac function
**Then, we have
, provided (28) holds.**Thus, by the Dulac criterion, system (4) does not have a nontrivial periodic solution. Therefore, the local asymptotically stability of system (4) around and the boundedness of the solutions of system (4) ensure its globally asymptotically stability around the coexistence equilibrium point.*

##### 3.7. Hopf Bifurcation

Characteristic equation of system (4) at is given by where with defined as in (23).

Theorem 9. *System (31) possesses a Hopf bifurcation around when passes through provided that **Proof. At , the characteristic of equation (4) can be written as , which has roots and . Thus, there exist a pair of purely complex eigenvalues. Also and are the smooth functions of . So for neighborhood of , the roots have the form and , where are real functions for . Next, we verify the transversality condition.
**Now, substituting into (31), derivating it along , and comparing real and imaginary parts, respectively, we get
**Simplifying it and it becomes
where
**Thus, we have
**At , we get
**Hence, the transversality condition holds. Therefore, the Hopf Bifurcation occurs at .*

#### 4. Numerical Simulation

In the present study, the rate of prey harvesting , the rate of predator harvesting , the fear factor coefficient , and hunting cooperation coefficient are the key parameters which will be taken as control parameters. The numerical simulation is carried out using the MATLAB software package for the set of parameter values given in Table 2.

##### 4.1. Effect of Varying the Prey Harvesting Rate

Let us fix the variables in Table 2 as .

It is observed that system (4) approaches to the trivial equilibrium point for , to the predator-free equilibrium point (axial equilibrium point) for , and to the coexistence equilibrium point for .The illustration of the coexistence equilibrium point exists and is locally asymptotically stable if , when lies in the range of ; then, the coexistence equilibrium point annihilates with the predator-free equilibrium point (axial equilibrium point) , and changes its stability and becomes locally asymptotically stable. Further increasing lies in the range of which leads to the predator-free equilibrium point that annihilates with the trivial equilibrium point , and changes its stability and becomes locally asymptotically stable.

Figures 1(a) and 1(b), respectively, show that the time series and parametric plot of the trajectory of system (4), i.e., it approaches locally asymptotically to the coexistence equilibrium point for . Figure 2(a) shows that the trajectories of system (4) approach locally asymptotically to the predator-free equilibrium point (axial equilibrium point) for , i.e., one can see that the prey species can coexist whereas the predator species goes to extinction as the value of . Figure 2(b) shows that the trajectories of system (4) approach locally asymptotically to the trivial equilibrium point stability of for , i.e., one can see that both species are extinct as the value of . Figure 3(a) shows that an increase in the prey harvesting rate leads to an increase in the prey population, whereas Figure 3(b) illustrates that as the rate of harvesting increases, the density of the predator population decreases.

**(a)**

**(b)**

**(a)**

**(b)**

**(a)**

**(b)**

##### 4.2. Effect of Predator Harvesting Rate

For the parametric values as in Table 2 with , the predator-free equilibrium point (axial equilibrium point) and the coexistence equilibrium point exists for and , respectively. The illustration of the coexistence equilibrium point exists and is locally asymptotically stable in . And lies between in the range ; the coexistence equilibrium point annihilates with the predator-free equilibrium point (axial equilibrium point) , and changes its stability and becomes locally asymptotically stable.

Figures 4(a) and 4(b) show the stability of the time series and parametric plot of system (4), respectively, i.e., it approaches to the coexistence equilibrium point for . Figure 5(a) shows the stability of system (4) approaches to locally asymptotically to the predator-free equilibrium point (axial equilibrium point) for . Furthermore, one can see that the prey species coexist whereas the predator species goes to extinction whenever . From Figures 6(a) and 6(b), it can be observed that an increase in the harvesting rate of predator leads to increase in the prey population density and a decrease predator population density, respectively.

**(a)**

**(b)**

**(a)**

**(b)**

##### 4.3. Effect of Hunting Cooperation

The figures show that the population size of both prey and predator oscillate, but they are not periodic. The amplitudes of the oscillation of both populations approach the equilibrium point of the system.

Figures 7(a) and 7(b) show the stability of the time series and parametric plot of system (4), respectively, i.e., it approaches asymptotically to the coexistence equilibrium point for . Figures 8(a) and 8(b) show the instability of time series and parametric plot of the system, respectively, i.e., it does not approach asymptotically to the equilibrium point for . From Figures 9(a) and 9(b), it can be observed that an increase in the hunting cooperation of predator leads to a decrease in the prey population density and predator population density. It is observed that for the parameter values given in Table 2 with , the trajectory of system (4) approaches the coexistence equilibrium point for and the Hopf bifurcation occurs at . Otherwise, it is unstable.

**(a)**

**(b)**

**(a)**

**(b)**

**(a)**

**(b)**

##### 4.4. Effect of Fear

Now, to investigate the impact of fear effect on the dynamics on system (4), for the parametric values with , the coexistence equilibrium point exists for ; otherwise, the predator-free equilibrium point exists.

The coexistence equilibrium point exists and is locally asymptotically stable in . And lies within the range ; the coexistence equilibrium point annihilates with the predator-free equilibrium point (axial equilibrium point) , and changes its stability and becomes locally asymptotically stable.

Figures 10(a) and 10(b) show the stability of the time series and parametric plot of system (4), i.e., it approaches asymptotically to the coexistence equilibrium point for . Figure 11(a) shows that the stability of system (4) approaches asymptotically to the predator-free equilibrium point (axial equilibrium point) for . Furthermore, one can see that the prey species coexist whereas the predator species goes to extinction as the value of is out of the range . From Figures 12(a) and 12(b), it can be observed that an increase in the fear strength leads to an increase in the prey population density and a decrease in the predator population density.

**(a)**

**(b)**

**(a)**

**(b)**

Here, we choose different values of ; system (2) shows stable focuses and periodic oscillation. If we incorporate with fear strength, a fear phenomenon into the system for a fixed hunting cooperation , the unstable system becomes stable for a fixed values of other than . Figures 13(a) and 13(b) show the instability of system (4) about the coexistence equilibrium point whenever for a fixed value of . Figure 14 shows the stability of system (4) about the coexistence equilibrium point whenever for a fixed value of . We can observe that fear factor () has a stabilizing effect on the dynamics of the system (c.f. Figures 13(a) and 13(b) and Figures 14(a) and 14(b)).

**(a)**

**(b)**

**(a)**

**(b)**

Further, we investigate the global stability of system (4) for with different initial conditions , which shows the global converges of the trajectories of system (2) to the coexistence equilibrium point (c.f. Figure 15). This shows the globally asymptotically stability of the interior equilibrium point .

On the analysis of system (4), one can see that, for the parametric values as in Table 2 and , the endemic equilibrium point is locally asymptotically stable for and unstable for . The existence of the Hopf bifurcation occurs at (see Figure 16).

#### 5. Conclusion

In this paper, it can be conclude that the formulated model is mathematically meaningfully valid and biologically well posed by providing boundedness positivity and existence of the solution of the model. Trivial, axial, and coexistence equilibrium points are investigated. Moreover, it is observed that in our model, trivial equilibrium point is always locally asymptotically stable for . Axial equilibrium point is locally asymptotically stable if and only if the variable satisfies the condition and . The local and global stability condition for the coexistence equilibrium point of system (4) are obtained. The prey population will be extinct if the prey growth rate is less than the prey harvesting rate (see Theorem 3) and the predator population extincts (Theorem 4 holds). A bifurcation analysis of system 4 undergoes a Hopf bifurcation under a certain condition.

Numerical simulations verify our analytical findings. The following main results are obtained for the parametric values given in Table 2.

Increasing prey harvesting lies in the range which leads to the extinction of predator and it is locally asymptotically stable to the axial equilibrium point (predator-free equilibrium point) (see Figure 2(a)). Again, further increase to the range results in the extinction of both the prey and predator populations and it is locally asymptotically stable to the trivial equilibrium point (see Figure 2(b)). Finally increase lies in the range and leads to increase in prey population and decreased in predator population; then, system 4 is locally asymptotically stable to the coexistence equilibrium point (see Figures 3(a) and 3(b)), respectively.

Increase predator harvesting lies in the range between which leads to the extinction of predator population and it is locally asymptotically stable to the axial equilibrium point (predator-free equilibrium point) (see Figure 5(a)). Again, increase lies in the range which increases the prey population and decreases the predator population and it is locally asymptotically stable to the coexistence equilibrium point (see Figures 6(a) and 6(b)).

Increase hunting cooperation of predator population lies in the range which leads to decrease both the prey and predator populations and it is locally asymptotically stable to the coexistence equilibrium point (see Figures 9(a) and 9(b)). The coexistence equilibrium point loses its stability, and the Hopf bifurcation occurs at ; otherwise, it is unstable.

Increase fear effect of prey population lies in the range out of which leads to the extinction of predator population and it is locally asymptotically stable to the axial equilibrium point (predator-free equilibrium point) (see Figure 11(a)). Again, further lies in the range which increases the prey population and decreases the predator population and it is locally asymptotically stable to the coexistence equilibrium point (see Figures 12(a) and 12(b)). Finally, we can observe that fear factor has a stabilizer effect on the dynamics of system 4 (see Figures 13(a) and 13(b) and Figures 14(a) and 14(b)).

#### Data Availability

The data used to support the findings of this study are included within the article.

#### Conflicts of Interest

The authors declare that there is no conflict of interest regarding the publication of this paper.