Hopf Bifurcation , Hopf-Hopf Bifurcation , and Period-Doubling Bifurcation in a Four-Species Food Web

Complex dynamics of a four-species food web with two preys, one middle predator, and one top predator are investigated. Via the method of Jacobian matrix, the stability of coexisting equilibrium for all populations is determined. Based on this equilibrium, three bifurcations, i.e., Hopf bifurcation, Hopf-Hopf bifurcation, and period-doubling bifurcation, are analyzed by center manifold theorem, bifurcation theorem, and numerical simulations. We reveal that, influenced by the three bifurcations, the food web can exhibit very complex dynamical behaviors, including limit cycles, quasiperiodic behaviors, chaotic attractors, route to chaos, perioddoubling cascade in orbits of period 2, 4, and 8 and period 3, 6, and 12, periodic windows, intermittent period, and chaos crisis. However, the complex dynamics may disappear with the extinction of one of the four populations, which may also lead to collapse of the food web. It suggests that the dynamical complexity and food web stability are determined by the food web structure and existing populations.


Introduction
The concept and research method of the food web was firstly proposed as early as in 1930s by Elton, who initiated exploring community structure in ecosystems from the viewpoint of network [1].After that, the food webs have been widely studied in various aspects of characteristics, involving stability, complexity, dynamics, etc. Via analyzing the characteristics of food webs, many research works have played an important role in understanding the properties of natural ecosystems [2][3][4].
Food web systems generally are composed of different combinations of relationships between interacting species, among which predation relationship is a basic one [5,6].Understanding the interactions between prey and predator is a significant topic of longstanding interest in biology.The prey-predator models with various functional responses, such as Holling I, II, III, and IV, Leslie-Gower, Michaelis-Menten, ratio-dependent, Beddington-DeAngelis, and so on, have been extensively discussed [7][8][9][10].Many other ecological mechanisms have been also included and investigated in the prey-predator models.For example, Raw et al. proposed and analyzed a tri-trophic prey-predator model with two predators and one prey exhibiting group defense mechanism and found occurrence of complex dynamical behaviors [7].
Over past decades, mathematical models have been investigated to study the dynamics of three or more interacting species in food webs.Early studies mainly focused on threespecies food web models, demonstrating the complexity of the dynamics for the food webs [6,11].The food web systems may converge not only to stable equilibria and periodic orbits, but also to more complex behaviors such as quasiperiodic oscillations and chaos induced by different bifurcations [11].As studied by Gakkhar and Gupta [6], generalized Hopf bifurcation can exhibit in a food web system with commensalism and Bogdanov-Takens bifurcation point is detected only for a certain value of commensal coefficient.Gupta and Chandra [11] investigated a prey-predator-scavenger model with quadratic harvesting and found that all three species have unpredictable and irregular behavior which may be controlled with the help of an appropriate number of harvesting so that the system can be dynamically balanced.

Mathematical Model and Stability Analysis
In 2002, Fussmann and Gerd [14] have studied the food web complexity and chaotic population dynamics, with proposition of architecture of food web models with two, three, and four to five trophic levels.Their results suggested that natural food webs possess architectural properties that may inherently lower the maybe of chaotic community dynamics [14].Based on their research work, we choose to study the architecture of food web models with three trophic levels including two preys, a middle predator, and a top predator.With the application of Hasting-Powell [16] model and Holling II functional response, a four-species food web model can be described by a dynamical system of four ordinary differential equations: with In the above equations, , , , and  are population densities of the four interacting species and  is time;  1 and  2 represent the intrinsic growth rates of species  and ;  1 and  2 are the carrying capacities;   ( = 1, 2, . . ., 5) represents the conversion rate of prey to predator;  1 and  2 are the death rates of species  and , respectively.  () represents a functional response of the Holling-II type;   and   ( = 1, 2, . . ., 5) are parameters used in the functional response:   describes the effects of capture rate and handling time;   determines how fast per capita feeding rate approaches its saturation rate.From the viewpoint of ecology, parameters  1 ,  2 ,  1 ,  2 ,   ,   ,   and   should all be positive constants.It should also be noticed that population densities are nonnegative.Hence, the state space of the food web system (1) must be strictly restricted in the region  4 + = {(, , , ) ∈  4 :  ≥ 0,  ≥ 0,  ≥ 0,  ≥ 0}.In order to reduce the number of parameters in (1), and to determine which combinations of parameters control the behaviors of the system, system (1) is changed with the application of the following expressions: Using (3) to replace the corresponding variables in (1), then we obtain a simpler form of the food web system: in which The equilibriums of system (4) can be obtained by soling equations / = 0, / = 0, / = 0, and / = 0. Direct calculation yields that system (4) has following six kinds of equilibriums: (i) Trivial equilibrium  0 = (0, 0, 0, 0) always exists, describing the absence of populations.
(iii) When one predator and one prey are absent simultaneously, system (4) degenerates to be a normal Holling II predator-prey system.In such cases, we can find the equilibriums for coexistence of predator and prey, including   ,   ,   , and   .
(iv) With the consideration of reducing one predator population, system (4) degenerates to be a two-trophic-level food web system with two preys and one predator.The equilibriums of such type of food web have been studied by Misra et al. [17].
(v) The equilibriums for the case with absence of one prey population have been investigated by Hsu et al. [18].One can refer to that literature for detailed calculations of equilibriums.
(vi) This research focuses on the positive interior equilibrium  * = ( * ,  * ,  * ,  * ), which represents the coexistence of the four species in the food web.From system (4), it is easy to know that ( * ,  * ,  * ,  * ) is a positive root of the following nonlinear algebraic equations: To determine the local asymptotical stability of the equilibrium  * , the method of Jacobian matrix is utilized.The Jacobian matrix associated to system (4) at any point (, , , ) can be described as follows: The characteristic equation of the Jacobian matrix () at  * can be written as where and   is the element of Jacobian matrix ( * ) at the th row and th volume.Define According to the Routh-Hurwitz criterion, when the following conditions 11) are satisfied simultaneously, the equilibrium  * is locally and asymptotically stable.Otherwise, this equilibrium becomes unstable.
A numerical example is provided here to illustrate the emergence of stable coexistent equilibrium of the food web system.We set parametric conditions of system (4) as follows:  Figure 1 shows the dynamics of the food web when the parameter values are given in (12).In the graph of time series (Figure 1(a)), variables , , ,  oscillate at the beginning with decreasing amplitude and become steady in long-term behaviors.In the -- and -- phase diagrams (Figures 1(b) and 1(c)), it is indicated that the coexistent equilibrium  * is a stable focus.It can be verified that From determinative condition (11), we also know that the coexistent equilibrium  * is locally and asymptotically stable.

Hopf Bifurcation
In this section, Hopf bifurcation theorem and center manifold theorem [19] are used to determine the occurrence conditions of Hopf bifurcation for the four-species food web model.To study the Hopf bifurcation, bifurcation parameter should be chosen at first.Among all parameters of system (4), parameter  is basic and represents the growth rate of prey , controlling one of energy inputs of the food web and therefore we choose it as the bifurcation parameter.
According to the Hopf bifurcation theorem, when system (4) undergoes Hopf bifurcation, the characteristic equation () = 0 should have a pair of conjugate pure imaginary roots, i.e., Substituting the eigenvalues ( 14) into the characteristic equation (8), and separating the real and imaginary parts, then we have Reducing the unknown quantity , the above equations ( 15) and ( 16) turn to be Equation ( 17) determines the relationship of the parameters which exactly leads to the emergence of Hopf bifurcation.
With the other parameter values given, we can calculate the threshold value of the bifurcation parameter .In the following, we denote this threshold value of Hopf bifurcation point as  =  0 .When the value of parameter  changes around the Hopf bifurcation point, the two conjugate imaginary eigenvalues ± will change to be complex eigenvalues, written as  =  +  ω.Substituting  =  + ω into (8), then we obtain the following [12]: Through (19), we find ω can be expressed by .Replacing such expression of ω into (18) results in where () = (3 3  2 +3 2 +4 3 + 1 )/( 3 +4).Differentiating  with respect to  and putting  =  0 , we get If system (4) undergoes Hopf bifurcation around the interior equilibrium  * , we should have (/)| = 0 ̸ = 0. System (4) undergoing Hopf bifurcation at  =  0 should be also under the satisfaction of the last determinative value.We consider two conjugate imaginary eigenvalues  1,2 = ± and two conjugate complex eigenvalues  3,4 =  ±  for the characteristic equation (8).To calculate the determinative value, we further transform the equations of system (4) by applying the following translation: Substituting ( 22) into the equations of system (4) and separating the linear part and nonlinear part.Then system (4) can be rewritten as where in which ((| ) is a fourth-order polynomial functions in regard to independent variable Via introducing an invertible transformation, i.e., where ) , where the elements   1  2 ( 1 ,  2 = 1, ⋅ ⋅ ⋅ , 4) of the matrix P are provided in Appendix, then the normal form of system (23) can be given by where Next, the center manifold   (0, 0, 0, 0) of system (28) at the origin is calculated.On the basis of center manifold theorem, a center manifold   (0, 0, 0, 0) exists and can be described as follows: where ℎ( 1 ,  2 ) and ( 1 ,  2 ) are assumed as By computation, we get Accordingly, the dynamics restricted to the center manifold   (0, 0, 0, 0) can now be considered.Combining with (31), system (28) is simplified as where in which We introduce the following partial derivative notation: in which subscript  1 demonstrates the partial derivative of the first variable; subscript  2 is the partial derivative of the second variable.On account of (34), we have Consequently, the final determination for the system (4) experiencing Hopf bifurcation is calculated by the following [20]: Direct calculation yields Summarizing the above analysis, if condition (38) is satisfied, then system (4) undergoes Hopf bifurcation at equilibrium ( *  0 ,  *  0 ,  *  0 ,  *  0 ).Furthermore, when (/)| = 0 < 0 establishes, if   0 < 0 and 0 <  −  0 << 1, the system experiences supercritical bifurcation, and a stable orbit generates with radius as about √|  0 /  0 |; if   0 > 0 and 0 <  0 −  << 1, subcritical Hopf bifurcation occurs and results in the emergence of an unstable orbit with radius as about √|  0 /  0 |.
For illustrating the occurrence of Hopf bifurcation in the food web system, an example is provided here.First, parameter  is still set as the bifurcation parameter and the other parameter values are given in (12).Under such parametric conditions, the bifurcation critical value of  for occurrence of Hopf bifurcation in system (4) is determined as  0 = 1.1621.Simultaneously, the coexistent equilibrium at the bifurcation critical point is  * (0.2412, 0.1422, 0.8133, 0.1575) and the corresponding eigenvalues are  1,2 = ±0.1642, 3 = −0.1260,and  4 = −0.0557,which show two pairs of imaginary numbers.Moreover, we can derive (/)| = 0 = −0.190021< 0 and the determinative value  0 is calculated as   0 = −5.62119< 0. Consequently, we know that system (4) experiences supercritical Hopf bifurcation at the coexistence equilibrium  * , i.e., when  <  0 , the dynamics of the food web converg to the coexistent equilibrium; when  =  0 , the coexistent equilibrium is weakly attracting; when  >  0 , the dynamics of the food web are attracted to a limit cycle (Figures 2 and 3).
Figure 2 is plotted to show the dynamics of the food web when  = 1.2.By Routh-Hurwitz criterion, we obtain suggesting that the coexistent equilibrium is unstable.In this case, , , , and  will keep periodically oscillating with long-term behavior (Figure 2(a)).In the phase diagram, the periodically oscillating behavior is clearly manifested as a limit cycle (Figures 2(b) and 2(c)).
Figure 3 further demonstrates the Hopf bifurcation diagrams with the change of parameter .As the value of parameter  increases, the system dynamics vary from stable equilibrium to limit cycle, forming bowl-shaped surfaces in Figure 3, which vividly and clearly demonstrate the occurrence of the Hopf bifurcation in the food web system.

Hopf-Hopf Bifurcation
Beyond the Hopf bifurcation, the food web system can still further undergo Hopf-Hopf bifurcation with parameter variation.Hopf-Hopf bifurcation is a complicated codimensiontwo bifurcation, which only appears in four-species food web system or higher-species food web systems.Hopf-Hopf bifurcation can lead to the emergence of a torus solution, appearing like tyre.Hopf-Hopf bifurcation in physical engineering field was fully studied [19,21,22].However, the research on Hopf-Hopf bifurcations of food web system is still not documented in the literature available to the authors.For studying the Hopf-Hopf bifurcation, we choose two parameters  and  as the bifurcation parameters.Depending on the variation of bifurcation parameters, the food web system may undergo Hopf-Hopf bifurcation.At the bifurcation point, equilibrium  * has two pairs of different purely imaginary eigenvalues, denoted as where  1 >  2 > 0.
Then (8) can be rewritten as Balancing the right side of both ( 8) and (41), we easily obtain (42) From the conditions  1 (, ) = 0 and  3 (, ) = 0, we can get the threshold values of the two bifurcation parameters as  =  0 and  =  0 .Moreover, if the conditions  2 ( 0 ,  0 ) > 0,  0 ( 0 ,  0 ) > 0, and  2 2 ( 0 ,  0 ) − 4 0 ( 0 ,  0 ) > 0 are satisfied, the values of  1 and  2 can be also obtained, i.e., In this critical case, the food web system undergoes Hopf-Hopf bifurcation at the equilibrium  * .To demonstrate the Hopf-Hopf bifurcation numerically, we set parameter values of the food web system as follows: Under the parametric conditions of (44), the threshold values of bifurcation parameters at the Hopf-Hopf bifurcation point are determined as  0 = 1.1985,  0 = 1.171.Meanwhile, the coexisting equilibrium is obtained as  * = (0.2746, 0.2474, 0.6053, 0.4371) and its eigenvalues are  1,2 = ±0.4465and  3,4 = ±0.0790,which verify the occurrence of the Hopf-Hopf bifurcation in the four-species food web system.As shown in Figure 4, the dynamics of the system converge to a stable invariant torus in such case, which represents quasiperiodic dynamics of the food web system.
Figure 5 demonstrates the dynamic transition around the Hopf-Hopf bifurcation point with the change of parameter .As the value of parameter  increases, the system dynamics vary from limit cycle to torus (Figure 5(a)).In the Poincaré map, this phenomenon reflects as a dynamic variation from a stable fixed point to a limit cycle, as shown in Figure 5(b).
Figure 6 further exhibits a bifurcation diagram when the food web system experiences Hopf-Hopf bifurcation.The bifurcation diagram is plotted with the successive maximum of  as the parameter  varies in the range [1.08, 1.32].As shown in Figure 6, with the  value increasing around the critical value  0 , the system dynamics experience transition from a stable limit cycle to a torus.After that the food web system comes into the zone of quasiperiodic dynamics, in which a few periodic windows can still take place, such as periodic windows with period 7, 8 orbits (Figures 7(b) and  7(d)).These periodic dynamics are intermittent and appear and disappear abruptly, resulting in alteration of periodic and quasiperiodic dynamics with the increase of bifurcation parameter.At the end of the quasiperiodic zone, the system dynamics turn to be periodic, as demonstrated by period 4 orbit in Figure 7(f).

Period-Doubling Bifurcation
Based on the Hopf bifurcation, we find that the food web can further undergo period-doubling bifurcation.Under the influence of the period-doubling bifurcation, the food web system may experience period-doubling cascade, eventually leading to the occurrence of complex periodic behaviors and chaotic dynamics.
In order to numerically display the period-doubling bifurcation and corresponding system dynamics, we give the parametric conditions as follows:  Here, we also choose parameter  as the bifurcation parameter, and the initial condition is set as (0.5, 0.5, 0.5, 0.5).The numerical simulation results are shown in Figures 7-11.With the change of parameter , the dynamics of all species exhibit period-doubling cascade, leading to dynamical transition from period to chaos.When  < 1.16, the system dynamics exhibit simple periodic oscillation, which appears as a limit cycle in phase portrait (Figure 8(a)) and a single point in Poincaré map (Figure 9(a)).When the value of parameter  increases to  ∈ (1.16, 1.316), the period of the periodic oscillation doubles, i.e., period 2 orbit.This reflects as two cycles in Figure 8  To further determine the chaotic behaviors of the food web, sensitivity analysis is performed via applying two close initial conditions and then observing the system dynamics in long-term behaviors.In Figure 10, red and blue waves show the dynamical behaviors with reference to two close initial conditions, namely, x 0 = 0.5, y 0 = 0.5, z 0 = 0.5,  0 = 0.5 and x 0 = 0. 5001, y 0 = 0.5, z 0 = 0.5,  0 = 0.5.We find that a tiny change (0.0001 in magnitude of x 0 ) in initial conditions can result in a large difference of z and  dynamics as the time progresses.It indicates that the solutions under such parametric conditions are highly sensitive to initial conditions, suggesting the emergence of chaotic dynamics for the food web.
To explicitly demonstrate the dynamic transition between nonchaos and chaos with the increase of bifurcation parameter, bifurcation diagrams are plotted, as shown in Figure 11.Firstly, we find that the period-doubling bifurcation starts a route to chaos (Figure 11(a)).Secondly, when the food web goes into the chaotic zone, periodic windows can repeatedly take place.In the periodic windows, period-doubling cascade also emerges and leads to multiperiodic behaviors (Figure 11(b)).Figure 12 displays the system dynamics in a periodic window with period 3 orbit and period-doubling cascade in orbits of period 3, 6, and 12. Thirdly, when the periodic windows end, the food web may experience chaos crisis, which leads to an abrupt transition from chaotic subsets to a chaotic band (Figure 11(c)).Lastly, the chaos disappears when parameter  increases to higher values (such as  > 1.75), and the food web dynamics return to be periodic.

Discussion and Conclusion
In this research, we have studied the dynamics of a fourspecies food web with Holling II functional response.Three bifurcations, Hopf bifurcation, Hopf-Hopf bifurcation, and period-doubling bifurcation, are found.Various dynamical behaviors are observed, including stable equilibrium, periodic oscillations, quasiperiodic oscillations, period-doubling cascade, periodic window, chaotic interior crisis, and chaotic attractors.With the variation of bifurcation parameter, the food web system may even experience the transitions between equilibrium and period, period and quasiperiod, or nonchaos and chaos.
The complex dynamics occurring in the food web result from the interactions among all populations.If one of the populations is extinct, possibly due to catastrophic shifts of environmental conditions, the system dynamics may experience dramatic change.In order to exhibit the contribution of populations on the dynamics of the food web system, we here choose three cases of typical dynamical behaviors, period (Figure 2), quasiperiod (Figure 4), and chaos (Figure 8(f)) and then remove one population and observe the Mathematical Problems in Engineering dynamics change.Figure 13 shows the situation where the prey population  is extinct.We find that with the absence of one prey, the middle predator or top predator is going to be extinct in the three cases of Figure 13 and that the food web system degenerates a predator-prey system.Similar phenomenon also occurs in the situation when prey  or predator  is extinct.As shown in Figure 14, it is discovered that, with the absence of predator , the prey population  always tends to extinction and that populations  and  keep periodically oscillating in long-term behavior for all three cases.This result demonstrates that the extinction of one population can lead to instability and collapse of the food web.
When we consider that the middle predator  is extinct, the periodic, quasiperiodic, and chaotic behaviors of the food web system all change to be a stable equilibrium (Figure 15).Although the left three populations are still coexisting, the dynamical complexity of the food web is reduced in these cases.The results in Figures 14 and 15 suggest that the omnivorous population can balance the population abundance and      enhance the stability of the food web system.Moreover, Figure 15 also implies that the food web has unpredictable and irregular dynamic behaviors which may be controlled with the help of appropriate number of predation so that the food web can be dynamical equilibrium.

Figure 1 :
Figure 1: (a) Time series of variables , , , ; (b) and (c) dynamics of the food web shown in -- and -- phase diagrams.The parameter values are provided in (12).

Figure 2 :
Figure 2: (a) Time series of variables , , , ; (b) and (c) dynamics of the food web in -- and -- diagrams. = 1.2, and the other parameter values are provided in (12).

Figure 3 :
Figure 3: Hopf bifurcation diagrams corresponding to the variation of parameter .The parametric conditions are given in (12).

Figure 5 :
Figure 5: Dynamic transition around the Hopf-Hopf bifurcation point corresponding to the variation of parameter .The parametric conditions are provided in (44).

Figure 6 :
Figure 6: Bifurcation diagram of the food web system corresponding to the variation of parameter  when Hopf-Hopf bifurcation occurs. = 1.28, and the other parameter values are given in equation (44).

Figure 7 :
Figure 7: Periodic and quasiperiodic dynamics in Poincaré maps, corresponding to the dynamic variation in the bifurcation diagram of Figure 6.

Figure 8 :
Figure 8: Period-doubling cascade and chaotic dynamics of the food web system as shown in the -- phase portraits.The other parameter values are provided in (45).

Figure 9 :
Figure 9: Periodic and chaotic dynamics corresponding to Figure 8 as shown in the Poincaré maps.

Figure 10 :
Figure 10: Sensitivity analysis on the dynamics of food web corresponding to the parametric conditions of Figure 8(f).

Figure 11 :
Figure 11: Bifurcation diagrams of the food web system.(a) Maxima for the population density of middle predator ; (b) and (c) magnifications of (a) when the r value ranges in [1.41, 1.435] and [1.51, 1.57].The parameter values are provided in (45).

Figure 12 :
Figure 12: Period-doubling cascade in orbits of period 3, 6, and 12, as shown in the Poincaré maps, corresponding to the periodic window in Figure 11(b).

Figure 13 :
Figure13: System dynamics in the situation with the absence of prey population , corresponding to the three cases in Figures2, 4, and 8(f).

Figure 14 :
Figure14: System dynamics in the situation with the absence of predator population , corresponding to the three cases in Figures2, 4, and 8(f).

Figure 15 :
Figure15: System dynamics in the situation with the absence of predator population , corresponding to the three cases in Figures2, 4, and 8(f).