Periodic solutions of some ecological models

I. theor. Biol. (1977) 69, 143-152

Periodic Solutionsof Some Ecological Models? FRED BRAUER Department of Mathematics, University of Wisconsin, Madison, Wisconsin 53706, U.S.A. and Department of Applied Mathematics, Weizmann Institute of Science, Rehovot, Israel (Received 19 January 1977, and in revised form

3 1 May 1977)

A result is given on the existence of an asymptotically stable periodic solution of a class of systems of periodic ordinary differential equations. Tbe result is applied to a lake eutrophication model with seasonal effects, and some suggestions are made for the solution of such models.

1. Introduction Since the work of Lotka (1924) and Volterra (193 1) on predator-prey systems, oscillatory solutions of mathematical models for ecological systems have been of interest. Such solutions have been shown to result from various types of models. For example, delay terms may lead to oscillations (May, 1973; Kaplan & Yorke, 1975). For predator-prey systems modelled by a pair of ordinary differential equations, there can be an unstable equilibrium and a stable limit cycle (Kolmogorov, 1936; Rosenzweig, 1971; Brauer, 1975). Recently it has been shown that ordinary differential equation models for more than two species in competition can lead to periodic solutions (Gilpin, 1975), although this can not happen if only two species are in competition. The purpose of this note is to show how the inclusion of periodic terms, such as seasonal variations, in ecological systems can give rise to periodic solutions. I concentrate on the situation where the system has an asymptotically stable equilibrium if the periodic terms are neglected, leaving for future investigations the more difficult question of the effect on an oscillatory t This work was supported in part by the Wisconsin Alumni ResearchFoundation. 143

144

F.

RRAUER

system of the inclusion of periodic terms. My work was motivated by discussions with Professor A. Nir and S. Lewis on models for lake eutrophication problems, and we indicate how our results may be applied to such questions.

2. The Existence of Asymptotically We wish to study an unperturbed

Stable Periodic Solutions

system of ordinary differential equations : x’ =f(x)

(1)

which has a periodic solution pa(t) of period w (not necessarily of least period w), and a perturbed system: Y’ = f(Y) + &

Y),

(2)

with g(t, y) periodic with period w in t. Here x, y, f, and g are n-dimensional column vectors, primes denote differentiation with respect to the scalar independent variable t, and E is a scalar parameter. We are interested in the existence and asymptotic stability of periodic solutions of the perturbed system (2). The basic result is the following: Theorem 1: Let&) and g(t, v) be continuously differentiable with respect to the components of y and let g(t, v) be periodic in t with period w for each y. Suppose that the system (1) has a periodic solution PO(t) of period w and suppose that the characteristic exponents of the variational system of (1) with respect to the solution PO(t), namely:

2’ = fyrPo(olY~

(3)

all have negative real part. Then the perturbed system (2) has an asymptotically stable periodic solution p(t, E) with period w, and lim p(t, E) = pe(t) E-+0

(4)

for all sufficiently small E. Theorem 1 is classical, and may be found, for example, in (Coddington & Levinson, 1955, pp. 348-350) or (Halanay, 1966). Here, fy[pO(t)] is the n x II matrix with the partial derivative of the ith component of f(v) with respect to the jth component of y evaluated at y = PO(t) in the ith row, jth column. The variational system (3) is thus a linear system with periodic coefficients, and its characteristic exponents are defined through the Floquet theorem for such systems (Coddington & Levinson, 1955, pp. 78-80). They are the eigenvalues of the matrix R, where a fundamental matrix of the system (3) has the form P(t) etR with P(t) periodic. The difficulty in applying Theorem 1 arises from the fact that these characteristic exponents can not be calculated explicitly in general.

PERIODIC

SOLUTIONS

OF

ECOLOGICAL

MODELS

145

In the special case of Theorem 1 where the periodic solution p,,(t) of (1) is a constant solution corresponding to an asymptotically stable critical point x* of (l), the variational system (3) is the linear system with constant coefficients : z’ = fy[%&, (5) and its characteristic exponents are just the eigenvalues of the matrix &[x,J, which is the community matrix of the system (1) at the equilibrium x,. Thus in this case, Theorem 1 reduces to the following result which has not previously been formulated explicitly. Theorem 2: Let f(y) and g(t, y) be continuously differentiable with respect to the components of y and let g(t, y) be periodic in t with period w for each y. Let x, be a critical point of (l), i.e. a solution of the equation f(x) = 0, which is asymptotically stable in the strong sense that the eigenvalues of the matrixf,[x,J all have negative real part. Then the perturbed system (2) has an asymptotically stable periodic solution ~(t, c) of the same period w for all sufficiently small E, with: lim p(t, &) = x,.

(6)

E-+0

It is possible for a critical point x, of (1) to be asymptotically stable but for the matrixjYY[xco] to have eigenvalues with real part zero. This situation, which is structurally unstable in the sense that it may be altered by small perturbations, is not covered by the above theorem. If this possibility is excluded, the content of the theorem is essentially that if an autonomous system with an asymptotically stable equilibrium is disturbed by a sufficiently small periodic perturbation, the result is an asymptotically stable periodic oscillation whose amplitude tends to zero with the magnitude of the perturbation. The purpose of this note is to indicate how Theorem 2 may be applied to ecological systems. Since Theorem 2 establishes the existence of a periodic solution only for sufficiently small E, it will be necessary to give operational procedures for judging whether there is a periodic solution even though the existence of a periodic solution can not be established rigorously. This may be done by numerical simulations with reasonable confidence but not absolute certainty. 3. Application to a Phytopladcton-Nutrient

Model

The growth of phytoplankton in a lake may be modelled by a predatorprey system with the nutrient (phosphorus) as prey and the phytoplankton as predator. While such models usually separate different layers of the lake T.R. 10

146

F.

BRAUER

(lmboden, 1973, 1974), it is possible to consider the hypolimnion as a subsystem whose outputs may be viewed as parameters in a pair of equations which describe the epilimnion. The advantage of doing this is that the hypolimnion equations are partial differential equations, and the theoretical analysis of the full system would be much more difficult; i.e. we simplify the system by dealing only with the epilimnion. This makes the model amenable to a theoretical analysis at the cost of a loss in realism. Conversations with Professor A. Nir and S. Lewis (at the Weizmann Institute of Science, Rehovot, Israel), have suggested the model:

p2l-V

@=

~ -ps+cj

+hn+F3-p3+4

nt = pa4 -+p,II-PJI+p3+4

&IV) cos g &IV

(t-90), (7)

cos

g5

(t -

90).

p3++

where 4 is the phosphorus concentration, and lI is the phytoplankton concentration. Here P,, P,, Pq, Pg, F and E are positive constants with E < P2. The system (7) may be viewed as a predator-prey system, but it is also a nutrient recycling system. The term involving P6 represents nutrient leaking from the system and lost, while the F3 term is an input rate of phosphorus to compensate for this leakage. In system (7), the first three terms of each equation correspond to f(y) in equation (2), while the last term corresponds to the periodic perturbation .sg(t, y) in (2). For E = 0, the system (7) is autonomous, and has an equilibrium at: dJ=

P,(P, +pcd p2--(p4+W

n = 5 ) PC5

provided P, > P,+P,. Linearization of the system (7) about this equilibrium gives the coefficient matrix (corresponding to fY[xoc] in Theorem 2)

whose eigenvalues have negative real parts. Thus Theorem 2 shows that for sufficiently small E the system (7) has a periodic solution which tends to the equilibrium as E + 0. Since this periodic solution is asymptotically stable, (it can in fact be proved that it is globally asymptotically stable), every solution tends to it as t -+ co. Thus the periodic solution can be approximated numerically by solving the system (7) numerically with any initial conditions.

PERIODIC

SOLUTIONS

OF

ECOLOGICAL

MODELS

147

Solutions tend to the periodic solution quite rapidly, and a computer simulation over a time period of one or two years (cycles) gives a very good approximation to the periodic solution. The fact that Theorem 2 gives no estimate of the range of E for which (7) has a periodic solution remains a serious problem. From a practical point of view, we suggest that a computer simulation beginning with E = 0 and gradually increasing E to and beyond the desired value should indicate whether the periodic solution breaks down or remains. For this purpose, a simulation which gives a phase-plane portrait rather than 4 and II as functions oft is convenient. It also gives useful estimates of the amplitude of the oscillation of various values of E. In (7) the terms containing Pz are the same in both equations, which means that 100% efficiency of transfer is assumed. It might be more realistic to replace (7) by a system:

qy= -- p9-W +P4n+F3-p3 p,+$J l-I’=

&W

-kos+$

CPJW -cP411-P,II+p,++ p,++

g5 (t - 90), (8)

CdIC#J

cos g5 (r-90),

with 0 < c < 1, but this would cause no essential change in the results. The point is that Theorem 2 gives a method of treating a large class of systems; the examples are given to show the approach and the type of conclusions which may be drawn. More complicated problems and more realistic examples may be studied by the same approach. A refinement of the model (7), currently being studied by Professor A. Nir and S. Lewis (Nir & Lewis, 1977) 4 = -[P,+~cos~~(f--120)]~~+~~~~

+4)+P411+F,, 3

1

Is =

n4

z(a +BWUJ3 + 4)

(9) - Pqrl-

PJI,

with Z, a, j? positive constants and E < P,. Just as for the model (7), it is not difficult to calculate that if E = 0 there is an asymptotically stable equilibrium at:

4=

+P,)W,

zP3(P4

+BF3)

P2P, - Z(P, + P6)(aP6 +j3F3) ’

provided : P2P6

-

z(P4

+P6)(aP6 + /?I;;) > 0.

Theorem 2 implies that for sufficiently small E, the system (9) has a periodic

148

F.

BRAUER

solution. A simulation has been performed on the IBM 370/168 at the Weizmann Institute of Science with CI= 0.2, /I = O*Ol, Z = 4.8, P, = 5, P, = 100, P4 = 0.2, P6 = O-3/4*8, F3 = 1 and various values of E, all with initial values 4(O) = 440, II(O) = 40. The results are shown in Figs l--7, which give the relation between II and 4, with the direction of increasing t indicated by an arrow. The following deductions can be made from examination of these figures. 500

FIG.

4

1.

12

8

16

20

9 FIG.

2.

(1) For 0 I E 5 6 there is a periodic solution with period one year which is approached quickly (even for the initial conditions used which are far from the periodic solution, the solution is indistinguishable from the periodic orbit after five months for each E). The speed of the approach to this periodic

PERIODIC

SOLUTIONS

OF

ECOLOGICAL

MODELS

16

4

12

8

16

20

+ FIG. 3.

32

24 l-l 16

8

16

24

32

60

60

9 FIG. 4.

60 II 40

20

40

FIG. 5.

I

149

150

i-.

BRAULR

9 Fro.

6.

60

n

40

80

120

160

200

9 FIG.

7.

orbit means that the periodic orbit may be used for predictions, even over time periods less than a full cycle. It also means that sudden changes, such as spring upwellings with resultant large changes in phosphorus concentration can be modelled by discontinuities in the solution, which can in turn be described by changes in initial values followed by rapid return to the same periodic orbit. (2) The periodic orbit for each E > 0 has the equilibrium point for E = 0, namely 4 = 9,977, II = 16, in its interior. For 0 _< E I 3, the amplitude appears to vary approximately linearly with E, but for E > 3 the amplitude increases much more rapidly. (3) For E 2 5, the simulation indicates that II = 0 over part of each cycle. In fact, ifs I Pz, it is easy to see that II’ 2 - (P4+ P&I, which implies that

PERIODIC

SOLUTIONS

OF

ECOLOGICAL

MODELS

151

while II could decrease exponentially, II can not reach zero in finite time. Thus if E < P2, there is no danger of extinction. If e < P2, II is bounded away from zero. The essential feature of the model is that it has a periodic solution which is approached rapidly by every solution. Thus knowledge of the periodic solution, easily obtained by a simulation, makes possible predictions of behavior of the system independent of initial conditions, and even allows for predictions of the effects of discontinuous changes after a short settling time. The model described here takes into account seasonal changes, with a period of one year. To build in daily effects it should be possible to consider the model over a time period over which the annual terms are almost constant, and then consider daily effects as a perturbation. This would give a periodic solution with a period of one day over such time periods, suggesting that there should be an annual cycle with daily cycles super-imposed on it.

4. Conclusions

Our main result says that if an autonomous system has an asymptotically stable equilibrium, then the effect of a periodic perturbation is to produce a periodic solution with the same period. The fact that this periodic solution is asymptotically stable makes it easy to approximate numerically. If the periodic orbit is approached rapidly by other solutions, then knowledge of the periodic orbit makes it easy to predict the behavior of the system. A more difficult problem arises when the unperturbed system has a nonconstant periodic solution. In this case, the results of Levinson (1950) together with the theory of differential equations on a torus (see, for example, Coddington & Levinson, 1955, Chapter 17) show that there is a rotation number ~1,continuous in E, and that if /J is rational, say ,u = p/q, then the system (2) has a periodic solution of period qw for sufficiently small E while if p is irrational there is no periodic solution. If /J is irrational the orbits are oscillatory in nature, even though there is no periodic solution. However, this suggests that the prediction of the effect of seasonal variations on oscillatory ecological systems is difficult. It is not true in general that a seasonal variation will produce a direct seasonal effect. We are now aware that periodic solutions may arise in ecological models either because of the existence of an unstable critical point or because of a periodic perturbation, but the effect of combination of both these factors may not be a periodic solution or may be a periodic solution whose period is different from that of the perturbation. Another factor which can produce

152

F.

BRAUER

periodic solutions is the presence of delay terms (Kaplan & Yorke, 1975; May, 1973); it would be interesting to determine whether a delay term together with a periodic perturbation necessarily produces a periodic solution. REFERENCES (1975).Znt. J. Cont. 23, 541. E. A. & LEVINSON, N. (1955). Theory of Ordinary Differential Z@ations. New York: McGraw-Hill. GILPIN, M. E. (1975).Am. Nat. 109, 51. HALANAY, A. (1966). Diferential &uations: Stability, OscilIations, Time Lugs, pp, 253-5, 263-4. New York: Academic Press. BRAUER, F. CODDINGTON,

IMB~DEN, D. M. (1973). Schweiz Z. Hydrof. 35, 29. IMBODEN, D. M. (1974). Limnol. Oceanogr. 19, 297. KAPLAN, J. & YORKE, J. A. (1975). Siam J. Math. Anal. KOLM~CXIROV, A. N. (1936). G. Inst. ital. Att. 7, 74. LEVINSON, N. (1950). A. Math. 52, 727. LOTKA, A. J. (1924). Elements of Physical Biology.

6, 268.

Baltimore: Williams & Wilkins. [Reissuedas EIements of Mathematical Biology, New York: Dover (1956)]. MAY, R. M. (1973). Ecology 54, 315. NIR, A. & LEWIS, S. (1977). To appear. ROSENZWEIG, M. L. (1971). Science 171, 385. VOLTERRA, V. (1931).Lecoons sur la Thkorie Mathkmatique de la Lutte pour la Vie. Paris: Gauthier-Villars.