Impact of Fear and Habitat Complexity in a Predator-Prey System with Two Different Shaped Functional Responses: A Comparative Study

Habitat complexity or the structural complexity of habitat reduces the available space for interacting species, and subsequently, the encounter rate between the prey and predator is decreased significantly. Different experimental shreds of evidence validate that the presence of the predator strongly affects the physiological behaviour of prey individuals and dramatically reduces their reproduction rate. In this study, we investigate the interplay between the level of fear and the degree of habitat complexity in a predator-prey model with two different shaped functional responses. We, therefore, develop the functional response using the timescale separation method, and the shape of the resulting functional response depends upon the monotonous property of catch rate, g(N) where N is the prey biomass. Whenever g(N) increases strictly, a saturating functional response occurs, but for nonmonotonic g(N), a dome-shaped functional response arises. For saturating case, it has been revealed that both prey and predator biomass may oscillate for lower levels of fear and a lower degree of habitat complexity. To stabilize this oscillatory behaviour to a coexistence state, we have to adequately increase the level of fear or degree of habitat complexity. However, for dome-shaped case, more complicated dynamics are observed. In this case, coexistence steady state, if exists, may be locally asymptotically stable for a lower degree of habitat complexity, but for intermediate values, the system is capable of producing multiple coexistence steady states with a bistable phenomenon between predator-free steady state and a coexistence steady state. Moreover, if the level of fear is sufficiently low, the system may experience a supercritical or/and subcritical Hopf bifurcation. In the dynamics of parametric disturbance for the degree of habitat complexity parameter, dome-shaped functional response predicts that disturbance may trap the system into a nearest attractor (either a large amplitude stable limit cycle or predator-free steady state); this can be overcome only by a larger alteration, or sometimes it is impossible to overcome (hysteresis phenomena), whereas the saturating-shaped functional response predicts a system resilience. For both the functional responses, a higher degree of habitat complexity always increases the extinction possibility of the predator, and no level of fear can compensate this biodiversity loss.


Introduction
In ecology, food is one of the basic substances essential for maintaining life and growth of a species. For the purpose of collecting food, members of one species may destroy the members of the other, which establish the most common predator-prey relationship in ecology. From the pioneering invention of Lotka-Volterra model [1,2], which is a concrete mathematical elucidation of a predator-prey relationship, various rich and prosperous ecological systems have been studied extensively to understand and interpret the dynamical consequences of interacting populations [3][4][5][6][7][8].
Most of these researches have emphasized qualitative analyses and dynamical behaviour together with stability properties, bifurcation analysis, creation of limit cycle, and permanence. ese investigations reveal that the qualitative dynamics of a predator-prey model critically depends on the way of predator's response towards the changing biomass of its prey [3]. Functional responses or numerical responses are widely used to describe the prey consumption rate as a function of prey density (biomass). In general, C. S.  classifies the functional responses in three types [9,10]. In type I, consumption rate is linearly increasing with prey density, while in type II, predator is restricted by its capacity to process food, and the number of consumed prey saturates at high levels of prey density. Type III is characterized by a sigmoidal relation in which the consumption rate is positively density-dependent over some regions of prey density. Also, an extensive range of other functional responses, such as Crowley-Martin, Monod-Haldane, and Beddington-DeAngelis, are derived as well [11][12][13][14]. Moreover, it has been observed that similar shaped functional responses may not necessarily predict the similar model dynamics [15]. For instance, Aldebert et al. [16] explained the different model predictions like equilibrium versus oscillations, single attractor versus multiple attractors, and resilience versus hysteresis for two similar shaped functional responses: Holling type II and Ivlev's functional responses. is phenomenon is known as structural sensitivity, which is one of the significant impediment to enhancing the predictive power of ecological models.
In spite of being one of the most instinctive concepts in ecology, understanding the role of habitat complexity, which is critically significant for maintaining biodiversity, community persistence, and ecosystem stability, is still an expanding branch in modern research. Populations survive in natural habitats which are spatially or environmentally heterogeneous and could be split up into small patches with various attributes like different levels of resources and environmental capacity. In this article, we define the heterogeneity of habitat as habitat complexity, which refers to the morphological characteristics within a structure itself or the heterogeneity in the arrangement of objects [17]. Habitat complexity is observed in both terrestrial or aquatic ecological systems [18,19]. For example, in marine habitat, complexity is increased in the presence of oyster and coral reefs, sea grass beds, mangroves, and salt marshes [20]. In lakes, habitat heterogeneity is most commonly present in the form of littoral zone vegetation or depth-gradient diversity [21]. MacArthur and MacArthur (1961) first performed the quantitative analysis of the impact of habitat complexity on species diversity, as the repercussions of which help to explore the local effects of habitat complexity related to increased richness and abundance in more complex habitats [22]. Lots of experimental studies reveal that habitat complexity may decrease the predation rates by reducing the encounter rates or attack rates between predator and prey [23,24]. In a series of laboratory experiments, Manatunge et al. [25] observed the behaviour of planktivorous fish Pseudorasbora parva (Cyprinidae) and their zooplankton prey Daphnia pulex with artificial vegetation at densities of 0, 350, 700, 1400, 2100, and 2800 stems.m − 2 and recorded the foraging rates of the fish at different prey densities for all stem densities. ey noticed that the foraging efficiency of P. parva decreased significantly with increasing habitat complexity. In a recent study, Bartholomew et al. (2000) proposed two nondimensional indices C t /A t and Sp/Pr that directly measure how the structural complexity interferes with the foraging ability of predator [26]. C t /A t estimates the amount of cover available within a habitat that interferes with predator's capability to see, and Sp/Pr measures the extent to which the structure interferes with predator's potentiality to move through the habitat in search of prey. Here, C t is the total covered area within a habitat, A t is the total area of the habitat, Sp is the average interstructural space size, and Pr is the size of the predator. e field experiment with the fish (Fundulus heteroclitus) as predator and amphipods as prey validate that prey survivorship will increase with increasing C t /A t and that survivorship will decrease with increasing Sp/Pr. In mathematical ecology, researchers are incorporating the effect of habitat complexity in a predator-prey system by modifying predator response functions, since this complexity can decrease the encounter coefficient of predators [17,27,28]. eoretical ecologists and evolutionary biologists have clarified that the presence of predator may affect the physiological behaviour of prey species to such an extent that it can be more effective in reducing prey biomass than direct predation [29,30]. However, due to the lack of experimental evidences, this behavioural change of prey species has not been incorporated in a mathematical model. Recently, a field experiment on song sparrows (Melospiza melodia) validated the effect of fear on prey population [31]. In that experiment, it is observed that the number of offspring produced per year is reduced by around 40% due to perception of predation risk. is reduction occurred due to the antipredator activities of birds such as high levels of perceived predation risk which can cause prey to decamp from their habitat temporarily; prey species are worried and frightened to come to an open habitat, so they do not have any free atmosphere to do their daily works including coupling. erefore, the costs of fear may additionally decrease prey fecundity and survival, and the overall reduction in prey numbers may exceed compared to direct predation alone. In another field experiment on Drosophila melanogaster as prey species and mantid as their predator species, Elliott et al. [32] observed that the presence of mantids decreased the reproductive performance of drosophila in their breeding as well as nonbreeding seasons. erefore, all these experimental pieces of evidence confirm that the perception of predation risk is persuasive enough to control the community structure and should thus be given greater attention in ecological conservation and management. ese experimental findings motivate the researchers to improve the predator-prey model by introducing the fear effects in prey population [33][34][35][36]. Wang et al. [37] considered a predator-prey model together with the cost of fear into the prey reproduction, and their analyses exhibit that a high level of fear stabilizes the system by ruling out the periodic oscillations. It has also been proved that the growing birth rate of prey or growing death rate of predator reduces the sensitivity of prey species in perceived predation risk. Wang and Zou [38] proposed a stage-structured predator-prey model with fear effect and adaptive avoidance of predators. Mathematical inspections 2 Discrete Dynamics in Nature and Society analyzed the role of fear and the maturation delay between juvenile prey and adult prey in determining the long-term population dynamics. ey proved that both high level of fear and strong adaptation of adult prey have a destabilizing effect, while large population of predators has a stabilizing effect on the predator-prey interactions. Most of the previous studies have considered that for high level of fear or larger predator biomass, fear function may affect the reproduction of prey to such an extent that the number of offspring decreases to zero. However, this is not the real scenario observed in various experimental evidences. erefore, in this work, we assume that even if the level of fear or size of predator population increases to a larger quantity, the prey population will not stop reproduction. It is biologically significant since the prey population become conscious (or aware of ) and habituated to the presence of predator.
In the above, we have discussed various biological aspects that can make the dynamics of a predator-prey model more complex and rich. Arguably, it is observed that the complexity of a predator-prey model is controlled by the choice of functional response. In this study, we construct a functional response utilizing the quasi-steadystate assumptions in which the dynamical model is simplified based on the idea that the predation processes happen in a much faster timescale than the birth and death processes of predator. We shall show that the monotonic property of catch rate shapes the resulting functional response in different forms. For example, the functional response will be saturating (dome-shaped) if the catch rate is strictly increasing (initially increases and then decreases) with prey density. Further, we incorporate the fear effect in prey reproduction and modify the encounter rate of predator by introducing the habitat complexity. Now, the questions we arise in this paper are the following. (i) How does the choice of two different shaped functional responses affect the qualitative behaviour of a system? (ii) How much the structural sensitivity can influence the dynamics of an ecological model? (iii) How do the habitat complexity and cost of fear jointly affect the population dynamics?
Inspection of the predator-prey interaction can be exploited using two major pathways. e most familiar pathway involves numerous methodical steps to characterize the nature of a specific predator-prey relationship by split up the interaction into various fundamental components. Although this approach helps to grasp the insights completely of a specific interaction, it is not fruitful to make that consequence in general statements. For that reason, the second pathway, that is, that of generalization, has been endeavoured. In this paper, we initially analyze the predatorprey model analytically by considering a general fear function and a prey-dependent catch rate and later examine whether this generalization is appropriate for a specific interaction. In recent times, a lot of species are facing the threat of extinction due to rapid climate change around the globe, pollution, or habitat loss, which may affect the ecological balance and biodiversity. So, coexistence in the ecosystem and maintaining biodiversity are the primary intentions of our investigation. In this paper, the study of multiple steady states receives special attention in the context of ecological conservation [39]. It explores the sensitivity of the system behaviour when that system faces some disruptions.
We arrange this work in the following manner. In Section 2, we upgrade the classical predator-prey model by incorporating the cost of fear (in prey reproduction) and habitat complexity that affects the encounter rate of predator. In this model, instead of using a particular functional response, a general response function has been employed whose graphical shape depends on the monotonic property of the catch rate. All the basic dynamical properties of the model, including positivity and boundedness of the solutions, permanence of the system, and local stability of the possible equilibrium points, are studied in this section. e entire discussion of Section 3 is divided into two sections. Section 3.1 is devoted to analyze the model dynamics with a particular fear function and a monotone increasing catch rate, whereas in Section 3.2, catch rate is taken as a domeshaped form. In both sections, the feasibility of steady states, the stability analysis, bistability phenomenon (if exists), local bifurcations, and existence of limit cycles are discussed in detail. Our model analysis also includes the effect of fear and habitat complexity in forming the population demography. All of the theoretical findings are convinced through numerical simulations, which show some potential roles of fear effect and habitat complexity that may significantly control the dynamics of predator-prey interplay. In Section 4, we discuss and summarize the biological significance of our analytical findings.

The Model
Our model is constructed by ameliorating the following classical predator-prey system with logistic growth in prey population: with initial conditions N(0) � N 0 > 0 and P(0) � P 0 > 0. Here, N(t) and P(t) are the respective population size of prey and predator at time t. Other parameters bear the following biological interpretations: r, d, and a are the reproduction rate, natural death rate, and intraspecific competition rate, respectively, of the prey population. e term f(N) represents a functional response that quantifies the variation of captured prey per unit time by one predator on average with changing densities; θ(0 < θ < 1) is the conversion efficiency and θf(N)P measures the biomass production of P due to predation; m is the mortality rate of predator. e most conventional approach to formulate a functional response, f(N), involves various time budgets of predator for searching prey and handling the victim prey. For an alternative approach, we assume that the whole predator community is divided into Discrete Dynamics in Nature and Society two different states: searching (S) and handling (H). erefore, P � S + H. Also, note that the transformation of searching individuals into handling individuals or vice versa occurs on a much faster timescale than the birth and death processes of predator. Hence, the dynamics of these alternations of states are given by with initial conditions S(0) � S 0 ≥ 0 and H(0) � H 0 ≥ 0. But both S(0) and H(0) are not simultaneously zero, since P(0) � S(0) + H(0) > 0. Biologically, the term g(N) represents the successful catch rate per searching predator, whereas c represents the attack rate or encounter rate. It is biologically feasible to consider g(0) � 0, that is, in the absence of prey, predator has nothing to catch. e time required for killing, consuming, and digesting prey (or in together handling prey) is symbolized as h � c − 1 , where c is the handling rate. After this handling time, handling individuals again turn back into searching individuals. Now, in this predator-prey system, we incorporate the effect of habitat complexity which is more likely to affect the attack coefficient c [23,24]. Hence, we replace the attack coefficient is a dimensionless parameter that quantifies the reduction of encounter rate due to habitat complexity and is interpreted as the degree or strength of habitat complexity. Hence, system (2) is modified as Implementing the timescale separation and holding P � S + H � constant for this faster timescale, a quasi-stationary solution for the searching predator is obtained as .
We also presume that the predation depends only on available searching individuals, and then the total number of victim prey is given by which introduces the functional response, f(N), as .
It should be mentioned that the monotonic property of g(N) plays a crucial role in characterizing different shapes of the functional response. For strictly increasing g(N), the derived functional response will be a saturating functional response. Otherwise, if g(N) is growing initially and after some certain prey density, due to group defense, g(N) will decrease and in that case, a domeshaped functional response will occur. From the practical perspective, derivation of the functional response in this manner and not to include it directly into the model has significant advantages. For many experiments, one may have to retrieve the catch rate from some artificial laboratory experiments and for that case, it is sufficient to initiate the experiment with an entirely searching (not satiated) predator population into a prey population of different sizes and extract the catch rate depending on the prey population.
Various empirical research studies on the predator-prey system reveal that prey's reproduction rate significantly decreases due to scare of being caught by the predator. Incorporating this biological effect, we modify the model by multiplying the reproduction term by a function ϕ(k, η, P) that accounts for the cost of fear, where k represents the level of fear and η ∈ (0, 1) describes the maximum cost of fear. It is biologically relevant to consider the following conditions on ϕ(k, η, P): ϕ(0, η, P) � ϕ(k, η, 0) � 1 means that in the absence of fear or in the absence of predator, the fear function does not affect the growth of prey population; (zϕ/zk) < 0, (zϕ/zP) < 0 signify that increase of fear level or enlarging the predator species impacts the negative effect on the reproduction of prey population; lim k⟶∞ ϕ(k, η, P) � lim P⟶∞ ϕ(k, η, P) � η implies that even if the level of fear or size of predator population increases infinitely large, the prey population will stress under certain level of fear due to physiological impact when the prey populations are accustomed to fear from predator species. In the experiment of Zanette et al. [31], although the reproduction of offspring songbirds reduces 40% per year due to fear of predator, they have not completely stopped the reproduction and that motivates to consider lim k⟶∞ ϕ(k, η, P) � lim P⟶∞ ϕ (k, η, P) � η instead of lim k⟶∞ ϕ(k, η, P) � lim P⟶∞ ϕ (k, η, P) � 0. It is biologically significant since the prey population is conscious and habituated to the presence of predator. Combining all of these aspects and incorporating these in (1), our resulting model is given by the following system of nonlinear ordinary differential equations: , 2.1. Positivity and Boundedness. To examine the well-posedness of system (7), we prove that every positive solution remains positive and will not grow abruptly for long time interval. Biologically, positivity and uniform boundedness assert that interacting species in our model are ecologically well behaved, and abundance of each species is restricted due to limited resource. (7) is positively invariant and uniformly bounded for any time t ≥ 0.

Theorem 1. Each solution of system
Proof. It can be easily proved that the RHS of system (7) is continuous and locally Lipschitzian in the domain R 2 + and hence the solution (N(t), P(t)) of the system exists uniquely on [0, ξ), where 0 < ξ ≤ + ∞. Further, it is assumed that lim N⟶0+ (g(N)/N) � g ′ (0) ≠ 0. en, from system (7) with initial conditions N(0) > 0 and P(0) > 0, we have is proves that the system is positively invariant for any time t > 0.
To prove the boundedness of the solution, we consider the function W � N + (1/θ)P. en, For each W(t) > 0, we get � Ω(say).

(10)
Using the theory of differential inequality [40] for W(t), Taking limit when t ⟶ ∞, we obtain 0 < W(t) ≤ (Ω/λ). Hence, every positive solution in R 2 + confined in the region B � (N, P) ∈ R 2 + : 0 < N + (P/θ) ≤ (Ω/λ) . is shows that the solution of the system is bounded. □ 2.2. Uniform Persistence. From mathematical point of view, persistence of a system assures that the solutions of that system are always away from zero. Ecologically, it means the sustainability of all the interacting species for a long time interval. Here, we prove the permanence of system (7) directly by utilizing average Lyapunov function [41].
Proof. For (N, P) ∈ R 2 + , we consider the following average Lyapunov function: where σ 1 and σ 2 are some positive constants. Now, differentiating with respect to t, we have Discrete Dynamics in Nature and Society e system will be uniformly persistent if Δ(N, P) > 0 at the boundary equilibrium points (0, 0) and ((r − d)/a, 0) for some σ 1 > 0 and σ 2 > 0. Now, evaluating Δ at the boundary equilibrium points, we obtain erefore, for some positive □ Note that the cost of fear in prey reproduction does not affect the persistent of the model system under consideration.

Equilibria and Its Stability.
Ecological equilibrium is a dynamical state at which a community of organisms ceases to grow. In predator-prey system, equilibrium occurs at the intersection of their respective nullclines in the nonnegative quadrant. Hence, all the steady states satisfy (dN/dt) � 0 and (dP/dt) � 0 simultaneously. Two boundary steady states of this model are given by Note that higher mortality rate is harmful to any species and may lead that species to extinction, which is biologically intuitive, so throughout this paper, it is assumed that r > d. Due to biodiversity loss, these steady states are harmful to the effective functioning of the ecosystem. So, coexistence in the ecosystem and maintaining biodiversity are the primary intention of our investigation. Existence of nontrivial steady state E 2 � (N * , P * ) depends on the form of g(N) and ϕ(k, η, P), and this takes the form of By definition, g(N) is a catch rate and hence, g(N * ) > 0. en, for existence of positive N * , (15) yields θ > mh. (17) If g(N) is bounded, then another additional requirement is From the biological viewpoint, (17) signifies that the mortality rate of predator needs to be smaller than the product of conversion efficiency θ and handling rate 1/h, which are both related to predation capabilities. However, we presume that searching and handling prey happen on a much smaller timescale than birth and death processes; hence, this inequality is feasible and likely holds. It is worth mentioning that increase of searching rate c fails to compensate for lower handling rate concerning the existence of the interior steady state. Since Z � ϕ(k, η, P) is decreasing in P with ϕ(k, η, 0) � 1 and Z � (d/r) + (aN * /r)+ (m/rθN * )P is a straight line with Z intercepts (d/r) + (aN * /r); they intersect at exactly one point for each N * , provided So, the number of coexistence steady states depends on the form of g(N) only: if g(N) is strictly increasing, then precisely one E 2 (N * , P * ) exists; otherwise, for dome-shaped g(N) at-most two E 2 , namely, E 2 1 (N * 1 , P * 1 ) and E 2 2 (N * 2 , P * 2 ) occur. For multiple cases, if we assume N * 1 < N * 2 , then g(N * 1 ) must lie in the increasing branch, and g(N * 2 ) must lie in the decreasing branch of g(N) (see Figure 1). Also, it should be noted that the number of coexistence steady states is not affected by the level of fear and from (16), we obtain which means the predator biomass of coexistence equilibrium point decreases due to increased level of fear. For the purpose of inspecting local stability around the steady states, we consider the following Jacobian matrix: Eigenvalues at the trivial equilibrium point E 0 are evaluated as (r − d)( > 0) and − m( < 0). Hence, E 0 is always a saddle point. At the semitrivial equilibrium point E 1 , Jacobian matrix has the eigenvalues λ 1 If θ < mh, no coexistence steady states exist and the semitrivial solution is locally asymptotically stable (LAS), but this inequality hardly holds. Hence, for stability of E 1 , we must have erefore, we have the following theorems.

Theorem 4.
e axial equilibrium point Further, we assume that the coexistence is possible and so N * < ((r − d)/a). If g(N) is strictly increasing in N, then exactly one coexistence steady state E 2 (N * , P * ) exists and g( , which establishes that existence of E 2 forces E 1 to be an unstable steady state. However, for a domeshaped predation rate, multiple coexistence steady states, namely, E 2 1 (N * 1 , P * 1 ) and E 2 2 (N * 2 , P * 2 ), with N * 1 < N * 2 may exist. In that multiple case, N * 2 < ((r − d)/a)⇒g((r− d)/a) < (m/c(1 − α)(θ − mh)) and hence E 1 exhibits stable behaviour; that is, bistability may occur between coexistence and semitrivial or axial equilibrium point. Figure 1: Geometrical explanation for existence of coexistence steady state(s). Here, N * is the intersection of the curve F 1 (N) � g(N) and the horizontal line F 1 (N) � (m/c(1 − α)(θ − mh)). For this particular N * , P * is the meet point of the curve F 2 (P) � ϕ(k, η, P) and the line F 2 (P) � (d/r) + (aN * /r) + (m/rθN * )P. (a, b) e horizontal line F 1 (N) � (m/c(1 − α)(θ − mh)) intersects the increasing g(N) at exactly one point provided conditions (17) and (18) are fulfilled, and the corresponding P * will be possible if (19) holds. (c, d) Multiple coexistences may occur when g(N) is a nonmonotonic function. Particular forms of g(N) and ϕ(k, η, P) are mentioned in the figure. Parameters values are taken as Values in the first bracket corresponding to c and m are used for plotting the last two figures.
Since the explicit form of E 2 (N * , P * ) cannot be assessed from (15) and (16), we address the local stability of E 2 using the sign of trace and determinant of the Jacobian matrix evaluated at E 2 . After some algebraic simplifications using (15) and (16), we obtain Det Since (zϕ/zP) < 0, for positive determinant, g ′ (N * ) > 0 must holds. For strictly increasing catch rate, this is always true and hence, for this case, stability of coexistence steady state (exactly one E 2 (N * , P * ) exists) depends on the sign of Tr(J E 2 ). Coexistence steady state will be locally asymptoti- For the either case, that is, for unstable coexistence state, there must be at least one limit cycle that occured [42]. Next, for domeshaped g(N), multiple coexistence steady states E 2 1 (N * 1 , P * 1 ) and E 2 2 (N * 2 , P * 2 ) with N * 1 < N * 2 may occur, and the corresponding trace and determinant can be easily calculated from (23) and (24). In that case, g(N * 2 ) always lies in the decreasing branch of g(N), and then g ′ (N * 2 ) < 0, and this implies that the determinant at E 2 2 is less than zero.
erefore, E 2 2 is always a saddle point. Stability behaviour of E 2 1 is the same as E 2 in monotone increasing case. Also, for this multiple case, mh)), and hence, axial steady state is locally asymptotically stable. erefore, for dome-shaped functional response, the system is capable of producing a bistable phenomenon between E 2 1 and E 1 .

Theorem 5
(i) If g(N) is strictly increasing, precisely one coexistence steady state, E 2 (N * , P * ) exists and will be locally asymptotically stable provided Tr( We summarize all of these parametric conditions related to the existence of steady states and their stability in Table 1. en, (dP/dt) � 0 only when P � 0, and therefore no coexistence equilibrium point is possible for θ � mh. Further, we have (dP/dt) < 0 for nonzero P as 0 < α < 1 and g(N) ≥ 0. erefore, the predator biomass decreases with time and the predator-free equilibrium point E 1 will be locally asymptotically stable. From the eigenvalue analysis, we can observe that both eigenvalues of the Jacobian matrix evaluated at E 1 are negative for θ � mh, which implies that E 1 is locally asymptotically stable.
implies both eigenvalues of the Jacobian matrix evaluated at E 1 are negative. Hence, E 1 is locally asymptotically stable. mh)), the model system experiences a saddle-node bifurcation at these parametric values. Detailed analysis of this bifurcation has been discussed later by taking a particular form of g(N) as g(N) � N/(1 + N 2 ).

Model with a Given Fear Function and Catch Rate
is section executes some extensive analysis of model (7) with some particular form of fear function and catch rate. Previously, we have already mentioned that the reproduction term of prey is modified by a fear function ϕ(k, η, P), where k is the level of fear and η ∈ (0, 1) is the maximum cost of fear. Since prey population is conscious and habituated with the presence of predator, we take Here, we assume the following particular form of ϕ is which agrees with all the conditions of ϕ as indicated earlier.
Diverse biological phenomena, such as saturation and group defense, may arise from the monotonic property of g(N), which characterize the different shapes of the functional response. Here, we investigate model (7) for two different forms of g(N): one is increasing-shaped, and the other one is dome-shaped.

Strictly Increasing Catch Rate.
Assuming g(N) � N for which the resulting functional response, f(N), becomes a saturating functional response and taking ϕ(k, η, P) as (26), we particularise model (7) as with initial population size N(0) � N 0 and P(0) � P 0 . It has been already proved that the model (27) is positively invariant and uniformly bounded for any time t ≥ 0. Two boundary equilibrium points that always exist are given by (i) the trivial equilibrium E 0 � (0, 0) and (ii) the predator-free equilibrium 8 Discrete Dynamics in Nature and Society � N is strictly increasing and unbounded, exactly one coexistence equilibrium, E 2 (N * , P * ), exists provided θ > mh, and N * < (r − d)/a. From (15) and (16), we obtain N * � (m/c(1 − α)(θ − mh)) and P * is the positive root of the following quadratic equation: e condition N * < (r − d)/a assures that the above equation has exactly one positive root. On the other hand, if hence, no positive P * exists. From the existence conditions of E 2 , we can define an existence region A 1 in α − θ parametric plane: where α � 1 − (am/c(r − d)(θ − mh)) and θ � mh + (am/c(r − d)).
e trivial equilibrium point E 0 is always an unstable saddle. Predator-free equilibrium point E 1 will be locally asymptotically stable if (r − d)/a < (m/c(1 − α)(θ − mh)), which can be deduced directly from (22). Interestingly, stability of E 1 ensures the nonexistence of E 2 .
us, the stability region of E 1 in α − θ parametric plane is given by In the previous section, it has been established that stability of coexistence steady state for strictly increasing catch rate exclusively depends on the sign of Tr(J E 2 ), where J E 2 is the Jacobian matrix for E 2 . Substituting g(N * ) � N * in (23) and using the equation of nontrivial predator nullcline, we obtain Hence, E 2 is locally asymptotically stable if Tr(J E 2 ) < 0. When k � 0, one can easily derive that system (27) will be Exactly one P * exists 1 No P * 1 and P * 2 exist 0 Only P * 1 exists, P * 2 does not exist 1 Both P * 1 and P * 2 exist 2 E 0 : unstable E 1 : LAS E 2 1 : LAS/ limit cycle E 2 2 : unstable Discrete Dynamics in Nature and Society locally asymptotically stable around E 2 in the parametric region: where . In presence of fear, that is, when k > 0, we can find α * such that α * < α ′ and E 2 is locally asymptotically stable for α * < α ′ < α < α. Hence, fear in prey promotes the possibility of coexistence in environment. (27). Here, we closely observe the change of system behaviour through local bifurcations with respect to the parameters α and k in which various dynamical consequences, such as stability exchange, occurrence of new steady state, and appearance of limit cycle, may happen. Such local bifurcations can also be analyzed from the positions of nontrivial nullclines.

Local Bifurcations: System
(a) Transcritical Bifurcation. In system (27), coexistence steady state E 2 exchanges its stability with axial equilibrium point E 1 with respect to the bifurcating parameter α; that is, a transcritical bifurcation between E 1 and E 2 occurs. Earlier, it has been already proved that when α * < α < α, the coexistence state E 2 is locally asymptotically stable, but the axial equilibrium E 1 is unstable. Whenever α � α, two equilibria coincide and exchange their stability. For α > α, E 2 becomes unstable (in this case, E 2 is biologically infeasible) and E 1 becomes stable. Here, we prove this bifurcation scenario using Sotomayor's theorem [43]. For α � α, the system exhibits only one axial equilibrium point E 1 . At E 1 , the Jacobian matrix is evaluated as One of the eigenvalues of J(E 1 ; α � α) is zero. e eigenvectors of J(E 1 ; α � α) and [J(E 1 ; α � α)] T corresponding to the zero eigenvalue are, respectively, given by )). Using the expression for V and W and adopting the similar notations from Perko (2013) [43], one can easily calculate the transversality conditions for transcritical bifurcation as follows: us, from Sotomayor's theorem [43], it can be concluded that the model system undergoes a transcritical bifurcation as parameter α passes through the threshold value α.
(b) Hopf Bifurcation. We have already observed that coexistence steady state E 2 (if exists) may lose their stability when the trace of the corresponding Jacobian matrix J E 2 becomes positive. We show that this stability loss occurs via Hopf bifurcation. Taking α as a variable parameter, we assume the characteristic equation of J E 2 is en, the characteristic equation has purely imaginary roots, which are λ 1 � i ����� � D(α * ) and λ 2 � − i ����� � D(α * ). Also, it can be shown that the transversality condition for Hopf bifurcation is fulfilled; that is, (d/dα)(Re(λ))| α * ≠ 0.
For quantitative study of the system and to explore the effect of fear and habitat complexity on population dynamics, we execute some numerical study. In this regard, we therefore consider k and α as variable parameters. Other parametric values are hypothetically fixed at Figures 2(a)-2(f) are the bifurcation diagrams of system (27) for the bifurcation parameter α. In Figures 2(a) and 2(b), we observe that in absence of fear, both prey and predator oscillate for lower degree of habitat complexity. But as α crosses the lower threshold value α [h] � α * � 0.515, a stable behaviour is observed around the coexistence steady state and a Hopf bifurcation occurs at α * � 0.515. Numerically, it can be shown that the first Lyapunov coefficient is negative, and hence the Hopf bifurcation is supercritical and the bifurcating limit cycle is stable. For the intermediate values of α, precisely, if 0.515 < α < 0.9332, the coexistence steady state remains stable; in fact, it is globally asymptotically stable. A transcritical bifurcation occurs at the higher threshold value α [tc] � α � 0.9332, and beyond this value, the predator population goes extinct, and the prey population settles to its carrying capacity (r − d)/a. Now, if we slightly increase the value of fear level to k � 0.1, the Hopf bifurcation point drops to α [h] � α * � 0.165 (Figures 2(c) and  Discrete Dynamics in Nature and Society 2(d)). Moreover, for higher level of fear, Hopf bifurcation disappears; that is, no oscillatory behaviour occurs (Figures 2(e) and 2(f)). erefore, we can conclude that the higher level of fear promotes the coexistence in environment.
In Figure 3, we have plotted the stability regions of different steady states emerging from system (27) in α − k parametric plane. Here, E 2 is unstable and exhibits a stable limit cycle below the Hopf bifurcation curve (redcoloured curve), and it is asymptotically stable in between Hopf bifurcation curve and a transcritical bifurcation curve (blue-coloured curve). Both equilibrium points, E 0 and E 1 , are saddle points below the transcritical bifurcation curve, α � α � 1 − (am/c(r − d)(θ − mh)), which is independent of fear parameter, k. Hence, if the degree of habitat complexity increases to α, whatever the value of k is, predator species always goes extinct, and the prey species survives to its environmental carrying capacity, (r − d)/a. Also, note that for α > α, E 2 is biologically infeasible and E 0 is a unstable saddle, and that suggests that E 1 is globally asymptotically stable. From each region, we take a typical values of (α, k) and plot the phase trajectories to show the representative behaviour of the solutions of system (27) in that region (Figure 4). Phase trajectories are drawn for different initial biomass (100, 20), (130, 10), (80, 60), and (150, 45) to validate that the observed results are global.

Impact of Fear and Habitat Complexity in Population
Biomass.
e prey biomass of the coexistence steady state E 2 (N * , P * ) is independent of fear parameter k, and so change of fear level does not affect the prey biomass. On the other hand, predator biomass in coexistence steady state that satisfies (28) explicitly depends on k. Also, which indicates that predator biomass decreases due to increased level of fear. Biologically, this means that due to the scare of being caught by predator, prey species forage less and that reduces the predation process and consequently predator biomass decreases because of starvation. Now, we differentiate N * and P * with respect to α to obtain the proper behavioural fluctuations of species biomass in respect of habitat complexity. Derivative of N * with respect to α is Hence, the prey population increases as the degree of habitat complexity increases. Now, the differentiation of P * with respect to α gives dP * dα � mP * /rθN * 2 − (a/r) erefore, variation of predator biomass depends on the sign of (mP * /rθN * 2 ) − (a/r). In Figures 2(a)-2(f), we have plotted the change of species biomass with respect to α and observed that prey population increases as α( > α [h] ) increases and attains to its carrying capacity at α � α [tc] . On the other hand, predator biomass initially increases for α > α [h] , but after a    Figure 3. In each figure, all the steady states, specified by small red-coloured filled circles, are the intersection of the prey nullclines (red-coloured curves) and the predator nullclines (green-coloured curves). Trajectories (blue-coloured curves) are drawn for different initial conditions and all the initial conditions are marked by small blue-coloured circles. e black-coloured orbit represents a stable limit cycle that occurs due to supercritical Hopf bifurcation.  Figure 3: Stability regions of different equilibrium points of the model system (27) in α − k parametric plane. For lower values of α and k, both prey and predator populations oscillate and settle to a stable limit cycle. As the degree of habitat complexity crosses the lower threshold curve "Hopf bifurcation curve," the system experiences a supercritical Hopf bifurcation. e coexistence steady state, E 2 , remains stable; in fact, globally asymptotically stability for the values of (α, k) belongs to the light blue shaded region. A transcritical bifurcation occurs as the degree of habitat complexity crosses the critical threshold α � α [tc] � 0.9332. Other set of parameter values is taken as certain value of α, predator biomass decreases steeply and disappears.
We can now identify whether the dynamical characteristics emerging from this particular type of system are applicable for model (7) with a general strictly increasing catch rate, g(N), that shapes the functional response saturating form. Suppose g (N) is unbounded increasing function, then the coexistence steady state occurs under the parametric conditions θ > mh and N * < (r − d)/a, but for bounded and increasing g(N) with sup g (N) � M, we must have another additional requirement (m/c(1 − α)(θ − mh)) < M for existence of coexistence steady state. erefore, one can observe different dynamical behaviour, in particular, existence versus nonexistence of interior equilibrium point (coexistence state) for two similar shaped functional responses. However, if the coexistence steady state occurs, all other qualitative behaviour, such as stability, bifurcations, occurrence of limit cycles, impact of fear, and habitat complexity, will be similar to this particular type of model.
e trivial equilibrium E 0 is always unstable saddle. e axial equilibrium point E 1 is locally asymptotically stable if condition (22) holds, which implies α > α. Note that whenever exactly one coexistence state occurs, E 1 becomes unstable. Now, we investigate the local stability of the coexistence steady state. We can easily get that which is a crucial expression for stability of the nontrivial equilibrium. In particular, a necessary condition for stability is that g ′ (N * i ) > 0, which is always true if exactly one coexistence steady state exists. If two coexistence steady states, say . Hence, bistability occurs between E 2 1 and E 1 . On the other hand, if only one coexistence steady state, namely, E 2 (N * , P * ), exists, then it will be asymptotically stable if Tr(J(E 2 )) < 0. In Figure 5(a), precisely one coexistence steady state E 2 (0.13, 1.11) exists and as Tr(J(E 2 )) < 0, that E 2 is locally asymptotically stable. For the initial size (0.1, 1.0) or (0.3, 1.3), population biomass approaches the coexistence state E 2 (0.13, 1.11) for long time interval. If we slightly change the initial size to (0.5, 1.5) or (2, 2), population biomass exhibits some limit cycle behaviour. is happens due to the occurrence of subcritical Hopf bifurcation, which will be discussed in the next section. In Figure 5(c), multiple coexistence states E 2 1 (0.38, 1.66) and E 2 2 (2.62, 0.53) occur, where E 2 1 is locally asymptotically stable and E 2 2 is unstable saddle. Also, the predator-free steady state E 1 (4.44, 0) is locally asymptotically stable. For initial size (1, 3) or (2, 4), population biomass approaches the coexistence steady state E 2 1 (0.38, 1.66), and for (1, 0.2) or (5, 0.3), predator population goes extinct; that is, (N(t), P(t)) tends to E 1 (4.44, 0) for long time interval. Hence, in this case, bistability between E 1 and E 2 1 occurs. In Figure 5(e), no interior steady state exists; hence, the predator-free steady state E 1 (4.44, 0) is globally asymptotically stable. Different trajectories for initial population biomass (0.1, 3.5), (0.5, 1.5), (3,4), and (4, 2) are drawn in Figure 5. Parameter values are described in Figure 5. (39). Now, we observe that system (39) undergoes various local bifurcations as the parameters k and α cross through critical thresholds.

Local Bifurcations: System
(a) Transcritical Bifurcation. We prove that the axial equilibrium point E 1 exchanges its stability with coexistence steady state E 2 as the bifurcating parameter α crosses the threshold value α � α; that is, a transcritical bifurcation occurs between E 1 and E 2 . For α � α, the Jacobian matrix of system (39) at E 1 is evaluated as follows: Clearly one of the eigenvalues is zero. e eigenvectors of J(E 1 ; α � α) and [J(E 1 ; α � α)] T corresponding to the zero eigenvalue are given by V � v 1 1 and W � 0 1 respec- )). Using the expressions for V and W, we compute the transversality conditions as follows: provided (r − d)/a ≠ 1. All the transversality conditions for transcritical bifurcation are satisfied, and hence by Sotomayor's theorem, we can conclude that system (39) exhibits transcritical bifurcation for α � α. In numerical demonstration, we call this α as α [tc] . We conclude this result in the following theorem.
(c) Hopf Bifurcation. We know that the coexistence steady state E 2 (if exists) lost its stability when the trace of the corresponding Jacobian matrix J E 2 becomes positive. Also, note that for multiple coexistence states, only E 2 1 may lose its stability as E 2 2 is always unstable saddle. We show that this instability occurs via Hopf bifurcation. Let the characteristic equation of J E 2 for a variable parameter α is given by For the suitable choice of other parameter values, we may find a critical value of α, say α * from the expression of T(α), such that T(α * ) � 0 and D(α * ) > 0. en, the characteristic equation has purely imaginary eigenvalues, which are λ 1 � i ����� � D(α * ) and λ 2 � − i ����� � D(α * ). e transversality condition for Hopf bifurcation (d/dα)(Re(λ))| α * ≠ 0 is also satisfied. erefore, system (39) undergoes a Hopf bifurcation, when α � α * .
For better visualization of the bifurcation analysis, we execute some numerical study and demonstrate how the fear level k and the strength of habitat complexity α control the system dynamics. In this regard, we consider the following set of parametric values hypothetically: Figures 6(a)-6(f ) illustrate the qualitative change in system dynamics for the variable parameter α in presence of different levels of fear. In Figures 6(a) and 6(b), we have plotted the change of population biomass and the corresponding bifurcation point for the bifurcation parameter α and the lower fear level, fixed at k � 0.1. We can observe that for sufficiently lower strength of habitat complexity, coexistence steady state is locally asymptotically stable. Further, if we increase the value of α to α * � α [h] � 0.6293, system (39) exhibits limit cycle behaviour and undergoes a Hopf bifurcation around the coexistence state. Also, we calculate the first Lyapunov coefficient numerically which is positive in sign. Hence, the Hopf bifurcation is subcritical, and the bifurcating limit cycle is unstable. Potentially this subcritical case is much more dramatic and dangerous in system biology. Once the strength of habitat complexity crosses the threshold value α [h] , the trajectories of this system must bounce to the nearest attractor (in this case, it must be either a limit cycle with high amplitude or predator-free steady state) which is far away from the coexistence steady state. Whenever this attractor is a limit cycle, it is impossible to get back the stable coexistence steady state by reversing the value of α (see Figure 7(a)). For either case, that is, when the attractor is the predator-free steady state, we have to reverse the α value adequately to get back the coexistence state (see Figure 7(b)). In other words, population biomass is either trapped in a high amplitude limit cycle through the subcritical Hopf bifurcation or exhibits a looping behaviour.
is phenomenon is known as hysteresis. Now at α � α [tc] � 0.7665, the system undergoes a transcritical bifurcation and the predator-free equilibrium becomes locally asymptotically stable. A saddle-node bifurcation also occurs at a higher degree of habitat complexity, α � α [sn] � 0.9, through which the coexistence steady states disappear. Multiple coexistences (E 2 1 and E 2 2 ) and bistability phenomenon occur in between α and α, which is specified by the grey-shaded region. Figure 6(a) demonstrates that prey equilibrium biomass increases continuously due to increasing value of habitat complexity, whereas Figure 6(b) shows that predator equilibrium biomass initially increases, but for higher value of habitat complexity, predator species goes extinct.
To observe the impact of fear in population dynamics, we analyze the model numerically for slightly increasing value of k. In Figures 6(c) and 6(d), we have plotted the bifurcation diagram and change of population biomass for the bifurcation parameter α and a fixed fear level at k � 0.8. In this case, the subcritical Hopf bifurcation occurs at α [h 1 ] � 0.84, which is larger than the previous one (α [h] � 0.6293). Also, another Hopf bifurcation occurs at α [h 2 ] � 0.89 through which the coexistence steady state E 2 1 becomes locally asymptotically stable. Numerically, it can be shown that the first Lyapunov coefficient is negative, and hence this Hopf bifurcation is supercritical. For higher fear level, k � 1.5, both Hopf bifurcations disappear (see Figures 6(e) and 6(f )) and only transcritical and saddle-node bifurcations remained unchanged. Hence, a higher level of fear promotes the coexistence in environment provided the habitat complexity is sufficiently low. Note that α � α [tc] � 1− (m/c(θ − mh))((r − d)/a + a/(r − d)) and α � α [sn] � 1− (2m/c(θ − mh)) are free from fear parameter, and, therefore, bistability phenomenon must occur between one coexistence steady state and the predator-free steady state. Also, for sufficiently higher degree of habitat complexity, coexistence in environment fails to exist and the system settles into the predator-free steady state.
(d) Bogdanov-Takens Bifurcation. For a suitable set of parameter values, system (39) exhibits a saddle-node bifurcation, as well as Hopf bifurcation. For example, we have already proved that for (r − d)/a > 1, two coexistence steady states E 2 1 and E 2 2 appear through saddle-node bifurcation below the threshold value of α � 1 − (2m/c(θ − mh)). Here, E 2 2 is always unstable saddle, and the stability of E 2 1 depends on the sign of Tr(J E 2 1 ). erefore, we can find another critical value of α � α * at which Tr(J E 2 1 ) � 0 and E 2 1 becomes unstable through Hopf bifurcation. Geometrically, Bogdanov-Takens (BT) bifurcation is a point in two-parametric bifurcation plane, where the saddle-node bifurcation curve meets the Hopf bifurcation curve, and at that steady state, the corresponding Jacobian matrix has a zero eigenvalue with algebraic multiplicity two.

(e) Bautin Bifurcation.
e coexistence equilibrium point E 2 1 , if it exists, may change its stability through a supercritical or a subcritical Hopf bifurcation. For a supercritical case, a stable limit cycle occurs through Hopf bifurcation and the first Lyapunov coefficient l 1 > 0, whereas for subcritical case, an unstable limit cycle appears and the first Lyapunov coefficient, l 1 < 0. Now, if l 1 � 0, system (39)    Multiple coexistence occurs in the grey-shaded region. Blue-coloured curve and red-coloured curve indicate the prey/predator biomass of E 1 � (N 1 , P 1 ) � ((r − d)/a, 0) and E 2 � (N * , P * ) (or E 2 1 (N * 1 , P * 1 ) and E 2 2 (N * 2 , P * 2 )), respectively. Solid lines and dashed lines refer to the stable and unstable behaviour of the system respectively. is bifurcation point separates the subcritical and supercritical Hopf bifurcations, and for neighbouring parametric values, the system exhibits one stable and one unstable limit cycles, which collide and vanish through a saddle-node bifurcation of periodic orbits.
In Figure 8, we have plotted the stability regions of various steady states emerging from model (39) in α − k parametric plane fixing all other parameter values at r � 0.6, η � 0.4, d � 0.2, a � 0.09, c � 0.5, h � 0.007, θ � 0.4, m � 0.01} and inspect how the system dynamics depends on the fear level, as well as the habitat complexity. e figure includes a transcritical bifurcation curve (magentacoloured), a saddle-node bifurcation curve (blue-coloured), and a Hopf bifurcation curve (black-coloured). Both transcritical curve, α � α � 1 − (m/c(θ − mh))((r − d)/a+ a/(r − d)), and saddle-node curve, α � α � 1− (2m/c(θ− mh)), are free from fear parameter k. Hence, the system always experiences these two local bifurcations for  Figure 7: Evidence of hysteresis phenomenon. Using the black arrows, we show how the predator biomass is either trapped in a high amplitude limit cycle through the subcritical Hopf bifurcation (a) or exhibits a hysteresis loop (b). In (a), we take k � 0.8, and in (b), k � 1.5.
Other parameter values are taken from Figure 6. Grey-shaded regions refer to the existence of multiple coexistence equilibrium points. Diagrams for prey populations are qualitatively similar.   Figure 8: Stability region of different equilibrium points and two-parametric bifurcation diagram of the model system (39) in α − k parametric plane. One transcritical bifurcation curve (magenta-coloured line), one saddle-node bifurcation curve (blue-coloured line), and one Hopf bifurcation curve (black-coloured curve) split up the whole region into five small subregions (R 1 − R 5 ). In R 1 , exactly one coexistence steady state (E 2 ) occurs which is locally asymptotically stable. In R 2 , that E 2 becomes unstable through subcritical Hopf bifurcation. Multiple coexistence (E 2 1 and E 2 2 ) occurs whenever α and k together belong to R 3 and R 4 . In R 3 , E 2 1 is locally asymptotically stable and E 2 2 is unstable saddle. Further, E 2 1 becomes unstable in R 4 through Hopf bifurcation. Also, in these two regions, predator-free equilibrium point E 1 is locally asymptotically stable. Hence, bistablity phenomenon occurs in R 3 and R 4 . In R 5 , no coexistence occurs and E 1 is globally asymptotically stable. A two-parametric bifurcation Bogdanov-Takens (BT) and generalized Hopf (GH) appear in this diagram.
Other parameter values are the same as in Figure 6. bifurcating parameter α without any restriction on k. One can observe that these bifurcation curves divide the whole parametric region into five small subregions (R 1 − R 5 ). Particularly in R 1 , unique coexistence steady state occurs which is locally asymptotically stable. at coexistence state eventually turns into an unstable state via Hopf bifurcation, if the parametric values (α, k) enter R 2 region. Numerically, one can prove that this Hopf bifurcation is subcritical. Now, if the degree of habitat complexity crosses the transcritical bifurcation curve and lies between α and α, multiple coexistence steady state occurs, of which one is unstable saddle and the other one either is locally asymptotically stable (in R 3 ) or exhibits stable/unstable limit cycle (in R 4 ). Further, in these parametric regions, the predator-free steady state is also locally asymptotically stable. Hence, bistability phenomenon occurs in both R 3 and R 4 . Further, for higher value of α, coexistence states disappear through a saddle-node bifurcation and all the trajectories settle into the predatorfree steady state (in R 5 ). At (α [GH] , k [GH] ) � (0.87, 0.9), two Hopf bifurcations-one is subcritical and the other one is supercritical-collide and disappear through saddle-node bifurcation of periodic orbits, for which a generalized Hopf bifurcation (GH) or Bautin bifurcation occurs. Also, a Bogdanov-Takens (BT) bifurcation occurs at (α [BT] , k [BT] ) � (0.9, 0.44), where the Hopf bifurcation curve touches the saddle-node bifurcation curve. Five phase portraits for five representative values of α and k, taken from each region (R 1 − R 5 ), are plotted in Figures 5(a)-5(e) and display the typical behaviour of the solution trajectories of system (39) in that region.

Impact of Fear and Habitat Complexity in Population
Biomass.
e prey biomass of the coexistence steady state E 2 (N * , P * ) is independent of fear parameter k, and so we only discuss here the effects of fear upon the predator biomass. Differentiating (41) with respect to k, we obtain which indicates that predator biomass decreases due to increased level of fear; that means due to the scare of being caught by predator, the interaction among predator and prey reduces significantly and consequently predator biomass decreases because of starvation. Now, to obtain the variation of species biomass, we differentiate N * and P * with respect to α. Here, Hence, if N * < 1, the prey population increases as the degree of habitat complexity increases, but for N * > 1, prey population decreases. In this particular model, if N * > 1, g(N * ) lies in the decreasing branch of g(N) and hence the corresponding coexistence steady state is unstable saddle, therefore, from stability view point, prey population increases for growing habitat complexity. Now, dP * dα � mP * /rθN * 2 − (a/r) erefore, variation of predator biomass depends on the sign of (mP * /rθN * 2 ) − (a/r). From Figures 6(a)-6(f ), one can observe the variation of species biomass with respect to α. Prey population increases as α increases and attains to its carrying capacity at α � α [tc] . On the other hand, predator biomass initially increases, but after a certain value of α, predator biomass goes extinct.
After some detailed analysis of this particular type of model, we can now recognize the principal characteristics which can be generalized in model (7), where g(N) is a general dome-shaped catch rate. Let, g(N) attain its maximum at N max and g(N max ) � M. e number of coexistence steady states depends on the different restrictions of α, which are described in Table 2, with α � 1− (m/cM(θ − mh)) and α � 1 − (m/c(θ − mh)g((r − d)/a)). Whenever exactly one coexistence steady state E 2 occurs, stability of that E 2 depends on the sign of Tr(J(E 2 )). e system experiences a Hopf bifurcation around E 2 , if Tr(J(E 2 )) becomes zero for some threshold value of α � α * or k � k * . For multiple cases, one coexistence state is either stable or becomes unstable via Hopf bifurcation and the other one is unstable saddle. e system also undergoes a transcritical bifurcation at α � α. If (r − d)/a > N max , a saddle-node bifurcation occurs at α � α. erefore, a bistability and hysteresis phenomenon appear in between α and α for any dome-shaped catch rate.

Discussion
Constructing a predator-prey interaction in modelling framework and identifying the crucial components that determine the ecological complexity are the main ecological research objectives. In particular, functional response is one of the fundamental components that determines the dynamics of a predator-prey interaction. In this study, we have derived a functional response based on timescale separation arguments, which includes the catch rate, g(N), and depending upon the monotonic property of g(N), various functional responses emerge. In particular, if g(N) increases strictly, a saturating functional response occurs; however, if g(N) increases initially and after some certain prey biomass due to group defense g(N) decreases, the resulting functional response will be in a dome-shaped form. Various experimental evidences suggest that increasing structural complexity of habitat reduces the foraging efficiency of predators, and consequently, the attack rate or encounter rate decreases significantly. It is also reasonably true that the presence of predator may affects the physiological behaviour of prey species to such an extent that it can be more effective in reducing prey biomass than direct predation. Incorporating these biological aspects, we have investigated the predator-prey model in terms of two different shaped catch rates and recognized how the choice of two different shaped functional responses affect the qualitative behaviour of a system.
Our proposed model is biologically feasible in the sense that any positive solution initiating in the first quadrant remains positively invariant and uniformly bounded for any time t ≥ 0. We have established the uniform persistent condition for the general model that ensures the existence of the positive equilibrium point of the system. e local asymptotic stability of different equilibrium points has been investigated. By taking a particular form of fear function, we have executed some extensive and comparative analyses of the model dynamics with two different shaped functional responses.
For a model with saturating functional response, only one interior equilibrium point can exist, but for domeshaped functional response, up to two coexistence equilibria can occur and, in general, the complexity of the dynamical behaviour is increased. e geometrical explanations for existence of steady state(s) for both saturating and dome-shaped functional response are discussed in Figure 1. We have implemented bifurcation analysis for the degree of habitat complexity (α) revealing that higher degree of habitat complexity increases the extinction possibility of the predator and the carrying capacity of the prey will be the only stable attractor, no matter how sensitive the prey population is to potential dangers the predator population. Also, interestingly, we cannot obtain a stable predator-free steady state by purely increasing the level of fear.
For the saturating case, lowering the degree of habitat complexity gives rise to stable coexistence and further deceasing may produce oscillatory behaviour (see Figures 2(a)-2(d)). However, this oscillations can be abolished if the level of fear increases sufficiently (see Figures 2(e)-2(f )). But the whole dynamics will collapse, if the coexistence steady state fails to exit, and this dynamical discrepancy may appear for two different but similar saturating-shaped functional responses. erefore, one can observe the structural sensitivity in this type of model systems. Although the structural sensitivity fails to predict the accurate dynamical behaviour, from our investigations, we can anticipate the overall possibilities.
On the other hand, for dome-shaped functional response, coexistence steady state will be locally asymptotically stable, even if the degree of habitat complexity decreases to zero. In that case, the coexistence state is surrounded by an unstable limit cycle which is emerged due to subcritical Hopf bifurcation (see Figures 6(a)-6(d)). Once the system undergoes a subcritical Hopf bifurcation for the bifurcation parameter α, the trajectories of that system must bounce to the nearest attractor which is far away from the coexistence steady state, and this dramatic behaviour can be overcome through a larger alteration of α, or sometimes it is impossible to overcome (hysteresis phenomena). For the intermediate values of α, the system exhibits multiple coexistence steady states and a bistability occurs between predator-free steady state and a coexistence steady state (see Figures 6(a)-6(d)). Now, if we adequately increase the level of fear, oscillatory behaviour will be disappeared which agree with dynamics for saturating case except the bistable phenomenon (see Figures 6(e) and 6(f )). us the present study describes different possible dynamical behaviour of the predator-prey system with a comparison between different shaped function responses. Hence, our study will be helpful to understand the interplay between predator and prey with more ecological impacts.

Data Availability
e data used to support the findings of this study are included within the article.

Conflicts of Interest
e authors declare that they have no conflicts of interest regarding this work.