Fear effect in a three-species food chain model with generalist predator

Soumitra Pal; Pankaj Kumar Tiwari; Arvind Kumar Misra; Hao Wang; Soumitra Pal; Pankaj Kumar Tiwari; Arvind Kumar Misra; Hao Wang

doi:10.3934/mbe.2024001

Mathematical Biosciences and Engineering

2024, Volume 21, Issue 1: 1-33. doi: 10.3934/mbe.2024001

Previous Article Next Article

Research article

Fear effect in a three-species food chain model with generalist predator

1.
Department of Mathematics, Institute of Science, Banaras Hindu University, Varanasi 221005, India
2.
Department of Basic Science and Humanities, Indian Institute of Information Technology, Bhagalpur 813210, India
3.
Department of Mathematical and Statistical Science, University of Alberta, Edmonton AB T6G 2G1, Canada

Received: 18 October 2023 Revised: 19 November 2023 Accepted: 27 November 2023 Published: 08 December 2023

Within the framework of a food web, the foraging behavior of meso-carnivorous species is influenced by fear responses elicited by higher trophic level species, consequently diminishing the fecundity of these species. In this study, we investigate a three-species food chain model comprising of prey, an intermediate predator, and a top predator. We assume that both the birth rate and intraspecies competition of prey are impacted by fear induced by the intermediate predator. Additionally, the foraging behavior of the intermediate predator is constrained due to the presence of the top predator. It is essential to note that the top predators exhibit a generalist feeding behavior, encompassing food sources beyond the intermediate predators. The study systematically determines all feasible equilibria of the proposed model and conducts a comprehensive stability analysis of these equilibria. The investigation reveals that the system undergoes Hopf bifurcation concerning various model parameters. Notably, when other food sources significantly contribute to the growth of the top predators, the system exhibits stable behavior around the interior equilibrium. Our findings indicate that the dynamic influence of fear plays a robust role in stabilizing the system. Furthermore, a cascading effect within the system, stemming from the fear instigated by top predators, is observed and analyzed. Overall, this research sheds light on the intricate dynamics of fear-induced responses in shaping the stability and behavior of multi-species food web systems, highlighting the profound cascading effects triggered by fear mechanisms in the ecosystem.

Keywords:

Citation: Soumitra Pal, Pankaj Kumar Tiwari, Arvind Kumar Misra, Hao Wang. Fear effect in a three-species food chain model with generalist predator[J]. Mathematical Biosciences and Engineering, 2024, 21(1): 1-33. doi: 10.3934/mbe.2024001

Related Papers:

[1]	Warin Rangubpit, Sunan Kitjaruwankul, Panisak Boonamnaj, Pornthep Sompornpisut, R.B. Pandey . Globular bundles and entangled network of proteins (CorA) by a coarse-grained Monte Carlo simulation. AIMS Biophysics, 2019, 6(2): 68-82. doi: 10.3934/biophy.2019.2.68
[2]	Ken Takahashi, Takayuki Oda, Keiji Naruse . Coarse-grained molecular dynamics simulations of biomolecules. AIMS Biophysics, 2014, 1(1): 1-15. doi: 10.3934/biophy.2014.1.1
[3]	Tomohiro Yanao, Sosuke Sano, Kenichi Yoshikawa . Chiral selection in wrapping, crossover, and braiding of DNA mediated by asymmetric bend-writhe elasticity. AIMS Biophysics, 2015, 2(4): 666-694. doi: 10.3934/biophy.2015.4.666
[4]	Warin Rangubpit, Pornthep Sompornpisut, R.B. Pandey . Thermal-induced unfolding-refolding of a nucleocapsid COVN protein. AIMS Biophysics, 2021, 8(1): 103-110. doi: 10.3934/biophy.2021007
[5]	James C. L. Chow . Recent progress in Monte Carlo simulation on gold nanoparticle radiosensitization. AIMS Biophysics, 2018, 5(4): 231-244. doi: 10.3934/biophy.2018.4.231
[6]	Masoud Jabbary, Hossein Rajabi . High-resolution phantom of a nephron for radiation nephrotoxicity evaluation in biophysical simulations. AIMS Biophysics, 2022, 9(2): 147-160. doi: 10.3934/biophy.2022013
[7]	Nily Dan . Bilayer degradation in reactive environments. AIMS Biophysics, 2017, 4(1): 33-42. doi: 10.3934/biophy.2017.1.33
[8]	Zaynah Sheeraz, James C.L. Chow . Evaluation of dose enhancement with gold nanoparticles in kilovoltage radiotherapy using the new EGS geometry library in Monte Carlo simulation. AIMS Biophysics, 2021, 8(4): 337-345. doi: 10.3934/biophy.2021027
[9]	Carl S. Helrich . Studies of cholesterol structures in phospholipid bilayers. AIMS Biophysics, 2017, 4(3): 415-437. doi: 10.3934/biophy.2017.3.415
[10]	Wei Zhang, Sheng Cao, Jessica L. Martin, Joachim D. Mueller, Louis M. Mansky . Morphology and ultrastructure of retrovirus particles. AIMS Biophysics, 2015, 2(3): 343-369. doi: 10.3934/biophy.2015.3.343

Abstract

1. Introduction

The outbreak of COVID 19 pandemic has caused an urgent need to understand a wide range of issues, both basic and applied, from epidemiology to diverse multi-scale molecular structures of its constitutive elements of the novel coronavirus (CoV-2) and its interaction with surroundings [1]–[4]. The RNA genome of CoV-2 is believed to encode only four major structural, sixteen non-structural and eleven accessory proteins. The four structural proteins including spike (S) protein, nucleocapsid (N) protein, membrane (M) protein, and the envelope (E) protein are key components for viral transmission and reproduction [4]–[9]. The CoVE protein is the smallest with 76 amino acid residues (¹M²Y...⁷⁶V) (for sequence see https://www.uniprot.org/uniprot/P59637) [10] but plays an essential role in replication cycle of the virus such as formation of envelope as an integral membrane protein in viral assembly, budding and pathogenesis [4]–[16].

Very recently, Kuzmin et al. [17] have examined the conformational variability of CoVE protein in presence of a membrane environment using an all-atom MD along with a coarse-grained simulation. Because of ‘prohibitively long time’ with the atomistic simulation in ‘exhaustive sampling of the conformational space’ utility of coarse-grained approach implemented in this study is very well emphasized. They find that CoVE induces curvature in membrane and ‘it has a preferable position’ in already curved membranes with ‘the C-termini localize to the convex region’. Apart from locating specific segments of the protein from the membrane sites, there is no estimate of the overall size of the protein or the underlying dynamics. For example, how the segments of protein perform its stochastic moves and how the overall dynamics (diffusive, non-diffusive, drift etc.) of the protein (i.e. that of its center of mass) sets in as its segments settles in the curved region. How its conformation and dynamics depends on temperature. It may be useful to explore the structure and dynamics of a protein monomer free from membrane as a function of temperature first. Then incorporate the environment such as membrane and other constitutive elements such as solvent to assess its effect on the structure and dynamics. We would like to focus on the folding dynamics of a free CoVE protein, a purely theoretical exercise as an idealized model system, the characteristic of which may be useful in hypothesizing its role in interpreting the laboratory data. Such a basic investigation has not been performed to our knowledge as presented here.

Sequence of CoVE protein and its specific characteristics such as formation of membrane ion channel, virulence intensity, etc. and other proteins are studied in depth [4]–[16]. Three distinct segments of the CoVE protein start with the hydrophilic amino (N) terminal with about 10 residues (¹M–¹⁰G), followed by the hydrophobic transmembrane (TM) segment of about 25–27 residues (¹¹T–³⁷L), and ends with a relatively large hydrophilic carboxyl terminal consisting of the remaining residues (³⁸R–⁷⁶V). The structural details of how these segments interact, compete and cooperate is not fully known. The objective of this article is to investigate the effect of temperature on the structural dynamics by analyzing some of the local and global physical properties by a coarse-grained (CG) Monte Carlo simulation. The model and methods section is described in brief next followed by results and discussion and a conclusion.

2. Model and methods

2.1. Model

The model adopted here is based on an effective and efficient method involving bond-fluctuating mechanism on a discrete lattice with ample degrees of freedom to investigate complex issues in probing long time properties in polymer systems [18] that cannot be accessed by traditional MC or MD in a continuum space. This coarse-grained approach has been further developed over the years [19]–[21] with options to fine-grain and incorporate the specificity of amino acids to model proteins. The relevant detail has been already described to examine structural dynamics of a number of proteins [19]–[21]. In coarse-grained description, we constructed a CG protein model of CoVE containing a chain of 76 amino acid residues [10] in a specific sequence on a discrete cubic lattice with ample degrees of freedom to describe essential motion and fluctuation for each residue in the protein. A residue is modeled using a single CG node with a cubic node of size (2a)³ where a refers to the lattice constant. Consecutive nodes are tethered together by fluctuating covalent bonds of which the bond length between consecutive nodes varies between 2 and √(10) in unit of lattice constant (a). A generalized Lennard-Jones potential (U_ij) is used to describe the interaction between each residue and surrounding residues within a range (r_c) as:

U_{i j} = [| ϵ_{i j} | {(\frac{σ}{γ_{i j}})}^{12} + ϵ_{i j} {(\frac{σ}{γ_{i j}})}^{6}], r_{i j} < r

$U_{i j} = [| ϵ_{i j} | {(\frac{σ}{γ_{i j}})}^{12} + ϵ_{i j} {(\frac{σ}{γ_{i j}})}^{6}], r_{i j} < r$ (1)

where r_ij is the distance between the amino acid residues at site i and j; r_c = √8 and σ = 1 in units of lattice constant. The potential strength e_ij in phenomenological interaction (1) is based on a knowledge-based [22] residue-residue contact matrix which has been developed over many years [23]–[27] from a huge (and rapidly growing) ensemble of protein structures in Protein Data Bank (PDB).

2.2. Methods

The protein chain is placed in a random configuration initially since one can reach any microstate from any other in the vast landscape of statistical ensemble in thermodynamic equilibrium—a fundamental principle of statistical physics. Therefore, starting configuration is irrelevant. Each residue undergoes its stochastic movement with the Metropolis algorithm. The protein configurations are thus generated by movements of each residue with a probability proportional to the Boltzmann factor exp(–ΔE/T) where ΔE is the change in energy between new and old configurations of the underlying residues subject to strictly imposed excluded volume constraints. Monte Carlo time step is defined by attempts to move each residue once. Simulation parameters i.e. length in unit of lattice constant are measured in arbitrary unit and the temperature T in reduced units of the Boltzmann constant.

It is important to point out that connecting the arbitrary units with the real sample size and temperature i.e. Kelvin, and time scale (micro-second) quantitatively is difficult due to phenomenological nature of the residue-residue interaction and the coarse-grained representation. Calibration of the parameters used in simulation requires measurements of some common physical quantities in both computer simulations as well as laboratory experiments which is not feasible. However, based on the approximate size of residue and its representation and response of the physical quantities with respect to changes in temperature, one can estimate the order of magnitude.

All simulations reported here are performed on a 150³ cubic lattice for a range of temperature T = 0.010–0.030 with 100 independent samples each with 10 million time steps; different sample sizes are also used to make sure that the qualitative results are independent of finite size effects. Both local and global physical quantities such as radius of gyration (R_g), root mean square displacement (RMSD) of the center of mass, structure factor (S), contact map, etc. are investigated as a function of reduced temperature. In the following ‘temperature’ refers to ‘reduced temperature’ and should not be confused with the absolute value.

3. Results and discussion

Snapshots of the protein (CoVE) at representative temperatures (T = 0.020, 0.024, 0.030) presented in Figure 1 illustrate the conformational changes as it denatures on raising the temperature. Obviously the protein conforms to a globular (compact) structure at low temperature (T = 0.020) and expands as it denatures (T = 0.030). The segmental globularity in the transmembrane segment of the protein is retained at least in part even at the high temperature (see below). Note that snapshots do not represent the average conformation; it is an illustrative mean that needs to be quantified. Therefore, some quantitative measures of segmental globularity are highly desirable. Average number (N_n) of residues within the range of interaction of each along the contour of the protein may provide some insight into the local restructuring as a function of temperature.

Figure 1. Snapshot of CoVE protein at the end of 10⁷ MCS time. Large black and grey sphere represents the first residue (¹M) and the last residue (⁷⁶V) of the protein, respectively. Smaller golden spheres represent the presence of other residues within the range of interaction excluding the consecutively connected residues. Backbone grey lines are covalent bonds connecting the residues. Snapshot at the left (T = 0.020), center (T = 0.024), right (T = 0.030) are at the representative temperatures to illustrate denaturing of the protein.

DownLoad: Full-Size Img PowerPoint

Figure 2. Residue contact profiles show average number (N_n) of surrounding residues around each residue of the protein within the range of interaction at a range of temperature (T = 0.012–0.024).

DownLoad: Full-Size Img PowerPoint

The contact profiles at representative temperatures are shown in Figure 2. At low temperatures (T = 0.012–0.020, in native phase), the coagulating centers are distributed sporadically particularly in N- and C-terminal regions along the backbone (e.g. ¹M²Y, ¹²L¹³I¹⁴V (in N-terminal), ⁵¹L⁵²V⁵³K, ⁵⁶V⁵⁷Y⁵⁸V⁵⁹Y, ⁷⁴L⁷⁵L⁷⁶V (in C-terminus)) with relatively large N_n. The transmembrane segment (¹⁸L–⁴⁴C), however, has the largest globular cluster with high number N_n of surrounding residues. Clearly, there is no significant change in contact profiles with the temperature in native phase. However, raising the temperature (T = 0.024) further leads to a substantial decrease in segmental globularity in both N- and C-terminals with low contact number N_n while retaining some globularity in TM segment with a relatively large N_n (see Figure 2). Thus, three distinct (N- and C-terminals separated by TM) segments can be identified easily from the evolution of contact profiles with the temperature.

The global structure and dynamics of the protein emerge from the stochastic movement of each residues and their configurational stabilities. Variation of the root mean square displacement of the center of mass (R_c) with the time step (t) at a range of temperature presented in Figure 3 show, how the asymptotic dynamics depends on the temperature. One can identify the type of dynamics by estimating the power-law exponent (ν) by fitting the data with a scaling, R_c ∝ t^ν where ν = ½ represents diffusion. Accordingly, the protein is found to move extremely slowly at low temperatures with ultra-sub-diffusive dynamics (i.e. ν = 0.05 at T = 0.020) and resumes diffusion at higher temperatures (i.e. ν = 0.47 at T = 0.024, see Figure 3). Dynamics do affect the conformation of the protein which can be quantified by analyzing radius of gyration and the structure factor (see below).

Figure 3. Average root mean square displacement (RMSD) of the center of mass (R_c) of the protein with the time step (t) versus temperature in range , T = 0.012–0.030 on a log-log scale. Slopes of the fitted data (solid black lines) in the asymptotic regime at a low (T = 0.012) and a high (T = 0.024) temperature are included with corresponding estimates.

DownLoad: Full-Size Img PowerPoint

How do we know that the protein chain has reached its steady-state equilibrium in 10 million time steps? It is important to check the convergence of physical quantities to make sure that system has reached steady-state within the range of simulation time steps in all such simulations. We have done here as in all our previous studies. For an illustration, variation of the radius of gyration with the time steps for a range of temperature is presented in Figure 4. We see that the radius of gyration has reached its steady-state equilibrium at almost all temperatures. It is not feasible to reach as good steady-state at low temperatures as is the case in most such interacting systems and even in laboratories. However, the probability of reaching the same asymptotic magnitude of R_g at all temperatures (T below 0.020) in sub-native phase (i.e. T = 0.018, 0.016 in Figure 4) is extremely low. Since it is not feasible to extend this simulation for such a long asymptotic time regime, it is worth issuing a caution to very cautious readers. The data on R_g in this sub-native phase (T = 0.010–0.020) in the following Figures should be therefore considered with caution.

Figure 4. Variation of the average radius of gyration (R_g) of the protein with the time step (t) for a temperature range, T = 0.016–0.030.

DownLoad: Full-Size Img PowerPoint

The radius of gyration (R_g) of the protein is a measure of its average size. Variation of the average (over configurations from 5 million to 10 million time steps each with 100 independent samples) radius of gyration (R_g) with the temperature (T) is presented in Figure 5. In contrast to denature phase (T = 0.020–0.030) where the radius of gyration (R_g) increases on raising the temperature, it decays on increasing the temperature in the native phase (T = 0.010–0.020). Why the protein contract on raising the temperature in sub-native phase? At very low temperatures, segmental interactions dominate over thermal noise. On raising the temperature, the frozen or frustrated trapped segments begin to explore and equilibrate better and therefore result into a more compact conformation. This process continues until a certain temperature beyond which thermal noise take over the segmental interactions with the temperature and protein begins to expand. Such an opposite thermal response in native and denature phases is also observed in other membrane proteins in recent years [20]. Based on a crude calibration of the measurements of the physical quantities with the same coarse-grained approach and an all-atom MD simulation [20], we believe that the range of temperature T = 0.020–0.030 may lie in the range of T = 300–400 K roughly. In order to claim such thermal response, a universal or sub-universal characteristics for the membrane proteins, more studies on many different membrane proteins are needed. The thermal response of the nucleocapsid protein (CoVN) of COVID 19 is also non-monotonic [21] but very different from that of CoVE. We hope this results will stimulate further investigations in this direction.

It is worth pointing out in passing that the validation of this model is based on the response of protein at extreme temperature where residue-residue interaction becomes irrelevant and conformation of a protein chain behaves similar to that of an athermal polymer chain. That is the case here as in all our previous work based on this coarse-grained model. Conformational dynamics of a range of proteins with various sizes are examined using this coarse-grained approach and observed different thermal response of proteins with comparable size. For example, proteins H3.1 and H2AX [19] are about the same size (136 and 140 residues) and exhibit very different thermal response to temperature: on raising the temperature, R_g of H3.1 increases continuously before saturation while R_g of H2AX increases to a maximum value followed by decay (a non-monotonic response different from the one observed here for CoVE). Thus, the specificity of residues and its distribution (i.e. sequence) in a protein pays a critical role in their collective response. Unfortunately, we do not have a data on a protein of the COVID related protein with the size comparable to CoVE to compare here. However, based on the analysis on many proteins (some cited in this study) with this approach, the conclusion on the thermal response of CoVE is reliable. A reader can always check it out by carry out simulation independently on an idealized model like the one presented here.

Figure 5. Average radius of gyration (R_g) of CoVE versus temperature.

DownLoad: Full-Size Img PowerPoint

Analysis of the structure factor can also provide insight into the global conformation of the protein. The structure factor (S(q)) is defined as:

S (q) = {⟨ \frac{1}{N} {| \sum_{j = 1}^{N} e^{- i \vec{q} \cdot r_{j}} |}^{2} ⟩}_{| \vec{q} |}

$S (q) = {⟨ \frac{1}{N} {| \sum_{j = 1}^{N} e^{- i \vec{q} \cdot r_{j}} |}^{2} ⟩}_{| \vec{q} |}$ (2)

where r_j is the position of each residue. The magnitude of the wave vector (q) of wavelength λ is given by |q| = 2π/λ. Using a power-law scaling for the structure factor with the wave vector, S(q) ∝ q^–1/γ, one can evaluate the exponent γ. Variations of the structure factor with the wavelength (comparable to radius of gyration) is presented in Figure 5 for representative temperatures.

In general, the radius of gyration (R_g) of a polymer chain shows a power-law scaling with number (N) of monomer (residues), R_g ∝ N^x with a well-defined exponent x. Since a protein is a hetero-polymer chain of residues, one may use the same scaling law here for the protein with the radius of gyration comparable to wavelength (R_g ~ λ). Then the power-law exponent x = γ and one may be able to estimate the globularity index g which is the effective dimension (D) for the spread of the protein (g ≅ D). Since, N ∝ R_g^D, D = 1/γ. Higher magnitude of D represents more compact structure with higher degree of globularity. Note that this scaling analysis is an estimate and should not be taken as a precise measurement. Obviously D cannot be higher than 3, however high value indicates more compact conformations. Clearly, the protein remains globular (three dimensional) at low temperatures in native phase and it denatures into a ramified random coil (D ~ 2) structures at high temperatures (see Figure 6). It is important to point out that the such a vast changes in conformation of the protein in sub-native i.e. low temperature (T = 0.010–0.020) to native-to-denature phase i.e. at higher temperature regime (T = 0.020–0.030) may pose the possibility of COVE to be an intrinsically disordered protein at least in some range of temperature and time scales. For example, the segmental interactions in the native phase may be somewhat similar to interaction of a protein with its underlying environment e.g. a specific surface [28] even in the absence of explicit external surface.

Figure 6. Structure factor (S) of the protein versus wave length (lambda) at representative low and high temperatures. Data points are selected with wavelength comparable to radius of gyration of the protein at corresponding temperatures with the estimates of slopes. The number along the axes are the range of the axes (starting and end points).

DownLoad: Full-Size Img PowerPoint

4. Conclusions

Monte Carlo simulations are performed with a coarse-grained model of envelope membrane protein CoVE of the human coronavirus to assess its structural variability as a function of temperature. The evolution of the protein's conformation as a function of temperature is clearly seen in visualization via snapshots and animations. Segmental organizing of residues is quantified by examining local physical quantities such as average contact and mobility profiles as a function of temperature. From the variations in contact profiles as a function of temperature, it is feasible to identify three distinct segments (N-terminal and C-terminal separated by a transmembrane segment) which is consistent with the previous studies on CoVE [4],[5],[11]–[16]. The radius of gyration exhibits a non-monotonic thermal response: the protein CoVE contracts on heating in native phase and expands continuously on raising the temperature before reaching a saturation in denatured phase. Such a unique non-monotonic thermal response has been recently observed in other membrane proteins (unrelated to coronavirus). This may lead to speculate that thermal response may be a unique characteristic to identify different classes of proteins—obviously more studies are needed to confirm this observation. Thermal response of CoVE is different from another structural protein CoVN (a nucleocapsid protein) of the coronavirus [21]. Variation of the structure factor with the wave vector is studied in detail for a wide range of temperature. Scaling analysis of the structure factor provides a quantitative estimate of the overall globularity with an effective dimension (D). CoVE is found to remain globular (D ~ 3) in its native phase and unfolds continuously to a random coil (D ~ 2) conformation in denature phase on raising the temperature.

References

[1]	A. J. Lotka, Elements of Physical Biology, Williams and Wilkins Company, Baltimore, 1925.
[2]	V. Volterra, Variations and fluctuations of the number of individuals in animal species living together, Anim. Ecol., (1926), 409–448.
[3]	P. Cong, M. Fan, X. Zou, Dynamic of three-species food chain model with fear effect, Commun. Nonlinear Sci. Numer. Simul., 99 (2021), 105809. https://doi.org/10.1016/j.cnsns.2021.105809 doi: 10.1016/j.cnsns.2021.105809
[4]	A. Erbach, F. Lutscher, G. Seo, Bistability and limit cycles in generalist predator-prey dynamics, Ecol. Complex., 14 (2013), 48–55. https://doi.org/10.1016/j.ecocom.2013.02.005 doi: 10.1016/j.ecocom.2013.02.005
[5]	C. Magal, C. G. Cosner, S. Ruan, J. Casas, Control of invasive hosts by generalist parasitoids Math. Med. Biol., 25 (2008), 1–20. https://doi.org/10.1093/imammb/dqm011 doi: 10.1093/imammb/dqm011
[6]	M. Verma, A. K. Misra, Modeling the effect of prey refuge on a ratio-dependent predator–prey system with the Allee effect, Bull. Math. Biol., 80 (2018), 626–656. https://doi.org/10.1007/s11538-018-0394-6 doi: 10.1007/s11538-018-0394-6
[7]	X. Wang, L. Zanette, X. Zou, Modelling the fear effect in predator-prey interactions, J. Math. Biol., 73 (2016), 1179–1204. https://doi.org/10.1007/s00285-016-0989-1 doi: 10.1007/s00285-016-0989-1
[8]	Y. Xu, A. L. Krause, R. A. Van Gorder, Generalist predator dynamics under kolmogorov versus non-Kolmogorov models, J. Theor. Biol., 486 (2020), 110060. https://doi.org/10.1016/j.jtbi.2019.110060 doi: 10.1016/j.jtbi.2019.110060
[9]	M. W. Sabelis, Predatory arthropods, in Natural Enemies: the Population Biology of Predators, Parasites and Diseases, (1992), 225–264.
[10]	M. G. Solomon, J. V. Cross, J. D. Fitzgerald, C. A. M. Campbell, R. L. Jolly, R. W. Olszak, et al., Biocontrol of pests of apples and pears in northern and central Europe-3. Predators, Biocontrol Sci. Technol., 10 (2000), 91–128. https://doi.org/10.1080/09583150029260 doi: 10.1080/09583150029260
[11]	H. Triltsch, Gut contents in field sampled adults of Coccinella septem-punctata (Col.: Coccinellidae), BioControl, 42 (1997), 125–131. https://doi.org/10.1007/BF02769889 doi: 10.1007/BF02769889
[12]	J. P. Harmon, A. R. Ives, J. E. Losey, A. C. Olson, K. S. Rauwald, Colemegilla maculata (Coleoptera: Coccinellidae) predation on pea aphids promoted by proximity to dandelions, Oecologia, 125 (2000), 543–548. https://doi.org/10.1007/s004420000476 doi: 10.1007/s004420000476
[13]	C. Roger, D. Coderre, G. Boivin, Differential prey utilization by the generalist predator Coleomegilla maculata lengi according to prey size and species, Entomol. Exp. Appl., 94 (2000), 3–13. https://doi.org/10.1046/j.1570-7458.2000.00598.x doi: 10.1046/j.1570-7458.2000.00598.x
[14]	B. Mondal, S. Sarkar, U. Ghosh, Complex dynamics of a generalist predator-prey model with hunting cooperation in predator, Eur. Phys. J. Plus, 137 (2022). https://doi.org/10.1140/epjp/s13360-021-02272-4 doi: 10.1140/epjp/s13360-021-02272-4
[15]	S. Creel, D. Christianson, Relationships between direct predation and risk effects, Trends Ecol. Evol., 23 (2008), 194–201. https://doi.org/10.1016/j.tree.2007.12.004 doi: 10.1016/j.tree.2007.12.004
[16]	B. L. Peckarsky, P. A. Abrams, D. I. Bolnick, L. M. Dill, J. H. Grabowski, B. Luttbeg, et al., Revisiting the classics: considering nonconsumptive effects in textbook examples of predator-prey interactions, Ecology, 89 (2008), 2416–2425.
[17]	J. Winnie, S. Creel, The many effects of carnivores on their prey and their implications for trophic cascades, and ecosystem structure and function, Food Webs, 12 (2017), 88–94. https://doi.org/10.1016/j.fooweb.2016.09.002 doi: 10.1016/j.fooweb.2016.09.002
[18]	L. Y. Zanette, A. F. White, M. C. Allen, M. Clinchy, Perceived predation risk reduces the number of offspring songbirds produce per year, Science, 334 (2011), 1398–1401. https://doi.org/10.1126/science.1210908 doi: 10.1126/science.1210908
[19]	M. J. Sheriff, C. J. Krebs, R. Boonstra, The sensitive hare: sublethal effects of predator stress on reproduction in snowshoe hares, J. Anim. Ecol., 78 (2009), 1249–1258. https://doi.org/10.1111/j.1365-2656.2009.01552.x doi: 10.1111/j.1365-2656.2009.01552.x
[20]	J. P. Suraci, M. Clinchy, L. M. Dill, D. Roberts, L. Y. Zanette, Fear of large carnivores causes a trophic cascade, Nat. Commun., 7 (2016), 10698. https://doi.org/10.1038/ncomms10698 doi: 10.1038/ncomms10698
[21]	W. Cresswell, Predation in bird populations, J. Ornithol., 152 (2011), 251–263.
[22]	E. L. Preisser, D. I. Bolnick, The many faces of fear: comparing the pathways and impacts of non-consumptive predator effects on prey populations, PLoS One, 3 (2008), e2465. https://doi.org/10.1371/journal.pone.0002465 doi: 10.1371/journal.pone.0002465
[23]	K. B. Altendrof, J. W. Laundre, C. A. Lopez Gonzalez, J. S. Brown, Assessing effects of predation risk on foraging behaviour of mule deer, J. Mammal., 82 (2001), 430–439. https://doi.org/10.1644/1545-1542(2001)082<0430:AEOPRO>2.0.CO;2 doi: 10.1644/1545-1542(2001)082<0430:AEOPRO>2.0.CO;2
[24]	S. Creel, D. Christianson, S. Liley, J. A. Winnie, Predation risk affects reproductive physiology and demography of elk, Science, 315 (2007), 960. https://doi.org/10.1126/science.1135918 doi: 10.1126/science.1135918
[25]	S. Creel, J. Winnie Jr, B. Maxwell, K. Hamlin, M. Creel, Elk alter habitat selection as an antipredator response to wolves, Ecology, 86 (2005), 3387–3397. https://doi.org/10.1890/05-0032 doi: 10.1890/05-0032
[26]	F. Hua, K. E. Sieving, R. J. Fletcher Jr, C. A. Wright, Increased perception of predation risk to adults and offspring alters avian reproductive strategy and performance, Behav. Ecol., 25 (2014), 509–519. https://doi.org/10.1093/beheco/aru017 doi: 10.1093/beheco/aru017
[27]	S. Pal, N. Pal, S. Samanta, J. Chattopadhyay, Fear effect in prey and hunting cooperation among predators in a Leslie–Gower model, Math. Biosci. Eng., 16 (2019), 5146–5179. https://doi.org/10.3934/mbe.2019258 doi: 10.3934/mbe.2019258
[28]	S. Pal, A. Gupta, A. K. Misra, B. Dubey, Chaotic dynamics of a stage–structured prey–predator system with hunting cooperation and fear in presence of two discrete delays, J. Biol. Syst., 31 (2023), 611–642. https://doi.org/10.1142/S0218339023500213 doi: 10.1142/S0218339023500213
[29]	S. Pal, P. Panday, N. Pal, A. K. Misra, J. Chattopadhyay, Dynamical behaviors of a constant prey refuge ratio–dependent prey–predator model with Allee and fear effects, Int. J. Biomath., 17 (2023), 2350010. https://doi.org/10.1142/S1793524523500109 doi: 10.1142/S1793524523500109
[30]	S. K. Sasmal, Population dynamics with multiple Allee effects induced by fear factors: a mathematical study on prey-predator interactions, Appl. Math. Model., 64 (2018), 1–14. https://doi.org/10.1016/j.apm.2018.07.021 doi: 10.1016/j.apm.2018.07.021
[31]	P. K. Tiwari, M. Verma, S. Pal, Y. Kang, A. K. Misra, A delay nonautonomous predator–prey model for the effects of fear, refuge and hunting cooperation, J. Biol. Syst., 29 (2020), 927–969. https://doi.org/10.1142/S0218339021500236 doi: 10.1142/S0218339021500236
[32]	R. K. Upadhyay, S. Mishra, Population dynamic consequences of fearful prey in a spatiotemporal predator–prey system, Math. Biosci. Eng., 16 (2018), 338–372. https://doi.org/10.3934/mbe.2019017 doi: 10.3934/mbe.2019017
[33]	J. Wang, Y. Cai, S. Fu, W. Wang, The effect of the fear factor on the dynamics of a predator–prey model incorporating the prey refuge, Chaos, 29 (2019), 083109. https://doi.org/10.1063/1.5111121 doi: 10.1063/1.5111121
[34]	H. Zhang, Y. Cai, S. Fu, W. Wang, Impact of the fear effect in a prey–predator model incorporating a prey refuge, Appl. Math. Comput., 356 (2019), 328–337. https://doi.org/10.1016/j.amc.2019.03.034 doi: 10.1016/j.amc.2019.03.034
[35]	D. Mukherjee, Role of fear in predator-prey system with intraspecific competition, Math. Comput. Simul., 177 (2020), 263–275. https://doi.org/10.1016/j.matcom.2020.04.025 doi: 10.1016/j.matcom.2020.04.025
[36]	B. Mondal, S. Roy, U. Ghosh, P. K. Tiwari, A systematic study of autonomous and nonautonomous predator-prey models for the combined effects of fear, refuge, cooperation and harvesting, Eur. Phys. J. Plus, 137 (2022). https://doi.org/10.1140/epjp/s13360-022-02915-0 doi: 10.1140/epjp/s13360-022-02915-0
[37]	P. Panday, N. Pal, S. Samanta, J. Chattopadhyay, A three species food chain model with fear induced trophic cascade, Int. J. Appl. Comput. Math., 5 (2019), 1–26. https://doi.org/10.1007/s40819-019-0688-x doi: 10.1007/s40819-019-0688-x
[38]	M. C. Allen, M. Clinchy, L. Y. Zanette, Fear of predators in free–living wildlife reduces population growth over generations, PNAS, 119 (2022), e2112404119. https://doi.org/10.1073/pnas.211240411 doi: 10.1073/pnas.211240411
[39]	M. Clinchy, M. J. Sheriff, L. Y. Zanette, Predator-induced stress and the ecology of fear, Funct. Ecol., 27 (2013), 56–65. https://doi.org/10.1111/1365-2435.12007 doi: 10.1111/1365-2435.12007
[40]	M. J. Sheriff, S. D. Peacor, D. Hawlena, M. Thaker, Non–consumptive predator effects on prey population size: A dearth of evidence, J. Anim. Ecol., 89 (2020), 1302–1316. https://doi.org/10.1111/1365-2656.13213 doi: 10.1111/1365-2656.13213
[41]	C. E. Gordon, A. Feit, J. Gruber, M. Lentic, Mesopredator suppression by an apex predator alleviates the risk of predation perceived by small prey, Proc. R. Soc. B, 282 (2015), 20142870. https://doi.org/10.1098/rspb.2014.2870 doi: 10.1098/rspb.2014.2870
[42]	W. J. Ripple, R. L. Beschta, Trophic cascades in Yellowstone: The first 15 years after wolf reintroduction, Biol. Conserv., 145 (2012), 205–213. https://doi.org/10.1016/j.biocon.2011.11.005 doi: 10.1016/j.biocon.2011.11.005
[43]	J. K. Hale, Theory of Functional Differential Equations, Springer, New York, 1977.
[44]	L. Perko, Differential Equation and Dynamical System, Springer, New York, 2001.
[45]	P. Liu, J. Shi, Y. Wang, Bifurcation from a degenerate simple eigenvalue, J. Funct. Anal., 264 (2013), 2269–2299. https://doi.org/10.1016/j.jfa.2013.02.010 doi: 10.1016/j.jfa.2013.02.010
[46]	B. Hassard, N. Kazarinoff, Y. Wan, Theory and Applications of Hopf Bifurcation, CUP Archive Cmbridge, 1981.
[47]	S. M. Blower, M. Dowlatabadi, Sensitivity and uncertainty analysis of complex models of disease transmission: an HIV model, as an example, Int. Stat. Rev., 62 (1994), 229–243. https://doi.org/10.2307/1403510 doi: 10.2307/1403510
[48]	S. Marino, I. B. Hogue, C. J. Ray, D. E. Kirschner, A methodology for performing global uncertainty and sensitivity analysis in systems biology, J. Theor. Biol., 254 (2008), 178–196. https://doi.org/10.1016/j.jtbi.2008.04.011 doi: 10.1016/j.jtbi.2008.04.011
[49]	P. Panday, N. Pal, S. Samanta, P. Tryjanowski, J. Chattopadhyay, Dynamics of a stage–structured predator–prey model: cost and benefit of fear–induced group defense, J. Theor. Biol., 528 (2021), 110846. https://doi.org/10.1016/j.jtbi.2021.110846 doi: 10.1016/j.jtbi.2021.110846

Reader Comments

Your name:*

Email:*
© 2024 the Author(s), licensee AIMS Press. This is an open access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0)