Algae-herbivore interactions with Allee effect and chemical defense

G Model

ECOCOM-549; No. of Pages 15 Ecological Complexity xxx (2015) xxx–xxx

Contents lists available at ScienceDirect

Ecological Complexity journal homepage: www.elsevier.com/locate/ecocom

Original Research Article

Algae-herbivore interactions with Allee effect and chemical defense Joydeb Bhattacharyya, Samares Pal * Department of Mathematics, University of Kalyani, Kalyani 741235, India

A R T I C L E I N F O

A B S T R A C T

Article history: Received 12 December 2014 Received in revised form 6 September 2015 Accepted 7 September 2015 Available online xxx

Macroalgae exhibit a variety of characteristics that provide a degree of protection from herbivores. One characteristic is the production of chemicals that are toxic to herbivores. The toxic effect of macroalgae on herbivorous reef fish is studied by means of a spatiotemporal model of population dynamics with a nonmonotonic toxin-determined functional response of herbivores. It is assumed that the growth rate of macroalgae is mediated by Allee effect. We see that under certain conditions the system is uniformly persistent. Conditions for local stability of the system is obtained with weak and strong Allee effects. We observe that in presence of Allee effect on macroalgae, the system exhibits complex dynamics including Hopf bifurcation and saddle-node bifurcation. The obtained results show that the spatiotemporal system does not exhibit diffusion-driven instability. Computer simulations have been carried out to illustrate different analytical results. ß 2015 Elsevier B.V. All rights reserved.

Keywords: Allee effect Macroalgal toxicity Hopf bifurcation Saddle-node bifurcation

1. Introduction Macroalgae-herbivore interactions are important in all ecosystems. Macroalgae can defend themselves against herbivores through a variety of chemicals that are toxic to herbivores (Appelhans et al., 2010). Toxin-mediated interactions between macroalgae and herbivores play an important role in marine ecology. Several algae species are known to produce chemical compounds that reduce the growth and reproduction of fishes. Dinoflagellates are a major marine phytoplankton group which are frequently found as epiphytes on macroalgae and corals. Ciguatoxin (CTX) originates in dinoflagellate species Gamberdiscus toxicus which reduces the growth and reproduction of herbivorous reef fish (Yasumoto et al., 1987). The adverse effects of CTX on fish embryos are studied by (Edmunds et al., 1999). On tropical coral reefs, herbivorous fish promote coral dominance by suppressing competing macroalgae (Rotjan and Lewis, 2006). Toxic-macroalgae reduce herbivory by satiating the detoxification system of herbivores (De Lara-Isassi et al., 2000). Interactions with toxicmacroalgae and herbivorous fish reduce herbivory by slowing the rate of digestion leading to reduced growth rates of herbivores and proliferation of macroalgae in coral reef ecosystem (Hay, 1997). In population dynamics an Allee effect describes the reduction in per capita population growth rate in a low density population level

(Mistro et al., 2012). Under an Allee effect for large populations, reproduction and survival rates are inversely proportional to increased population density (Gonza´lez-Olivares et al., 2011; Sen et al., 2012). A strong Allee effect introduces a population threshold, and the population must surpass this threshold to grow. In contrast, a population with a weak Allee effect does not have a threshold. An Allee effect may arise from difficulties in finding mates, reproductive facilitation and predation (Holt et al., 2004). Increase in herbivory in coral reef keeps macroalgae population in control. Field observations by (Momo, 1995; Wear et al., 1999) suggest that a threshold concentration of macroalgal population is required for its cell division. This means, under excessive herbivory, macroalgae can experience an Allee effect. We have considered a bi-trophic food chain model where toxicmacroalgae are subject to Allee effect and herbivorous fish are growing on toxic-macroalgae at the second trophic level. A toxindetermined functional response is considered to explore how toxic-macroalgae affect herbivore dynamics (Li et al., 2006). We assumed that herbivorous fish are harvested at a rate proportional to its population density. In the present paper the main emphasis will be put on studying dynamics of the system with an Allee effect. We have studied the model analytically as well as numerically. The proofs are all deferred to the Appendix. 2. The basic model

* Corresponding author. Tel.: þ91 33 25666571; fax: þ91 33 25828282. E-mail addresses: [email protected] (J. Bhattacharyya), [email protected] (S. Pal).

We consider a two compartmental system of differential equations in which the concentrations of toxic-macroalgae and

http://dx.doi.org/10.1016/j.ecocom.2015.09.002 1476-945X/ß 2015 Elsevier B.V. All rights reserved.

Please cite this article in press as: Bhattacharyya, J., Pal, S., Algae-herbivore interactions with Allee effect and chemical defense. Ecol. Complex. (2015), http://dx.doi.org/10.1016/j.ecocom.2015.09.002

G Model

ECOCOM-549; No. of Pages 15

K

u

0

K

(a)

u

Growth rate with strong Allee

Logistic growth rate without Allee effect

0

Growth rate with weak Allee

J. Bhattacharyya, S. Pal / Ecological Complexity xxx (2015) xxx–xxx

2

0

m

K

u

(c)

(b)

Fig. 1. (a) Logistic growth rate. (b) Growth rate with weak Allee effect (K < m  0). (c) Growth rate with strong Allee effect (0 < m < K).

herbivores are respectively u(t) and vðtÞ at time t. The growth rate of toxic-macroalgae in the absence of herbivorous fish is subject to   an Allee effect with per capita growth rate f 1 ðuÞ ¼ rðu  mÞ 1  Ku , where m is the Allee threshold i.e., the minimum biomass concentration of macroalgae required to begin a positive growth in absence of herbivorous fish, r is the intrinsic growth rate of macroalgae and K is the carrying capacity, K < m < K. When m > 0, the macroalgal-growth rate decreases if the population size is below the threshold level m and macroalgae goes to extinction, describing strong Allee effect. If K < m  0, it is said that the macroalgal-population is affected by a weak Allee effect (cf. Fig. 1). The toxin generated by macroalgae decreases the growth rate of herbivorous fish. The consumption rate of herbivorous fish is given   ðuÞ eu by f 2 ðuÞ ¼ f ðuÞ 1  f4G , where f ðuÞ ¼ 1þheu , e is the encounter rate per unit of macroalgae, h is the handling time per unit of macroalgae in the absence of toxin and G is the toxin-adjusted maximal amount of toxic-macroalgae a herbivorous fish can ingest per unit time and so, it is a measure of macroalgal-toxicity level (Feng et al., 2011).

G

Since 0  f 2 ðuÞ  f ðuÞ  1h, it follows that observed that f2(u) is monotonic increasing for 1
1 . 2h

1 < G < 1h. It 4h 1  G < 1h and 2h

is is

In the non-monotone case f2(u) reaches

2G with f2(um) = G (cf. Fig. 2). its maximum at um ¼ eð12GhÞ

Considering the drift of macroalgae on ocean currents and the movement of herbivorous fish, a reaction-diffusion system, under the assumption that macroalgae and herbivorous fish are diffusing according to Fick’s law in a rectangular domain V = [0, L1]  [0, L2]  R2 is as follows:

@u ¼ uf 1 ðuÞ  vf 2 ðuÞ þ d1 r2 u; ðx; y; tÞ 2 V  ð0; 1Þ @t

(1)

@v ¼ avf 2 ðuÞ  ðD þ HÞv þ d2 r2 v; ðx; y; tÞ 2 V  ð0; 1Þ @t with the initial conditions uðx; y; 0Þ  0; vðx; y; 0Þ  0, for all (x, y) 2 V and the zero-flux boundary conditions @@un ¼ @@vn ¼ 0 in @V  (0, 1), where n is the outward unit normal vector of the boundary @V, which is assumed to be smooth. Also, r2 is the Laplacian operator in two-dimensional space, which describes random movement. The zero-flux boundary conditions imply the

G

1/4h

2

2

Functional Response (f )

Functional Response (f )

1/2h

0

um

(a)

0

u

Fig. 2. (a) The functional response of herbivorous fish for monotonically increasing with the asymptote below G.

1 4h


1 2h

(b)

u

2G is unimodal and reaches maximum at um ¼ eð12GhÞ . (b) For

1 2h

 G<

1 h

the functional response is

Please cite this article in press as: Bhattacharyya, J., Pal, S., Algae-herbivore interactions with Allee effect and chemical defense. Ecol. Complex. (2015), http://dx.doi.org/10.1016/j.ecocom.2015.09.002

G Model

ECOCOM-549; No. of Pages 15 J. Bhattacharyya, S. Pal / Ecological Complexity xxx (2015) xxx–xxx

mobility of both macroalgae and herbivorous fish population within a bounded habitat. Here D and a are the death rate and growth efficiency (0 < a < 1) respectively of herbivorous fish, H is the rate of harvesting of herbivorous fish and d1, d2 are coefficients of self-diffusivity that measures the dispersal rate; all of these are positive quantities. 3. Adimensionalization of the system We non-dimensionalize the system (1) so that the scaled system contains a minimal number of parameters. Let us change the variables of the system (1) to non dimensional ones ev D by substituting u ¼ Ku ; v ¼ rK ; t ¼ rKt; m ¼ m K ; e ¼ eK; D ¼ rK ; H ¼ d2 H ea rK ; a ¼ r ; d ¼ d1 : Then 1 < m < 1. Under these substitutions, by dropping the bar’s the system (1) reduces to   @u uv eu ¼ uðu  mÞð1  uÞ  1 (2) þ r2 u 1 þ heu 4Gð1 þ heuÞ @t   @v auv eu ¼ 1  ðD þ HÞv þ dr2 v: 4Gð1 þ heuÞ @t 1 þ heu The initial and boundary conditions are uðx; y; 0Þ  0; vðx; y; 0Þ  0, for all (x, y) 2 V and @@un ¼ @@vn ¼ 0 in @V  (0, 1). The variables are nondimensional and the parameters have been scaled by the operating environment, determined by r and K. 4. Spatially homogeneous solutions Before analyzing the reaction-diffusion model, we investigate the spatially homogeneous system and the properties of its solutions that are useful for the analysis of the reaction-diffusion system. In this case the system (2) of partial differential equations is equivalent to the system of ordinary differential equations:   du uv eu ¼ uðu  mÞð1  uÞ  1 (3)  F1 dt 1 þ heu 4Gð1 þ heuÞ   dv auv eu ¼ 1  ðD þ HÞv  F 2 dt 1 þ heu 4Gð1 þ heuÞ with the initial conditions uð0Þ  0; vð0Þ  0. Obviously the right hand sides of system (3) are continuous smooth functions on R2þ ¼ fðu; vÞ : u; v  0g: n oi R h t

We

have R h t

vðtÞ ¼ vð0Þe

uðtÞ ¼ uð0Þe 0 n o

au

0 1þheu

eu 14Gð1þheuÞ

v eu ðumÞð1uÞ1þheu 14Gð1þheuÞ

i

dt

and

ðDþHÞ dt

. This implies, u(t)  0 and vð0Þ  0 whenever u(0) > 0 and vð0Þ > 0. Therefore, all solutions remain within the first quadrant of the u  v plane starting from an interior point of it. Hence, the interior of the positive octant of R2þ is an invariant region. 4.1. Boundedness and persistence Lemma 4.1.1. For all e > 0, there exists te > 0 such that all the  1 solutions of (3) enter into the set ðu; vÞ 2 R2 : uðtÞ þ vðtÞ < a  hðh  m þ D þ HÞ þ e whenever t  te, where h = max {1, u(0)}. DþH The system will be persistent if 0 < liminf t ! 1 xi ðtÞ  limsupt ! 1 xi ðtÞ < 1, for any initial value xi(0) in R2þ , where xi(t) is the population of the ith organism in the system

3

(Thieme, 2003). Persistence represents convergence on an interior attractor from any positive initial conditions and so it can be regarded as a strong form of coexistence (Ruan, 1993). From a biological point of view, persistence of a system ensures the survival of all the organisms in the long run.  Since limsupt ! 1 uðtÞ þ a1 vðtÞ < hðhmþDþHÞ , there exist posiDþH such that tive numbers M1, M2 with M1 < h and M 2 < hðhmþDþHÞ DþH u(t)  M1 and vðtÞ  M 2 for large values of t. The condition given in the following lemma rules out the possibility of extinction of any organism in the system. Lemma 4.1.2. For large values of t, if either M1 ðM1 1Þ < G < 1h h 4GðDþHÞ að4GeM1 Þ4GheðDþHÞ

holds,

there

exists

M2 1M1

u 1 ; v1 > 0

< m < 1 or with

u1 >

such that u1  uðtÞ  M 1 ; v1  vðtÞ  M 2 .

Under the conditions as stated in Lemma (4.1.2), each solution of the system (3) with positive initial values enters in the compact n o set B ¼ ðu; vÞ 2 R2 : u1  uðtÞ  M1 ; v1  vðtÞ  M 2 and remains in it. The minimum population thresholds of toxic-macroalgae and herbivorous fish are given by u1 and v1 respectively, whereas M1, M2 are the maximum threshold levels of toxic-macroalgae and herbivorous fish in the long run. 4.2. Equilibria and their stability In this section we determine biologically feasible equilibrium solutions of the model and investigate the dependence of their stability on several key parameters. The system (3) possesses the following equilibria: (i) Organism-free equilibrium E0 = (0, 0); (ii) herbivorous fish-free equilibria E1 = (1, 0) and E2 = (m, 0); (iii) interior equilibrium E ¼ ðu ; v Þ, given by the intersections of n o v eu au the nullclines ðu  mÞð1  uÞ ¼ 1þheu and 1þheu 1  4Gð1þheuÞ n o eu ¼ D þ H in the interior of the first quadrant, 1  4Gð1þheuÞ where u* is a positive root of the equation Au2 + Bu + C = 0, A = e{a(4Gh  1)  4Gh2e(D + H)}, B = 4G{a  2he(D + H)} and pffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2Gð2hG 1ÞÞ 2 GðGG Þ , of Au2 + C = 4G(D + H). The roots u1;2 ¼ ef4Ghð1hG Þ1g  Bu + C = 0 are real if G  G ¼ eðDþHÞ a . Also, if maxf0; mg < ui < 1 holds, then ðui ; vi Þ are the two possible interior equilibria of eu ðu mÞð1u Þ

i the system, where vi ¼ i i G ; ði ¼ 1; 2Þ. For ease of    notation, we denote u1 ¼ u and v1 ¼ v .

The equilibria E0 and E1 always exist, whereas E2 exists in the case of strong Allee effect only. Examples of mutual positions of the macroalgae and herbivorous-fish nullclines are presented in Fig. 3(a) and (b) with m = 0.1 and m = 0.5 respectively where r = 1, K = 1.6, e = 0.2, h = 0.3, G = 0.9, a = 0.7, D = 0.1, H = 0.1, showing that E2 does’t exist in the case of weak Allee effect and exists for strong Allee effect. The following lemma gives the conditions for existence of the unique interior equilibrium ðu ; v Þ by ruling out the possibility of existence of ðu2 ; v2 Þ under any circumstances. Lemma 4.2.1. For max {0, m} < u* < 1, if either one of the following two conditions holds, then the interior equilibrium ðu ; v Þ of the system exists uniquely: (a) h < (b) h <

1 1 and G > 4hð1hG  ; G Þ  1 1 and G  G < 4hð1hG   . 2G Þ

Please cite this article in press as: Bhattacharyya, J., Pal, S., Algae-herbivore interactions with Allee effect and chemical defense. Ecol. Complex. (2015), http://dx.doi.org/10.1016/j.ecocom.2015.09.002

G Model

ECOCOM-549; No. of Pages 15 J. Bhattacharyya, S. Pal / Ecological Complexity xxx (2015) xxx–xxx

4 0.4

0.2

0.35 0.15

0.3 0.25

0.1

0.15

v

v

0.2 0.05

0.1 0.05

0

0 −0.05

−0.05 −0.1

−0.1

−0.2

0

0.2

0.4

0.6

0.8

1

1.2

0

0.2

0.4

0.6

u

u

(a)

(b)

0.8

1

1.2

Fig. 3. (a) Mutual position of u-nullclines (red) and v-nullclines (blue) with weak Allee effect (for m = 0.1). E2 and E* do not exist. (b) Mutual position of u-nullclines (red) and v-nullclines (blue) with strong Allee effect (for m = 0.5). E* does not exist, whereas, E0, E1 and E2 exist. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

Under the conditions of Lemma (4.2.1), it follows that: (a) with low macroalgal toxicity level, the interior equilibrium exists if the handling time is less than G1 ; (b) with high macroalgal toxicity, the interior equilibrium exists if the handling time is less than 2G1  . Therefore, macroalgal handling time plays a critical role for the existence of the interior equilibrium. With weak Allee effect, the nullclines in Fig. 4(a) intersect at E* 1 in the first quadrant with h < G1 ¼ 3:5 and G > 4hð1hG  ¼ 0:9115 Þ satisfying Lemma 4.2.1(a), where r = 1, m = 0.2, K = 2, e = 0.2, h = 0.3, G = 1, a = 0.7, D = 0.1, H = 0.1. With strong Allee effect, the nullclines in Fig. 4(b) intersect at E* in the first quadrant with 1 h < 2G1  ¼ 1:75 and 0:2857 ¼ G < G < 4hð1hG  ¼ 0:9115 satisfying Þ Lemma 4.2.1(b), where r = 1, m = 0.5, K = 2, e = 0.2, h = 0.3, G = 0.9, a = 0.7, D = 0.1, H = 0.1. The linearized system of (3) about an equilibrium Eˆ is given by dX ˆ ˆ is the Jacobian matrix of the ¼ JðEÞX, where X ¼ ð u v ÞT and JðEÞ dt

ˆ system (3) evaluated at E.

We analyze the stability of system (3) by using eigenvalue analysis of the Jacobian matrix evaluated at the appropriate equilibrium. At E0 the eigenvalues of the Jacobian matrix of the system (3) are m and D  H. Therefore, for m > 0, all the eigenvalues of the Jacobian matrix J(E0) are negative, whereas for m < 0, the product of the eigenvalues of the Jacobian matrix J(E0) is negative. This gives the following lemma: Lemma 4.2.2. The system (3) is always locally asymptotically stable at E0 with strong Allee effect and is a saddle point at E0 with weak Allee effect and m < 0. Lemma 4.2.3. The system (3) is stable at E1 if G < and a > (1 + he)(D + H).

ea 4ð1þheÞfað1þheÞðDþHÞg

Since G is the maximal ingestion rate of toxic-macroalgae by herbivorous fish, for low values of G, the macroalgae becomes highly toxic and so in this case herbivorous fish will not survive in the system. Also, the conditions for stability of (3) at E1 is independent of m. Thus, the stability criterion is valid for both weak and strong Allee effects. Lemma 4.2.4. The system (3) is always unstable at E2 with strong Allee effect.

0.4 0.18

0.35

0.16

0.3

0.14 0.25

0.12

v

v

0.2

0.1 0.08

0.15

0.06

0.1

0.04 0.05

0.02

0 −0.05 −0.2

0 −0.02 0

0.2

0.4

0.6

0.8

1

0

0.1

0.2

0.3

0.4

0.5

u

u

(a)

(b)

0.6

0.7

0.8

0.9

1

Fig. 4. (a) Mutual position of u-nullclines (red) and v-nullclines (blue) with weak Allee effect (for m = 0.2); E2 does not exist, whereas, E0, E1, E* exist. (b) Mutual position of unullclines (red) and v-nullclines (blue) with strong Allee effect (for m = 0.5); E0, E1, E2 and E* exist. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

Please cite this article in press as: Bhattacharyya, J., Pal, S., Algae-herbivore interactions with Allee effect and chemical defense. Ecol. Complex. (2015), http://dx.doi.org/10.1016/j.ecocom.2015.09.002

G Model

ECOCOM-549; No. of Pages 15 J. Bhattacharyya, S. Pal / Ecological Complexity xxx (2015) xxx–xxx

5

The Jacobian J*  J(E*) of the system (3) evaluated at an interior equilibrium E* is 0 1 v f2ðD þ HÞð1 þ heu Þ  au g G     B ðu  mÞð1  2u Þ þ u ð1  u Þ  e C B C au ð1 þ heu Þ2 C J ¼ B B C v f2ðD þ HÞð1 þ heu Þ  au g @ A 0  2  u ð1 þ heu Þ

Lemma 4.2.5. Assume that the conditions of Lemma 4.2.1 are satisfied. If G 6¼ G* and Tr(J*) < 0, the system (3) is locally asymptotically stable at E*.

Appendix. A comprehensive numerical investigation at the bifurcation point signifies that the Hopf bifurcation is supercritical.

The results of local stability analysis of the steady states of the system (3) are summarized in Table 1.

4.3. Optimal harvesting policy

Lemma 4.2.6. If maxf0; mg < u < 12 and G 6¼ eab, the system (3) undergoes a saddle-node bifurcation at G = G*, where 







fv ð12heu Þþ2v1 u ð1þheu Þg b ¼ a4ð1þheu .  Þfv1 ð1þheu Þhev g

n o Lemma 4.2.7. If max 0; m; 23m < u < 1 and h <

1 , the G *

system (3)

undergoes a Hopf bifurcation when m crosses m , where m ¼ mu þ3u2 2u and mþ2u 1

 Þf2G þeu ð2hG 1Þg

m ¼ ð1u

G ð1þheu Þ

2

.

We investigate the orbital stability of Hopf-bifurcating periodic solution by using Poore’s sufficient condition (Poore, 1976). The supercritical and subcritical nature of Hopf-bifurcating periodic solution is determined by positive and negative sign of the real part of F, where F ¼ al Ful j um us b j bm bs þ 2al Ful j um b j ðJE1 Fr b b þ  Þ mr u p uq p q h i 1 al Fu julk b j ðJE  2iv0 Þ Fur p uq b p bq , the repeated indices within kr

each term imply a sum from 1 to 2, all the derivatives of Fl (l = 1, 2) are h i 1 evaluated at E* with u1 ¼ u; u2 ¼ v and ðJ  Þ denotes the mr

1

element in row m, column r of ðJ  Þ . Also, a = (a1, a2) and b = (b1, b2)T are left and right normalized eigenvectors of J* with respect to the eigenvalues iv0 at m = m* so that a b = 1. The detailed calculations for supercritical and subcritical Hopf bifurcation are given in the

Let us consider h = qE where E is external harvesting effort and q is the catchability coefficient. For optimal harvesting policy, we consider that the present value R of a continuous line stream of R1 revenues is given by R ¼ 0 edt Pðv; E; tÞdt, where P = pqvE  cE is the economic rent (net revenue) at time t, d is the instantaneous annual rate of discount, c is the harvesting cost per unit effort and p is the per unit biomass of the prey species. The problem is to maximize the objective functional, using Pontryagin’s maximal principle, subject to the state R1 equation R ¼ 0 edt ð pqv  cÞEðtÞdt, where E(t) is the control variable subject to the constraints 0  E  Emax which define the control set Ut = [0, Emax]. Here Emax denotes a feasible upper limit of E subject to the infrastructural support available for harvesting. Generally, Emax is expected to be a function of v and t. The Hamiltonian function for the problem is given by H ¼ edt ð pqv  cÞE þ l1 ðtÞ½uðu  mÞð1  uÞ n o h n o i uv eu auv eu  þ l2 ðtÞ 1þheu  ðD þ HÞv , 1  4Gð1þheuÞ 1  4Gð1þheuÞ 1þheu where l1(t) and l2(t) are adjoint variables. Suppose that E is the optimal control, and u and v are the corresponding responses. By the maximum principle, there exists dl dl adjoint variables l1, l2 for t  0, such that dt1 ¼  @@Hu and dt2 ¼  @@vH, where

(

)

(

@H v euv av eauv ¼ l1 ðtÞ 2uð1 þ mÞ  3u2  m  þ  þ l2 ðtÞ 2 3 2 3 @u 2Gð1 2Gð1 ð1 þ heuÞ þ heuÞ ð1 þ heuÞ      þ heuÞ @H u l ðtÞ eu a u eu 1 1 1 ¼ edt pqE   l2 ðtÞ  ðD þ qEÞ 1 þ heu 4Gð1 þ heuÞ 1 þ heu 4Gð1 þ heuÞ @v

)

Table 1 Stability analysis of the system (3) with Allee effect. Equilibrium

Criterion for existence Weak Allee

E0 = (0, 0) E1 = (1, 0)

Always exist Always exist

E2 = (m, 0)

Doesn’t exist (m < 0) 1 G 1 and G >  4hð1  hG Þ 1 u < 1; h <  2G 1 and G  G <  4hð1  hG Þ u < 1; h <







E ¼ ðu ; v Þ

pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi  Gð2hG  1Þ þ GðG  G Þ  ef4Ghð1  hG Þ  1g=2 eu ðu  mÞð1  u Þ  v ¼ G u ¼

Conditions for stability Strong Allee

Weak Allee

Strong Allee

Always exist

ae 4ð1 þ heÞfa  ð1 þ heÞðD þ HÞg and a > (1 + he)(D + H) – Saddle point Tr(J*) < 0 0
Always exist 1 m < u < 1; h <  1G and G >  4hð1  hG Þ 1 m < u < 1; h <  2G 1 and G  G <  4hð1  hG Þ

G 6¼ G* and Tr(J*) < 0

Please cite this article in press as: Bhattacharyya, J., Pal, S., Algae-herbivore interactions with Allee effect and chemical defense. Ecol. Complex. (2015), http://dx.doi.org/10.1016/j.ecocom.2015.09.002

G Model

ECOCOM-549; No. of Pages 15 6

J. Bhattacharyya, S. Pal / Ecological Complexity xxx (2015) xxx–xxx

and so (

dl1 dt dl2 dt

v

l1 ðtÞ 3u2  2uð1 þ mÞ þ m þ

¼

euv

)

( þ l2 ðtÞ

eauv

av

)

 3 3 2 ð1 þ heuÞ ð1 2Gð1  2Gð1 þ heuÞ   þ heuÞ  þ heuÞ u l ðtÞ eu a u eu 1 1 1 edt pqE þ þ l2 ðtÞ ðD þ qEÞ  : 1 þ heu 4Gð1 þ heuÞ 1 þ heu 4Gð1 þ heuÞ

¼

 2

(4)

We now consider an optimal solution in the equilibrium state of the above problem so that we may write

ˆ vˆ Þ, the linearized form of the system spatial model (2) by Eˆ ¼ ðu; about Eˆ can be expressed as

  v eu ðu  mÞð1  u Þ  ¼0  1  1 þ heu 4Gð1 þ heu Þ   au eu  D  qE: ¼  1  1 þ heu 4Gð1 þ heu Þ

ˆ wt ¼ LðwÞ  D0 r2 w þ Jw





where w ¼ ðU; V ÞT , D0 ¼ diag ð0; dÞ, Jˆ ¼

"

ˆ ð1þheuÞ

ˆ 4Gð1þheuÞ

#



edt pqE þ l1 ðtÞðu  mÞð1  u Þ:

¼

where

A1 ¼ 3u2  2u ð1 þ mÞ þ m þ v

 f2ðDþqEÞð1þheu Þau g

v f2ðDþqEÞð1þheu Þau gl2 ðtÞ , au ð1þheu Þ2

,

au ð1þheu Þ2

B1 ¼

2

 A1

d li  A2 B1 li ¼ M i edt ; ði ¼ 1; 2Þ dt

where M1 =  pqB1E and M2 = pqE(d + A1). Taking i = 1, the complete solution

It may be noted that (U, V) are small perturbations of ðu; vÞ about ˆ the equilibrium point E. Let 0 = m0< m1 < m2 < be the eigenvalues of the operator  r 2 on V with zero flux boundary conditions. Then, l is an eigenvalue of L if and only if l is an eigenvalue of the matrix Aij ¼ mi D0 þ Jˆ for some i  0. If all the eigenvalues of the operator L have negative real parts, then the system (2) is asymptotically stable at Eˆ and if there is an eigenvalue of L with positive real part, then Eˆ is unstable. Also, if all the eigenvalues of L

d

A2 = (u*  m)(1  u*) and C2 =  e tpqE. Eliminating l1 and l2 we get

dt

vˆ ˆ 2 ð1þheuÞ

v f2ðD þ qEÞð1 þ heu Þ  au g v f2ðD þ qEÞð1 þ heu Þ  au gl2 ðtÞ ¼ l1 ðtÞ 3u  2u ð1 þ mÞ þ m þ þ 2 au ð1 þ heu Þ au ð1 þ heu Þ2 2

This can also be rewritten as dl1 ¼ A1 l1 þ B1 l2 ðtÞ dt dl2 ¼ A2 l1 ðtÞ þ C 2 dt

2

n o euˆ euˆ uˆ ; aˆ 12 ¼  1þhe ; 1  2Gð1þhe 1  4Gð1þhe ˆ ˆ uÞ uÞ uˆ   n o avˆ euˆ auˆ euˆ ; aˆ 22 ¼ 1þhe 1  2Gð1þhe 1  D  H. 2 2 ˆ uÞ uˆ

ˆ  uÞ ˆ  uð1 aˆ 21 ¼

d li

aˆ 11 aˆ 21

(7)  aˆ 12 , u ¼ uˆ þ U; v ¼ aˆ 22

ˆ vˆ þ V; aˆ 11 ¼ ðuˆ  mÞð1  2uÞþ n o



At E ðu ; v Þ, (4) reduces to

dl1 dt dl2 dt

(5)

of

(5)

is

l1 ðtÞ ¼

M

a1 ea1 t þ a2 ea2 t þ N1 edt , where a1, a2 are arbitrary constants, a1, a2 are roots of the auxiliary equation m2  A1m  A2B1 = 0 and N = d2  A1d  A2B1 6¼ 0. It is clear that l1 is bounded iff a1, a2 < 0 or a1 = 0 = a2. M M Then, edt l1 ¼ N1 and similarly, we obtain edt l2 ¼ N2 . Assuming that the optimal equilibrium does not occur either at E = 0 or E = Emax, i.e. the control constants are not binding, we must have the singular control given by

@H @P dt ¼ 0 ) l2 qv ¼ edt ðpqv  cÞ ¼ e : @E @E

(6)

Therefore, the total user cost of harvest per unit effort (left hand side of (6)) must be equal to the discounted value of the future price at steady state effort level (right hand side (6)). d M Eliminating l2 between edt l2 ¼ N2 and l2qv = e t(pqv  c) we   M2 M 2 qEv get c ¼ p  N qv and P ¼ N . 5. Stability analysis in two spatial dimensions The equilibria E0, E1, E2 and E* for the spatially homogeneous system (3) are still the steady-states for the reaction-diffusion system (2). Denoting the general uniform steady state (USS) of the

have non-positive real parts while some eigenvalues have zero real parts, then the stability of (2) at Eˆ cannot be determined by the linearization process (Henry, 1981; Wang et al., 2011). We start by assuming solutions of the form

  n U 1 l tþiðkx xþky yÞ , where l > 0 is the frequency, ¼ e n2 V ni > 0 represents the amplitude (i = 1, 2) and kx, ky > 0 are wave numbers of the perturbations in time t. Then the system (7) becomes

@U ¼ ðaˆ 11  k2x  k2y ÞU þ aˆ 12 V @t n   o @V ¼ aˆ 21 U þ aˆ 22  k2x þ k2y d V @t

(8)

In the spatial model, the value of l depends on the sum of the square of wave numbers k2x þ k2y . Thus both wave numbers affect the eigenvalues. This reflects the fact that some Fourier modes will vanish in the long-term limit while others will amplify. To simplify the situation we can make use of l being rotational symmetric function on the kxky-plane. We substitute k2 ¼ k2x þ k2y and derive the results for the two-dimensional case from one-dimensional formulation. The corresponding characteristic equation of the   linearized system (8) is given by Det Jˆk  lI ¼ 0, where Jˆk ¼

 aˆ 12 aˆ 11  k2 and I stands for 2  2 identity matrix. 2 aˆ 21 aˆ 22  k d

Please cite this article in press as: Bhattacharyya, J., Pal, S., Algae-herbivore interactions with Allee effect and chemical defense. Ecol. Complex. (2015), http://dx.doi.org/10.1016/j.ecocom.2015.09.002

G Model

ECOCOM-549; No. of Pages 15 J. Bhattacharyya, S. Pal / Ecological Complexity xxx (2015) xxx–xxx

7

Table 2 Stability analysis of the system (2) with Allee effect. Equilibrium

Criterion for existence

E0

Always exist

E1

Always exist

E2

Doesn’t exist

Conditions for stability

Weak Allee

Strong Allee

Weak Allee

Stable ae 4ð1 þ heÞfa  ð1 þ heÞðD þ HÞg and a > (1 + he)(D + H) – aem2 0 mð1  mÞ a> m 0
Always exist

(m < 0) u < 1; h < E*

Strong Allee

m >  k2

1 G

1 G 1 and G >  4hð1  hG Þ 1  m < u < 1; h < 2G 1  and G  G <  4hð1  hG Þ m < u < 1; h <

1 and G >  4hð1  hG Þ 1  u < 1; h < 2G 1  and G  G <  4hð1  hG Þ

G 6¼ G* and either one of the two conditions (i)Tr(J*) < k2 (ii)Tr(J*) > k2 and

TrðJ  Þ  k2


k2

DetðJ  Þ

k2 fTrðJ Þ  k2 g

Table 3 Set of parameter values Original Parameters

Description of Parameters

Default value

Dimension

Reference

r m K e h D G

Intrinsic growth rate of toxic-macroalgae Allee threshold Carrying capacity of toxic-macroalgae Encounter rate per unit of toxic-macroalgae Handling time one unit of macroalgae in the absence of toxin Mortality rate of herbivorous fish Max consumption rate of toxic-macroalgae by herbivorous fish Growth efficiency of herbivorous fish Harvesting rate of herbivorous fish Coefficients of diffusivity

1.5 0.5 7 0.08 0.01 0.1 60 0.65 0.1 0.1

1/time Mass/volume Mass/volume 1/time time 1/time 1/time  1/time 

 Feng et al. (2009)   Yasumoto et al. (1987)  Holt et al. (2004), Yasumoto et al. (1987) Thieme (2003)  

a H d1, d2

The eigenvalues for the organism-free steady state E0, obtained from the characteristic equation of the linearized system (8) are m  k2 and D  H  k2d. Therefore, if m  0 holds, then the diffusive system (2) collapses, wiping out all the populations. Further, if m < 0 and k2 >  m hold, then the diffusive system (2) is locally asymptotically stable. This leads to the following result:

Lemma 5.1. The system (2) is always locally asymptotically stable at E0 with strong Allee effect and is locally asymptotically stable at E0 with weak Allee effect if m >  k2. Lemma 5.2. The stability of the spatially homogeneous system (3) at E1 implies the stability of the diffusive system (2) at E1. The spatially homogeneous system is always unstable at E2 with strong Allee effect. The following lemma gives the conditions

0.4

1.4

0.3 0.2

1.0

Herbivores (v)

Toxic−macroalgae (u)

1.2

0.8 0.6 50

0.4 40 0.2 0 300

30

200

150

10 100

Time (t)

50

0

0

0 0.5 0.4 0.3

Space (x)

50

0.2

40

0.1 0 300

20 250

0.1

30 20 250

200

150

10 100

50

0

0

Space (x)

Time (t)

Fig. 5. Biomass distribution of toxic-macroalgae and herbivorous fish over time and space of the system (2) for parameter values as given in Table 3.

Please cite this article in press as: Bhattacharyya, J., Pal, S., Algae-herbivore interactions with Allee effect and chemical defense. Ecol. Complex. (2015), http://dx.doi.org/10.1016/j.ecocom.2015.09.002

G Model

ECOCOM-549; No. of Pages 15 J. Bhattacharyya, S. Pal / Ecological Complexity xxx (2015) xxx–xxx

8

Fig. 6. Biomass distribution of toxic-macroalgae and herbivorous fish over time and space of the model (2) for m = 0.75 and other parameter values as given in Table 3.

Fig. 7. Biomass distribution of macroalgae and herbivores over time and space of the model (2) for m = 0.8 and other parameter values as given in Table 3.

for which diffusion driven stability occurs at E2 with strong Allee effect: k2 dÞ em2 a Lemma 5.3. If G < 4ð1þhemÞfamðDþHþ ; a > ð1þhemÞðDþHþ m k2 dÞð1þhemÞg pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi and k > mð1  mÞ hold, the diffusive system (2) with strong Allee effect is locally asymptotically stable at E2. *

At E the characteristic equation of the linearized system (8) is

l2 + Pl + Q = 0, where P =  a11 + k2(1 + d), Q = Det(J*)  dk2a11 + k4d and a11 = Tr(J*). Lemma 5.4. The stability of the spatially homogeneous system (3) at E* implies the stability of the diffusive system (2) at E*. Therefore system (2) does not exhibit diffusive instability around E0, E1, E2 and E*. In absence of diffusion-driven instability, Turing Hopf bifurcation does not occur in the system (2) and

consequently Turing pattern does not form (Banerjee and Banerjee, 2012). The following lemma gives the conditions for which diffusion can stabilize the unstable spatially homogeneous system at E*: Lemma 5.5. Diffusion driven stability of the system at E* occurs if either one of the following two conditions holds: (a)0 < Tr(J*) < k2; (b)Tr(J*) > k2 and

TrðJ  Þk2

k2



DetðJ Þ < d < k2 fTrðJ  Þk2 g.

The results of local stability analysis of the steady states of the system (2) are summarized in Table 2. 6. Numerical simulations In this section, we investigate the numerical approach as demonstrated in (Bhattacharyya and Pal, 2013; Chattopadhyay

Please cite this article in press as: Bhattacharyya, J., Pal, S., Algae-herbivore interactions with Allee effect and chemical defense. Ecol. Complex. (2015), http://dx.doi.org/10.1016/j.ecocom.2015.09.002

G Model

ECOCOM-549; No. of Pages 15 J. Bhattacharyya, S. Pal / Ecological Complexity xxx (2015) xxx–xxx

et al., 2013; Chaudhuri et al., 2012), the effect of the various parameters on the qualitative behavior of the system using parameter values given in Table 3 throughout, unless otherwise stated. The performed numerical simulations are carried out in a one-dimensional spatial domain. Considering the non-spatial and spatial system, we observe that under the set of parameter values as given in Table 3, the system is locally asymptotically stable at E* (cf. Fig. 5). Effect of Allee threshold parameter (m): Increasing the Allee threshold value (viz. m = 0.75), other parameter values as given in Table 3, the system (2) becomes oscillatory at E* (cf. Fig. 6). For further increase of m (viz. m = 0.8), the system (2) collapses, wiping out all the populations (cf. Fig. 7). Effect of macroalgal toxicity (G): For high macroalgal toxicity (viz. for G = 0.3 and h = 0.85, other parameter values as given in Table 3), the system stabilizes at E1 (cf. Fig. 8(a)).

9

In this case lowering the harvesting rate of herbivorous fish (viz. H = 0.05) stabilizes the system at E* (cf. Fig. 8(b)). Effect of carrying capacity (K): For low carrying capacity (viz. K = 1), other parameter values as given in Table 3, the system stabilizes at E0. Increasing the carrying capacity (viz. K = 1.5), the system becomes stable at E1. Further increase of carrying capacity (viz. K = 4.5) stabilizes the system at E* (cf. Fig. 9). For high carrying capacity (viz. K = 7.3), the system becomes oscillatory around E* (cf. Fig. 10(b)). Effect of harvesting (H): (i) It is observed that with high threshold value of Allee (viz. m = 0.75), other parameter values as given in Table 3, the system is oscillatory at E* (cf. Figs. 6 and 10(a)). In this case increasing the rate of harvesting of herbivorous fish (viz. H = 0.2) stabilizes the system at E* (cf. Fig. 10(a)). (ii) It is observed that with high carrying capacity (viz. K = 7.3), other parameter values as given in Table 3, the system is oscillatory

0.35 0.35 0.3

0.25

Herbivores (v)

Herbivores (v)

0.3

0.2 0.15

0.25 0.2 0.15

0.1

0.1

0.05

0.05

0

0.1

0.2

0.3

0.4

0.5

0.6

0.7

0.8

0.9

1

0

0.1

0.2

Toxic−macroalgae (u)

0.3

0.4

0.5

0.6

0.7

0.8

0.9

1

Toxic−macroalgae (u)

(a)

(b)

Fig. 8. (a) Phase plane diagram for G = 0.3, h = 0.85 and other parameter values as given in Table 3. The domain is characterized by bistability at E0 and E1 depending upon initial conditions. (b) Phase plane diagram for G = 0.3, h = 0.85, H = 0.05 and other parameter values as given in Table 3. The domain is characterized by bistability at E0 and E* depending upon initial conditions. Nullclines are shown as dashed lines.

E*=(0.8584,0.1063)

I = (0.5, 0.1) 1

0.1 0.09

Herbivores (v)

0.08 0.07 0.06 0.05 0.04 0.03 0.02 0.01 E0 0

0.1 E 2

0.2

0.3

0.4

0.5

0.6

0.7

0.8

0.9

1 E1

Toxic−macroalgae (u) Fig. 9. Phase plane diagram for K = 1 and other parameter values as given in Table 3. The system is LAS at E0 (solid blue). For K = 1.5, other parameter values in Table 3, the system is LAS at E1 (dotted blue). For K = 4.5, other parameter values in Table 3, the system is LAS at E* (dotted black). (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

Please cite this article in press as: Bhattacharyya, J., Pal, S., Algae-herbivore interactions with Allee effect and chemical defense. Ecol. Complex. (2015), http://dx.doi.org/10.1016/j.ecocom.2015.09.002

G Model

ECOCOM-549; No. of Pages 15 J. Bhattacharyya, S. Pal / Ecological Complexity xxx (2015) xxx–xxx Toxic−macroalgae (u)

Toxic−macroalgae (u)

10

0.65 0.6 0.55 0.5 0.45

0

100

200

300

400

500

600

700

800

900

0.8 0.7 0.6 0.5 0.4

1000

0

100

200

300

400

0.21 0.205 0.2 0.195 0.19 0.185 0.18

0

100

200

300

400

500

600

700

800

900

1000

600

700

800

900

1000

Time Herbivores (v)

Herbivores (v)

Time

500

600

700

800

900

0.23 0.22 0.21 0.2 0.19 0.18 0.17

1000

0

100

200

300

400

500

Time

Time

(a)

(b)

Fig. 10. (a) Time series analysis of the system with m = 0.75 and other parameter values as in Table 3. The system is oscillatory around E* (solid). For m = 0.75, H = 0.2 and other parameter values as in Table 3, the system is LAS at E* (dotted). (b) Time series analysis of the system with K = 7.5 and other parameter values as in Table 3. The system is oscillatory around E* (solid). For K = 7.3, H = 0.2 and other parameter values as in Table 3, the system is LAS at E* (dotted).

5

x 10

−3

Det (J*)

4 3 2 1 0

0.1

0.2

G*=0.2857

0.4

0.5

0.6

G**=0.7

0.8

0.9115

G

0

*

Tr (J )

−0.2 −0.4 −0.6

II

I

III

−0.8 0

0.1

0.2

G*=0.2857

0.4

0.5

0.6

G**=0.7

0.8

0.9115

G Fig. 11. Bifurcation diagram with r = 1, m = .1, K = 8, e = .08, h = .3, a = 0.7, D = 0.1, H = 0.1 and G as bifurcation parameter.

at E*. In this case increasing the rate of harvesting of herbivorous fish (viz. H = 0.2) stabilizes the system at E* (cf. Fig. 10(b)). Saddle node bifurcation: For G 2 (0, 0.7), r = 1, m = 0.1, h = 0.3, K = 8, a = 0.7 and other parameter values as given in Table 3, we obtain G* = 0.2857. From Fig. 11 we see that DetðJ  ÞjG < 0:2857 ¼ 0 and TrðJ ÞjG < 0:2857 < 0, implying that the Jacobian matrix J* has one zero eigenvalue and a negative real eigenvalue at G* = 0.2857. Also, for G < 0.2857, E* does not exist and consequently, J(E*) and Tr(J*) cannot be evaluated. For 0.2857 < G < 0.7, we have Det(J*) > 0 and Tr(E*) < 0, implying that the system is stable at E*. Thus, the crossing of G = 0.2857 involves in creation and destruction of the fixed point E* leading to a saddle node bifurcation at G = 0.2857. Fig. 12(a) and (b) gives all possible phase portraits corresponding to the parametric domains I and II in Fig. 11. Hopf bifurcation: (i) For G 2 (0.7, 0.9115), r = 1, m = 0.1, h = 0.3, a = 0.7 and other parameter values as given in Table 3, from Fig. 11 we see that Det(J*) > 0. Also, from Fig. 11 we see that TrðJ ÞjG < 0:7 < 0 and

TrðJ ÞjG > 0:7 > 0. Also, at G** = 0.7 we have Det(J*) > 0 and Tr(J*) = 0. Further, the tangent to the curve Tr(J*) at G = 0.7 is not parallel to G-axis. Therefore the system (3) undergoes a supercritical Hopf bifurcation when G crosses G** = 0.7. (ii) For m 2 (0.2, 1.4), other parameter values as given in Table 3, from Fig. 13 we see that Det(J*) > 0. Also from Fig. 13 we see that TrðJ Þjm < 0:714 < 0 and TrðJ  Þjm > 0:714 > 0. Also, at m* = 0.714 we have Det(J*) > 0 and Tr(J*) = 0. Further, the tangent to the curve Tr(J*) at m = 0.714 is not parallel to m-axis. Therefore the system (3) undergoes a supercritical Hopf bifurcation when m crosses m* = 0.714. Fig. 14(a) gives the phase portrait corresponding to the parametric domains III in Fig. 11 and Fig. 14(b) corresponds to the parametric domains IV and V in Fig. 13. 7. Discussion We have considered a food chain with toxic-macroalgae and herbivorous fish where macroalgal growth is under the influence of

Please cite this article in press as: Bhattacharyya, J., Pal, S., Algae-herbivore interactions with Allee effect and chemical defense. Ecol. Complex. (2015), http://dx.doi.org/10.1016/j.ecocom.2015.09.002

G Model

ECOCOM-549; No. of Pages 15 J. Bhattacharyya, S. Pal / Ecological Complexity xxx (2015) xxx–xxx

11

0.45 0.3

0.4

0.25

Herbivores (v)

Herbivores (v)

0.35 0.3 0.25 0.2 0.15

0.2 0.15 0.1

0.1 0.05 0.05 E

E

0

0

0

0.1

0.2

0.3

0.4

0.5

0.6

0.7

0.8

0.9

1

0

0.1

0.2

Toxic−macroalgae (u)

0.3

0.4

0.5

0.6

0.7

0.8

0.9

1

Toxic−macroalgae (u)

(a)

(b)

Fig. 12. (a) Phase plane diagram of the system with G = 0.23 and other parameter values as in Fig. 11. The domain I in Fig. 11 is characterized by bistability at E0 and E1 depending on initial conditions.(b) Phase plane diagram of the system with G = 0.5 and other parameter values as in Fig. 11. The domain II in Fig. 11 is characterized by bistability at E0 and E* depending upon initial conditions. Nullclines are shown as dashed lines.

Tr (J*)

Tr (J*)

0.04

Tr (J*) > 0

0.02 0

Tr (J ) = 0

*

Tr (J ) < 0 −0.02 −0.04 0.2

0.4

0.6

m*=0.714

0.8

1

1.2

m x 10

−3

Det (J*)

*

Det (J )

4 Stable at E* Tr (J*) < 0 & Det (J*) > 0

3.5

3 0.2

Untable at E* Tr (J*) > 0 & Det (J*) > 0

IV

0.4

V

0.6

m*=0.714

0.8

1

1.2

m

Fig. 13. Bifurcation diagram with r = 1.5, K = 7, e = .08, h = .01, a = 0.65, D = 0.1, H = 0.1, G = 60 and m as bifurcation parameter. Hopf bifurcation occurs when the curve Tr(J*) meets m-axis at m* = 0.714.

Allee effect and the functional response representing the growth of herbivorous fish is toxin-determined. The main objective of this paper is to study the dynamic behaviour of the system in presence of macroalgal toxicity and Allee effect. We have shown that solutions of the system are bounded in the long run. It is observed that limited macroalgal toxicity ensures the survival of the organisms in the system in the long run. Also, it is observed that for low handling time, if the minimum viable macroalgal population density crosses a certain critical value, the system enters into Hopf bifurcation that induces oscillation around the interior equilibrium. The stability as well as the direction of bifurcation near the interior equilibrium is obtained by applying the algorithm due to Poore. Also, under certain conditions, the system undergoes a saddle-node bifurcation when the maximum toxic-macroalgae consumption rate by herbivorous fish reaches some threshold. Finally, we investigate spatiotemporal dynamics with constant rate of diffusion for all the species. Our analysis confirms the non-existence of diffusion-driven instability at all the rest points. Further, from analytical and numerical observations we get the following conclusions:

(i) If the minimum viable macroalgal population density is high, the system becomes oscillatory around the interior equilibrium. In this case, increasing the rate of harvesting of herbivorous fish stabilizes the system at the interior equilibrium. For further increase of Allee threshold, the system collapses, wiping out all the populations. (ii) With high macroalgal toxicity, herbivorous fish wipes out from the system. Lowering the harvesting rate of herbivorous fish stabilizes the system at the interior equilibrium. (iii) With low carrying capacity, all the organisms become extinct from the system. Increase of carrying capacity restores the stability of the system at the interior equilibrium. With high carrying capacity, the system becomes oscillatory around the interior equilibrium. In this case, increasing the rate of harvesting of herbivorous fish stabilizes the system at the interior equilibrium. Throughout the article an attempt is made to search for a suitable way to control the growth of toxic-macroalgae and

Please cite this article in press as: Bhattacharyya, J., Pal, S., Algae-herbivore interactions with Allee effect and chemical defense. Ecol. Complex. (2015), http://dx.doi.org/10.1016/j.ecocom.2015.09.002

G Model

ECOCOM-549; No. of Pages 15 J. Bhattacharyya, S. Pal / Ecological Complexity xxx (2015) xxx–xxx

Toxic−macroalgae (u)

12

0.3

0.2

0.6 0.55 0.5 0.45

0

100

200

300

400

0.15

500

600

700

800

900

1000

600

700

800

900

1000

Time 0.23

0.1

Herbivore (v)

Herbivores (v)

0.25

0.65

0.05

E

0

0.22 0.21 0.2 0.19 0.18

0

0.1

0.2

0.3

0.4

0.5

0.6

0.7

0.8

0.9

1

0

100

200

Toxic−macroalgae (u)

300

400

500 Time

(a)

(b)

Fig. 14. (a) Phase plane diagram of the system with G = 0.85 and other parameter values as in Fig. 11. The domain III in Fig. 11 is characterized by dynamic instability at E* and stability at E0 depending upon initial conditions. Nullclines are shown as dashed lines. (b) Time series analysis of the system with m = 0.5 and other parameter values as in Fig. 13. The domain IV in Fig. 13 is stable at E* (dotted). For m = 0.75 and other parameter values as in Fig. 13, the domain V in Fig. 14 is oscillatory around E* (solid).

herbivorous fish for stable coexistence of all the species in the system. From analytical and numerical observations it is seen that Allee effect, macroalgal toxicity level and handling time play a key role in controlling the stability of the system at different rest points. From numerical simulations, it is seen that harvesting of herbivorous fish helps stabilizing the system at the coexistence equilibrium.

Acknowledgements The paper has improved considerably, both in presentation and scientific content, by the critical comments made by the reviewers. The Authors are thankful to them. The research is supported by Science and Engineering Research Board, Government of India, Grant No. SR/S4/MS:863/13.

Appendix A Proof of Lemma 4.1.1. Rt n o F ðu;v;t Þdt v eu , where F 1 ðu; v; tÞ ¼ ðu  mÞð1  uÞ  1þheu 1  4Gð1þheuÞ . We have uðtÞ ¼ uð0Þe 0 1 We consider the following cases: Case (i): Let 0 < u(0) < 1. If possible, let u(t)  1 for all t  0. Then there exists 0 < t1 < t2 such that u(t1) = 1 and u(t) > 1 for all t 2 (t1, t2). Now, for all t 2 (t1, t2), we have R t2 R t2 Rt R t1 F 1 ðu;v;t Þdt F 1 ðu;v;t Þdt F ðu;v;t Þdt F ðu;v;t Þdt uðtÞ ¼ uð0Þe 0 1 ¼ uð0Þe 0 1 e t1 ¼ uðt 1 Þe t1 . Since u(t) > 1 for all t 2 (t1, t2), it follows that F 1 ðu; v; tÞ < 0, and so u(t) < u(t1) = 1, a contradiction. Therefore, for 0 < u(0) < 1, we have u(t)  1, for all t  0. Case (ii): Let u(0) > 1. Then for u(t)  0 for all t  0, we have Rt F ðu;v;t Þdt < uð0Þ and so there is no equilibrium point in the region fðu; vÞ : u > 1; v  0g. uðtÞ ¼ uð0Þe 0 1 Hence combining the two cases, we can say that for any t  0, 0  u(t)  max {1, u(0)}. Let SðtÞ ¼ uðtÞ þ a1 vðtÞ. Then S0 (t)  h(h  m + D + H)  (D + H)S(t), where h = max {1, u(0)}. Let wðtÞ be a solution of S0 (t) = h(h  m + D + H)  (D + H)S(t) satisfying wð0Þ ¼ Sð0Þ. n o eðDþHÞt þ hðhmþDþHÞ and so for large values of t, we have uðtÞ þ a1 vðtÞ < hðhmþDþHÞ þ e, where e is an Then wðtÞ ¼ Sð0Þ  hðhmþDþHÞ DþH DþH DþH arbitrary positive number. Since u(t)  0 and vðtÞ  0, for all t  0, it follows that 0  uðtÞ þ a1 vðtÞ  hðhmþDþHÞ , for all t  0 and so the system (3) is dissipative. DþH Proof of Lemma 4.1.2. There exists T1 > 0 such that u(t)  M1 and vðtÞ  M 2 , for all t  T1. h i uðtÞ Also, for all t  T1, we have du  uðtÞ uðtÞð1  uðtÞÞ  mð1  M 1 Þ  M 2 þ 4M : Gh dt 1

This implies liminf t ! 1 uðtÞ  u1 , where u1 ¼ M If M1 < 1 and 1M2 < m < 1, then u1 > 0.

4M1 Ghfmð1M1 ÞM 2 g . 4M2 Ghð1M1 Þþ1

1

M ðM 1Þ

< G < 1h, then u1 > 0. If M1  1 and 1 h1 Therefore, there exists 0 < u1  M1 such that u1  u(t)  M1, for all t  T1.

Please cite this article in press as: Bhattacharyya, J., Pal, S., Algae-herbivore interactions with Allee effect and chemical defense. Ecol. Complex. (2015), http://dx.doi.org/10.1016/j.ecocom.2015.09.002

G Model

ECOCOM-549; No. of Pages 15 J. Bhattacharyya, S. Pal / Ecological Complexity xxx (2015) xxx–xxx

For t  T1 we have

dv dt

v

n

au1 1þheu1



1

eM 1 4G



13

o

4GðDþHÞ  D  H > 0 if u1 > að4GeM . 1 Þ4GheðDþHÞ

4GðDþHÞ > 0 for uðtÞ  u1 > að4GeM and so in this case there exists T2 > 0 and 0 < v1  M2 such that Therefore, for all t  T1, dv dt 1 Þ4GheðDþHÞ vðtÞ  v1 ; 8 t  T 2 . M M ðM 1Þ 4GðDþHÞ Therefore, for large values of t, if either 1M2 < m < 1 or 1 h1 < G < 1h holds, there exists u1 ; v1 > 0 with u1 > að4GeM such 1 1 Þ4GheðDþHÞ that u1  uðtÞ  M 1 ; v1  vðtÞ  M 2 . Proof of Lemma 4.2.1. 1 (a) Let h < G1 and G > 4hð1hG  hold. Then A > 0 and C < 0. Þ  2

ð12hG Þ 1 G  < 0 and  ¼  pffiffiffiffiffiffiffiffiffiffiffiffiffi4hð1hG Þ  4hðhGp1Þ ffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 2Gð2hG 1ÞÞþ2 GðGG Þ Bþ B 4AC  u ¼ ¼ . 2A ef4Ghð1hG Þ1g

Also,

so

B2  4AC.

This

implies,

Au2 + Bu + C = 0

has

a

unique

positive

real

root

Since max {0, m} < u* < 1, we have v > 0 and so E ðu ; v Þ exists uniquely. and G  G <

hold. Then A < 0, C < 0 and B2  4AC. pffiffiffiffiffiffiffiffiffiffiffiffiffi pffiffiffiffiffiffiffiffiffiffiffiffiffi B2 4AC B2 4AC This implies Au2 + Bu + C = 0 has two positive real roots u ¼ Bþ 2A and u2 ¼ B 2A .  *   Since max {0, m} < u < 1, it follows that E ¼ ðu ; v Þ exists. pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi Now, u < 1 ) 4eðD þ HÞð1  heÞð1  2GhÞ < B þ 2A < B2  4AC . 1 Since he < 1 and G < 2h , it follows that B + 2A > 0 pffiffiffiffiffiffiffiffiffiffiffiffiffi 2 eu ðu mÞð1u2 Þ < 0. This implies u2  1 ¼  Bþ2Aþ 2AB 4AC > 0 ) u2 > 1 and consequently, v2 ¼ 2 2 G (b) Let h <

1 2G

1 4hð1hG Þ

Therefore, if maxf0; mg < u < 1; h < 2G1  and G  G < Proof of Lemma 4.2.3. The Jacobian matrix of the system (3) at E1 is 0 B JðE1 Þ ¼ B @

m1

1 4hð1hG Þ

hold, then E* exists uniquely.

1

e  4Gð1 þ heÞ

C 4Gð1 þ heÞ C A ea  D  H 1 þ he 4Gð1 þ heÞ2 2

(9)

a

0

a  ea At E1 the eigenvalues of the Jacobian matrix are m  1; 1þhe  D  H. 4Gð1þheÞ2 Since 1 < m < 1, we have m  1 <0. a  ea ea Also, 1þhe 2  D  H < 0 if 0 < G < 4ð1þheÞfað1þheÞðDþHÞg and a > (1 + he)(D + H). 4Gð1þheÞ

Therefore, if G <

ea 4ð1þheÞfað1þheÞðDþHÞg

and a > (1 + he)(D + H) hold, the system (3) is locally asymptotically stable at E1.

Proof of Lemma 4.2.4. E2 exists when 0 < m < 1 and so, we analyze the stability of the system at E2 in case of strong Allee effect only. The Jacobian matrix of the system (3) evaluated at E2 is 0

B mð1  mÞ B JðE2 Þ ¼ B @ 0



4mG þ em2 ð4Gh  1Þ

1

C 2 4Gð1 þ hemÞ C C A ma em2 a  D  H 2 1 þ hem 4Gð1 þ hemÞ

ma  At E2 the eigenvalues of the Jacobian matrix are mð1  mÞ; 1þhem

em2 a 4Gð1þhemÞ2

 D  H.

With strong Allee effect we have m(1  m) > 0 and so the system (3) is always unstable at E2. Proof of Lemma 4.2.5. We consider the following cases: * 1 Case-1: Let maxf0; mg < u < 1; h < G1 and G > 4hð1hG  hold. Then E exists uniquely. Þ  







þeu ð2hG 1Þg then DetðJ Þ ¼ av G f2G > 0 and so E* is locally asymptotically stable if Tr(J*) < 0. 2 e2 u ð1þheu Þ pffiffiffiffiffiffiffiffiffiffiffiffiffi pffiffiffiffiffiffiffiffiffiffiffiffiffi  B2 4AC ðB B2 4AC Þ  2G  If h < 2G1  , then u þ 2C < 0 (since A, B > 0). This implies, u <  2C B ¼ B ¼ eð12hG Þ ) 2G þ eu ð2hG  2AB  * * 1Þ > 0 ) DetðJ Þ > 0 and so E is locally asymptotically stable if Tr(J ) < 0. * * 1 Therefore, if maxf0; mg < u < 1; h < G1 ; G > 4hð1hG  and Tr(J ) < 0, then E is locally asymptotically stable. Þ

If

1 2G


1 , G

and G < G <

hold. Then E* exists uniquely. pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 1 Since h < 2G1  and G < G < 4hð1hG B2  4AC .  , it follows that A < 0 and B > Þ pffiffiffiffiffiffiffiffiffiffiffiffiffi pffiffiffiffiffiffiffiffiffiffiffiffiffi B2 4AC ðB B2 4AC Þ * 2G This implies u þ 2C < 0 ) u < eð12hG  and so Det(J ) > 0. B ¼ 2AB Þ Case-2: Let maxf0; mg < u < 1; h <

1 2G

1 4hð1hG Þ

* * 1 Therefore, if maxf0; mg < u < 1; h < 2G1  ; G < G < 4hð1hG  and Tr(J ) < 0, then E is locally asymptotically stable. Þ Proof of Lemma 4.2.6. 1 2G For max {0, m} < u* < 1 and G < 2h , there is only one interior equilibrium E ðu ; v Þ of the system (3) where u ¼ eð12G  hÞ and

v ¼ eu

 ðu mÞð1u Þ G

; m < u < 1. At G = G*, we have



J ¼

ðu  mÞð1  2u Þ þ u ð1  u Þ 0

G e 0

!

Please cite this article in press as: Bhattacharyya, J., Pal, S., Algae-herbivore interactions with Allee effect and chemical defense. Ecol. Complex. (2015), http://dx.doi.org/10.1016/j.ecocom.2015.09.002

G Model

ECOCOM-549; No. of Pages 15 14

J. Bhattacharyya, S. Pal / Ecological Complexity xxx (2015) xxx–xxx

Moreover, for m < u < 12, the zero eigenvalue of the Jacobian matrix is simple. Let V and W be the eigenvectors corresponding to the zero eigenvalue for J* and J*T respectively. Then we obtain

v1 ÞT

V ¼ ð1 and

1 ÞT

W ¼ ð0









Þþu ð1u Þg where v1 ¼ efðu mÞð12u . G Let us express the system (3) in the form

X˙ ¼ Fðu; vÞ

(10)

where X ¼ ð u

v ÞT and

 1 uv eu B uðu  mÞð1  uÞ  1 þ heu 1  4Gð1 þ heuÞ C C   Fðu; v; GÞ ¼ B @ A auv eu 1  ðD þ HÞv 1 þ heu 4Gð1 þ heuÞ 0

Then W T F G ðu ; v ; G Þ ¼

eau2 v 2 4G"2 ð1þheu Þ 2

D2 Fðu ; v ; G ÞðV; VÞ ¼

> 0 and

!

@ Fðu; v; G Þ @2 Fðu; v; G Þ @2 Fðu; v; G Þ @2 Fðu; v; G Þ þ þ v1 þ v1 v1 2 @u@v @v@u @u @v2

a½2ð1 þ heu Þf2G v1 ð1 þ heu Þ  eu v1  2G hev g þ ev ð1  2heu Þ

¼

 4

2G ð1 þ heu Þ

# ðu;vÞ¼ðu ;v Þ

afv ð1  2heu Þ þ 2v1 u ð1 þ heu Þg 6¼ b ¼ 6 0 if ¼ .    e 4ð1 þ heu Þfv1 ð1 þ heu Þ  hev g a G

Therefore, if maxf0; mg < u < 12 and G 6¼ eab, the system (3) undergoes a saddle-node bifurcation at G = G*. Proof of Lemma 4.2.7. The characteristic equation of the Jacobian J* of the system (3) evaluated at an interior equilibrium E* is l2  lTr(J*) + Det(J*) = 0, where Tr(J*) = (u*  m)(m + 2u*  1)  u*(1  u*) and m ¼ ð1u

 Þf2G þeu ð2hG 1Þg

.

2

G ð1þheu Þ *



2



2u Since m + 3u  2 >0 it follows that m + 2u  1 >1  u > 0 and so mcr ¼ mumþ3u > 0. þ2u 1  * 1 Also, it is observed that for h < G , we have Det(J ) > 0 and ðTrðJ ÞÞjm¼m ¼ 0. Thus, there exists m = m* where the trace vanishes. To ensure the existence of a Hopf bifurcation we have to check the transversality condition. d Differentiating the expression for Tr(J*) with respect to m, we get dm ðTrðJ  ÞÞjm¼m ¼ 1  2u  m < 0 and so the transversality condition for a Hopf bifurcation is satisfied. n o Therefore, if max 0; m; 23m < u < 1 and h < G1 , the system (3) undergoes a Hopf bifurcation when m crosses m*. Also, at E* we obtain *

1



ðJ  2iv0 I2 Þ

*

0 1 @ 2iv0 ¼ Dv0 Fuj2  E

DþH

1

a

A

1 Fuj  2iv0  E

and

1

0 1 @

1

DþH

0

a A 2 1 Fuj F  ujE  E  2 2 1 DþH 2 where Dv0 ¼ DþH a FujE  4v0  2iv0 FujE and D0 ¼ a FujE .

ðJ Þ

¼

D0

The left and right normalized eigenvectors a and b of variational matrix J* are !T   iF 2 ujE DþH a ¼ a1 1; i av and b ¼ b1 1; v . 0

0

av2

0 Using a b = 1 we obtain a1 b1 ¼ j ¼ av2 þðDþhÞF . 2 0 ujE "(

Now al Ful j um us b j bm bs ¼ jjb1 j2 " 1 Fuuj E

þ

2 ðDþHÞFuj

F2 E uvjE 2 av0

E

E

2 Fuj

(

2 ðDþHÞFuuj

þi

av0

jjb j2

Fur p uq b p bq ¼ D 1 ðX 11 X 12 þ X 21 X 22 Þ v0 ! ( ðDþHÞF 2

F2 ujE uvjE

av0

1 0 Fuuj  E

þv

, X 21 ¼

ðDþHÞF 2

1 Fuuuj



F1 E uvjE v0

where

Fu1vj  E

F2 ujE uuvjE av20

þ

þ

)#

)

jjb j2

þ D1 0

n

1 X 11 ¼ Fuuj

iðDþHÞF 2 uvjE

av0

) 

(

Fu1vj

E



ðDþHÞF 2

uuujE av0

þi

E

2 þ iðDþHÞ av Fuvj 0

1 ðDþHÞFuj

2 F1 Fuuj E ujE

F2

E

av20



Fu2vj

F1 ujE uuvjE v0



E

o E

þi

2 Fuj F1 E uujE



)# Fr b b ¼ , al Ful j um b j ðJE1  Þ mr u p uq p q

1 ðFuj F2  uuj E

2 Fuj

E

Fu1vj

v0

E

E

2  Fuj F1 Þ  uuj  E

þ

2 4Fuj F2 E uvjE

2 ðDþHÞFuuj



av0

and

E

! E

X 12 ¼

, (

, X 22 ¼ 2i

ðF 2

ujE

v0

E

aD0

h i al Ful j u b j ðJE  2iv0 Þ1 k

2 ðDþHÞFuuj

a

2

Þ F1 uvjE

2 jjb1 j2 ðDþHÞFuuj

F1

þ

E

kr

2  4Fuj F1  2i E uvjE )

F2 F2 ujE ujE uvjE v20

2  v0 Fuuj

E

For these values of Ful i u j u ; Ful i u j ; a and b, we obtain the value of ðFÞm¼m . k

Please cite this article in press as: Bhattacharyya, J., Pal, S., Algae-herbivore interactions with Allee effect and chemical defense. Ecol. Complex. (2015), http://dx.doi.org/10.1016/j.ecocom.2015.09.002

.

G Model

ECOCOM-549; No. of Pages 15 J. Bhattacharyya, S. Pal / Ecological Complexity xxx (2015) xxx–xxx

15

If ðFÞm¼m > 0 then the system (3) undergoes a supercritical Hopf bifurcation as m is increased through m = m* so that the bifurcating periodic orbit is asymptotically orbitally stable. Proof of Lemma 5.2. Let the system (3) is stable at E1. n o ea a e Then G < 4ð1þheÞfað1þheÞðDþHÞg and a > (1 + he)(D + H) hold implies 1þhe 1  4Gð1þheÞ < D þ H. n o a e  D  H  k2 d. 1  4Gð1þheÞ At E1 the eigenvalues of the characteristic equation of the linearized system (8) are m  1  k2 and 1þhe 2 Now, as 1 < m < 1, it follows that m  1  k < 0. n o n o a e a e < D þ H ) 1þhe < D þ H þ k2 d. Also, 1þhe 1  4Gð1þheÞ 1  4Gð1þheÞ Therefore, the stability of the non-spatial system (3) implies the stability of the diffusive system (2). Proof of Lemma 5.3. At E2 the eigenvalues of the characteristic equation of the linearized system (8) are m(1  m)  k2 and n o am em  D  H  k2 d. 1  4Gð1þhemÞ 1þhem pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi Since 0 < m < 1, it follows that m(1  m)  k2 < 0 if k > mð1  mÞ. n o k2 dÞ am em em2 a  D  H  k2 d < 0 if G < 4ð1þhemÞfamðDþHþ 1  4Gð1þhemÞ and a > ð1þhemÞðDþHþ . Now, 1þhem m k2 dÞð1þhemÞg pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ð1þhemÞðDþHþk2 dÞ em2 a and k > mð1  mÞ hold, the diffusive system (2) with strong Allee Therefore, if G < 4ð1þhemÞfamðDþHþk2 dÞð1þhemÞg ; a > m effect is locally asymptotically stable at E2. Proof of Lemma 5.4. If the non-spatial system (3) is locally asymptotically stable at E*, then Tr(J*) < 0 and Det(J*) > 0. Now, Tr(J*) < 0 implies P > 0. Also, Tr(J*) < 0 and Det(J*) > 0 implies Q > 0. Hence, it follows that whenever Tr(J*) < 0 holds, the stability of the non-spatial system (3) at E* implies the stability of the diffusive system (2) at E*. Proof of Lemma 5.5. If E* exists, then we have Det(J*) > 0. Since Tr(J*) > 0, the spatially homogeneous system is unstable at E*. We consider the following cases: Case(a): Let 0 < Tr(J*) < k2. Then P = dk2 + k2  Tr(J*) > 0 and Q = Det(J*) + dk2{k2  Tr(J*)} > 0. Therefore, the spatiotemporal system is locally asymptotically stable at E*. Case(a): Let Tr(J*) > k2. Then, P > 0 if d > Therefore, if

TrðJ  Þk2

TrðJ  Þk2

k2

k2



DetðJ Þ and Q > 0 if d < k2 ðTrðJ  Þk2 Þ. 

DetðJ Þ * < d < k2 ðTrðJ  Þk2 Þ holds, then the spatiotemporal system is locally asymptotically stable at E .

References Appelhans, Y.S., Lenz, M., Medeiros, H.E., Perez da Gama, B.A., Pereira, R.C., Wahl, M., 2010. Stressed, but not defenceless: no obvious influence of irradiation levels on antifeeding and antifouling defences of tropical macroalgae. Mar. Biol. 157, 1151–1159. Banerjee, M., Banerjee, S., 2012. Turing instabilities and spatio-temporal chaos in ratio-dependent Holling–Tanner model. Math. Biosci. 236, 64–76. Bhattacharyya, J., Pal, S., 2013. The role of space in stage-structured cannibalism with harvesting of an adult predator. Comput. Math. Appl. 66, 339–355. Chattopadhyay, J., Venturino, E., Chatterjee, S., 2013. Aggregation of toxin-producing phytoplankton acts as a defence mechanism: a model-based study. Math. Comput. Model. Dyn. Syst. 19 (2), 159–174. Chaudhuri, S., Chattopadhyay, J., Venturino, E., 2012. Toxic phytoplankton-induced spatiotemporal patterns. J. Biol. Phys. 38, 331–348. De Lara-Isassi, G., lvarez-Herna´ndez, S., Collado-Vides, L., 2000. Ichtyotoxic activity of extracts from Mexican marine macroalgae. J. Appl. Phycol. 12 (1), 45–52. Edmunds, J.S.G., McCarthy, R.A., Ramsdell, J.S., 1999. Ciguatoxin reduces larval survivability in finfish. Toxicon 37 (12), 1827–1832. Feng, Z., Liu, R., DeAngelis, D.L., Bryant, J.P., Kielland, K., Chapin III, F.S., Swihart, R.K., 2009. Plant toxicity, adaptive herbivory, and plant community dynamics. Ecosystems 12, 534–547. Feng, Z., Qiu, Z., Liu, R., DeAngelis, D.L., 2011. Dynamics of a plant-herbivorepredator system with plant-toxicity. Math. Biosci. 229 (2), 190–204. ˜ ez, B., Mena-Lorca, J., RojasGonza´lez-Olivares, E., Meneses-Alcay, H., Gonza´lez-Yan Palma, A., Ramos-Jiliberto, R., 2011. Multiple stability and uniqueness of the limit cycle in a Gause-type predator-prey model considering the Allee effect on prey. Nonlinear Anal. Real World Appl. 12 (6), 2931–2942. Hay, M.E., 1997. The ecology and evolution of seaweed-herbivore interactions on coral reefs. Coral Reefs 16 (1), S67–S76.

Henry, D., 1981. Geometric Theory of Semilinear Parabolic Equations, volume 840 of Lecture Notes in Mathematics. Springer, Berlin, Germany. Holt, R.D., Knight, T.M., Barfield, M., 2004. Allee Effects, Immigration, and the Evolution of Species’ Niches. Am. Nat. 163 (2), 253–262. Li, Y., Feng, Z., Swihart, R., Bryant, J., Huntly, N., 2006. Modeling the impact of plant toxicity on plant-herbivore dynamics. J. Dyn. Differ. Equ. 18 (4), 1021–1042. Mistro, D.C., Rodrigues, L.A.D., Petrovskii, S., 2012. Spatiotemporal complexity of biological invasion in a space- and time-discrete predator-prey system with the strong Allee effect. Ecol. Complex. 9, 16–32. Momo, F.R., 1995. A new model for periphyton growth in running waters. Hydrobio 299, 215–218. Poore, A.B., 1976. On the theory and applications of Hopf-Friedrichs bifurcation theory. Arch. Rat. Mech. Anal. 60, 371–393. Rotjan, R.D., Lewis, S.M., 2006. Parrotfish abundance and selective corallivory on a Belizean coral reef. J. Exp. Mar. Biol. Ecol. 335, 292–301. Ruan, S., 1993. Persistence and coexistence in zooplankton-phytoplankton-nutrient models with instantaneous nutrient recycling. J. Math. Biol. 31, 633–654. Sen, M., Banerjee, M., Morozov, A., 2012. Bifurcation analysis of a ratio-dependent prey-predator model with the Allee effect. Ecol. Complex. 11, 12–27. Thieme, H.R., 2003. Mathematics in Population Biology, Princeton Series in Theoretical and Computational Biology. Princeton. Wang, J., Shi, J., Wei, J., 2011. Dynamics and pattern formation in a diffusive predator-prey system with strong Allee effect in prey. J. Differ. Equ. 251 (45), 1276–1304. Wear, D.J., Sullivan, M.J., Moore, A.D., Millie, D.F., 1999. Effects of water-column enrichment on the production dynamics of three seagrass species and their epiphytic algae. Mar. Ecol. Prog. Ser. 179, 201–213. Yasumoto, T., Seino, N., Murakami, Y., Murata, M., 1987. Toxins produced by benthic dinoflagellales. Biol. Bull. 172, 128–131.

Please cite this article in press as: Bhattacharyya, J., Pal, S., Algae-herbivore interactions with Allee effect and chemical defense. Ecol. Complex. (2015), http://dx.doi.org/10.1016/j.ecocom.2015.09.002