Stochastic modelling for determining zooplankton abundance

Stochastic modelling for determining zooplankton abundance

e c o l o g i c a l m o d e l l i n g 2 0 4 ( 2 0 0 7 ) 372–378 available at www.sciencedirect.com journal homepage: www.elsevier.com/locate/ecolmod...

332KB Sizes 0 Downloads 49 Views

e c o l o g i c a l m o d e l l i n g 2 0 4 ( 2 0 0 7 ) 372–378

available at www.sciencedirect.com

journal homepage: www.elsevier.com/locate/ecolmodel

Stochastic modelling for determining zooplankton abundance Oren Karouby a,b , Athanassios Iliadis a,∗ , Jean-Pierre Durbec b , Virginie Riandey b , Franc¸ois Carlotti b a b

Department of Pharmacokinetics, Faculty of Pharmacy, 27, bld Jean Moulin, 13385 Marseilles Cedex 5, France L.M.G.E.M. Centre d’Oc´eanologie de Marseille, Campus de Luminy, Case 901, 13288 Marseilles, France

a r t i c l e

i n f o

a b s t r a c t

Article history:

The most reliable method for studying the dynamic behaviour of a population of organ-

Received 8 March 2006

isms in a marine pelagic environment is to analyse the biomass or the abundance spectra.

Received in revised form

The spectra presented in semi logarithmic scale features the measurement variability as

3 January 2007

a reversed cone embodying the recorded abundances. This distribution in data variability

Accepted 19 January 2007

motivated us to propose the existence of a process uncertainty.

Published on line 21 March 2007

To analyse abundance spectra in stochastic context, we propose a compartmental configuration in which, at a given time, marine organisms are located according to their size

Keywords:

in the compartments of a compartmental chain. The movement of organisms across com-

Stochastic modelling

partments is described in probabilistic terms expressing transfer probabilities. The model

Abundance spectra

considered at steady state allows the expression of the process uncertainty and the con-

Marine population

struction of the confidence corridors embodying the abundance spectra data. The developed model may be used for comparison of several collections of data, for analysis of environmental factors influencing the marine populations, or for new developments involving longitudinal data. © 2007 Elsevier B.V. All rights reserved.

1.

Introduction

Probably the most reliable method for studying the dynamic behaviour of a population, such as that of organisms in a marine pelagic environment, is to classify these organisms according to the biomass or to the abundance spectra. Such spectra, usually called abundance graphics, describe the observed frequency of organisms classified according to their measured sizes, weights, or volumes. New instruments for counting and measuring plankton like optical plankton counter will increase the use of such spectra (Herman et al., 2004). For instance, the abundances graphics of zooplankton have higher values for small organisms and progressively lower val-



Corresponding author. Tel.: +33 491835645; fax: +33 491835667. E-mail address: [email protected] (A. Iliadis). 0304-3800/$ – see front matter © 2007 Elsevier B.V. All rights reserved. doi:10.1016/j.ecolmodel.2007.01.011

ues for larger organisms (Sheldon et al., 1972; Rodriguez and Mullin, 1986). Although the working hypotheses to conduct experiments, to gather measurements, and to perform such analyses are well established and accepted by several authors (Platt and Denman, 1977, 1978; Heath, 1995; Zhou and Huntley, 1997), there is no general standardised procedure for analysing the abundance spectra (Zhou, 2006). Platt and Denman (1977, 1978) were inspired from observations reported by Sheldon et al. (1972) to propose the theory of biomass spectrum in a pelagic ecosystem. The theory is based on the assumption that the system is at the steady state without immigration or emigration of organisms. Mechanisms for organisms’ renewal and elimination from the system were considered the biomass flux associated with the respiration

e c o l o g i c a l m o d e l l i n g 2 0 4 ( 2 0 0 7 ) 372–378

and the growth rates, and the death rate due to predation, respectively. In that model, biomass flux propagation is conceived as moving from small organisms towards larger ones. The biomass flux is estimated from the slope of the regression line that fits the normalised biomass spectrum. The introduction of the OPC instrumentation in the early 1980s renewed the interest in analysing biomass or abundance spectra. Heath’s model (Heath, 1995) was based on the Platt and Denman theory (Platt and Denman, 1977, 1978) and it described the rate of variation of the organisms’ number as a function of the organisms’ weight. In that model, the biomass flux was associated with the death and the growth rates and it was estimated from a linear regression on the biomass data. Zhou and Huntley developed a subsequent global approach (Zhou and Huntley, 1997) for analysing the normalised biomass spectra by establishing a link between Heath’s assumptions and the Platt and Denman theory. More precisely, Zhou and Huntley included in the model the growth rates for each cell of organisms’ sizes, and they derived the differential equations describing the biomass flux propagation from small organisms towards larger ones. The obtained differential equations adequately predict the biomass variations and the growth rate for each cell of organisms’ sizes. Therefore, this model proved highly reliable to fit a function to the observed data in biomass spectrum. On the other hand, mathematical modelling is a valuable tool aimed the analysis of abundance spectra. The proposed models soundly depend on the set goal, the experimental procedure, and the observed data, e.g., the procedure used for data classification, the kind of dispersion diagram obtained, etc., all that before any analysis using mathematical models. Mathematical models are used not only to extract fundamental features from abundance spectra concerning the population of studied organisms, but also, given the degree of modelling reliability, to simulate and predict under various experimental conditions the subsequent behaviour. Nevertheless, experimental conditions, e.g., temperature, salinity, hydro dynamism, feeding, etc., highly influence the real process and consequently play an important part in the use of models in conditions far from those for which they were established. In the context of the present study, we propose to analyse the observed data directly by a stochastic model instead of the classical analysis of an abundance spectrum. The model regresses on the observed data.

2.

Material and methods

2.1.

Sampling site and measurements

The marine station, called SOFCOM (Station d’Observation ´ Frioul du Centre Oceanographique de Marseille), is located in the Bay of Marseilles (43◦ 14.3 N/5◦ 17.3 E) at a depth of 60 m inside a low polluted sector. During 4.5 months between 4 April 2001 and 7 July 2001, 11 marine samples were gathered at irregular time intervals. To establish the so-called abundance spectrum, the number of organisms was counted inside a fixed sampling volume for each size between 280 and 1500 ␮m. Recent progress in electronics and acoustics contributed to

373

studying the zooplankton abundance and distribution in situ. The optical plankton counter (OPC) is an operational instrument that detects, sizes, and counts individual particles based on measuring the attenuance, or diminution in intensity, of a collimated light beam intercepted by transiting particles (Herman, 1988). The OPC can determine the size of individual particles with effective diameters from 250 ␮m to 2 cm. Thus, nearly the full size spectrum of mesozooplankton can be recorded. OPC automatically provides abundance spectra as a function of organisms’ sizes and considerably shortens analysis and improves accuracy with respect to traditional approaches, e.g., the binocular optical instrumentation (Riandey et al., 2005; Sourisseau and Carlotti, 2006). Using allometry relations, other relevant indexes such as biomass and production fluxes can be directly obtained from the OPC instrumentation (Riandey et al., 2005; Sourisseau and Carlotti, 2006; Zhou, 2006).

2.2.

Data conversion

Presentation of the raw data spectra in semi logarithmic scale features the measurement variability as a reversed cone embodying the recorded abundances (Fig. 1). This distribution in data variability motivated us to propose the existence of a process uncertainty. Since the measured abundances ai are the ratio of the number ni of particles contained in the sample divided by a fixed volume V, and by assuming detection of at least one particle, i.e., min(ai ) = 1/V, we can convert ai to ni by setting ni = ai /min(ai ). In fact, all measured abundances are now converted into integer numbers with min(ni ) = 1.

2.3.

Stochastic model

The developed model associates a compartment with each size of organisms, which varies as integers between 280 and 1500 ␮m in radii. Thus, the model is represented by a directed chain of compartments i with i = 0, . . ., m. The input (i = 0) and output (i = m) ends correspond to the smallest (280 ␮m) and largest (1500 ␮m) observable sizes, respectively. Each

Fig. 1 – Distribution of the observed abundances ai as a function of particles’ sizes i. Data corresponds to the spectrum number 2.

374

e c o l o g i c a l m o d e l l i n g 2 0 4 ( 2 0 0 7 ) 372–378

for the growth constant rate i is the increasing exponential model: i = 0 exp(ˇi)

(2)

with parameters 0 and ˇ. To build in the global model, we evaluate pi (t + t) for a given compartment i by enumerating the stochastic events and writing up the associated probabilities: Fig. 2 – Diagram of the compartmental structure proposed for modelling abundance spectra.

compartment i accounts for the following transfers pictured in Fig. 2:

(1) organism reaches i from compartment S in t, ei Rt, (2) organism reaches i from i − 1 by growing in size, pi−1 (t)i−1 t, and (3) organism remains in i during t + t, pi (t)[1 − (i + )t]. By assuming that the above mentioned events are independent, we obtain: pi (t + t) = pi (t)[1 − (i + )t] + pi−1 (t)i−1 t + ei Rt

(1) inflow from the environment associated with the rate ei , (2) transfers from compartment i − 1 or toward compartment i + 1 with organisms growing in size, and (3) elimination because of predation or death.

This probabilistic transfer model becomes a differential equation at the limit t → 0:

The environment is depicted as a “source” compartment S receiving organisms with a rate R and dispatching them in compartments according to their sizes. In such a case, inflows ei can be assumed following either:

If the process is at steady state (functional hypothesis H0), i.e., dpi (t)/dt = 0, the differential equations over compartments i lead to the following recursive relations:

• A model inspired from the Poisson distribution law:

which quantify the probabilities i with i = 0, 1, . . . at the steady state. By introducing (1) and (2) in the above recurrent relations, we obtain:

ei =

i exp(−) i!

(1a)

(1b)

with parameter  > 1 and where () is the Riemann Zeta function, or • A simple boundary condition: e0 = 1

and

ek = 0,

for k = 1, 2, . . .

Re0 0 + 

e0 , +1

and

k =

qk =

Rek + k−1 k−1 k + 

for k = 1, 2, . . .

ek +  exp[ˇ(k − 1)]qk−1  exp(ˇk) + 1

for k = 1, 2, . . . (3)

−

(i + 1) ()

0 =

q0 =

with parameter  > 0, or • A model inspired from the Pareto distribution law: ei =

dpi (t) = Rei + i−1 pi−1 (t) − (i + )pi (t) dt

(1c)

The key variables materialising the movement of organisms along the compartmental chain are a marginal probability pi (t) that “an organism initially present in the compartment S belongs in the compartment i at time t” and two conditional probabilities expressing the movement of organisms out of a given compartment i:

and i = qi

for i = 0, 1, . . .

The above relations involve the positive dimensionless parameters ,  = 0 /, ˇ, and = R/. The newly introduced parameters  and do not allow discrimination between incoming and out coming mechanisms, but they provide a structurally identifiable model. The number of organisms Ni associated with the compartment i, is assumed random, following the binomial distribution (distributional hypothesis H1) Ni ∼ Bin(N,i ), where N is the total number of organisms implied in the process under steady-state conditions. Thus, the probability of obtaining ni tiles among N with prior probability i is: n

• the probability t that “a given organism resident in a compartment at time t will next transfer out of that compartment at time t + t” because of predation, and • the probability i t that “a given organism resident in a compartment i at time t will next transfer to compartment i + 1 at time t + t” because of growing size. A good candidate

(4)

i Pr(ni ; N, i ) = CN ni (i ) (1 − i )

N−ni

(5a)

with expectation and variance given by E[Ni ] = Ni respectively.

and

Var[Ni ] = Ni (1 − i )

(5b)

e c o l o g i c a l m o d e l l i n g 2 0 4 ( 2 0 0 7 ) 372–378

When N > 50 and i < 0.1, the binomial distribution is accurately approximated by the Poisson distribution (distributional hypothesis H2) Ni ∼ Poi( i ), with parameter i = Ni . In that case, the probability of obtaining ni tiles for the ith compartment is:

ini

Pr(ni ; i ) =

ni !

exp(− i )

(6a)

with expectation and variance given by E[Ni ] = Var[Ni ] = i .

2.4.

(6b)

Parameter estimation

Given the data ni i ∈ I where I is the set of indices of the M sampled compartments, the maximum likelihood criterion could be used for parameter estimation. This criterion J1 may be grounded on the distributional hypothesis H1 to compute the conditional probability from the binomial distribution. For the M observed data, J1 gives:

 J1 = ln





Pr(ni ; N, i )

i∈I

= (M + 1) ln (N + 1) −

+



ni ln i +



i∈I



ln (N − ni + 1) −

i∈I



ln (ni + 1)

i∈I

(N − ni ) ln(1 − i )

This last expression was obtained by using the relation ln(n!) = ln (n + 1). The criterion J1 depends on , , ˇ, and through the model i defined by Eqs. (1)–(4), and seemingly on parameter N. Although structurally identifiable in criterion J1 , parameters and N are numerically non identifiable (Norton, 1986). This fact can be verified under the distributional hypothesis H2 for which the log-likelihood is:

 J2 = ln

=



and the Hessian matrix was computed. The numerical precision on estimates is evaluated from the inverted Hessian.

3.

Results

The three models given by Eqs. (1a)–(1c) were identified and the best model was selected by using the Akaike Information Criterion (Ljung, 1999). The best model is that with the simple boundary condition (Eq. (1c)). The two others frequently lead to inconsistent estimates for parameter . As a fitting example, Fig. 3 presents for sample number 2 the observed data, the expectation function E[Ni ], and the confidence corridors E[Ni ] ± Var[Ni ]1/2 expressing the process uncertainty. Modelling the process uncertainty allows the confidence corridors to account for the observed variability. Table 1 summarises the estimated parameters and the associated numerical confidence intervals. We note that all parameters are correctly estimated. By using inequalities 7 and parameter estimates, the last two columns of Table 1 give the max and min values for parameters and N, respectively. From an estimation point of view, the weighted leastsquares criterion operating on a summation of the weighted discrepancies between observations ni and predictions E[Ni ] could be employed (Seber and Wild, 1989):

 ( n − E[N ] )2 i i i∈I

i∈I

375

i2

Expectation E[Ni ] allows the regression of the model function to the observed data. Variance Var[Ni ] is used to take into account the model uncertainty associated with the probabilistic movement of organisms independently of the measurement error (Matis and Tolley, 1979; Matis and Wehrly, 1979). Therefore, the weighting factor i2 = Var[Ni ] + { KE[Ni ]}2 accounts for both the process uncertainty Var[Ni ] and the





Pr(ni ; Ni )

i∈I

ni ln(Ni ) −



i∈I

Ni −

i∈I



ln (ni + 1)

i∈I

Since i is proportional to (Eq. (4)), the only structurally identifiable parameter is the product = N . Subsequently, four dimensionless parameters, , , ˇ, and will be estimated from observed data by maximising the criterion J2 . m  ≤ 1, from Eq. (4) we obtain: Since i=0 i

 ≤

m  i=0

−1 qi

and

N≥

m 

qi

(7)

i=0

Maximisation of J2 is achieved by nonlinear optimisation techniques implying the quasi-Newton algorithm available within MATLAB (2004). Sensitivity functions were evaluated

Fig. 3 – Model fitting on data from spectrum number 2. The graphic expresses the distribution of the particles as a function of particles’ sizes i. Open circles represent the observed number of particles ni , the solid line represents the expectation function E[Ni ], and dashed lines define the confidence corridor E[Ni ] ± Var[Ni ]1/2 expressing the process uncertainty.

376

e c o l o g i c a l m o d e l l i n g 2 0 4 ( 2 0 0 7 ) 372–378

Table 1 – Estimated parameters of the abundance spectra No.

1 2 3 4 5 6 7 8 9 10 11

−ln J

199.32 178.18 185.19 224.81 175.56 221.31 197.94 197.30 174.42 202.02 223.07



ˇ × 103



Esta

SErb

Esta

SErb

Esta

117.24 129.15 129.33 134.10 128.41 130.69 137.90 171.38 118.69 150.19 142.73

4.96 5.97 6.36 4.49 5.11 4.32 4.95 9.03 4.72 6.14 4.80

2.788 2.047 2.711 0.921 1.440 1.268 1.706 1.739 1.964 1.724 1.380

0.354 0.448 0.643 0.225 0.421 0.277 0.491 1.062 0.471 0.507 0.332

13,606 14,036 13,699 22,952 15,214 19,162 17,143 13,348 14,190 14,386 18,261

c

Nd

1.054 1.031 1.067 1.003 1.010 1.008 1.023 1.052 1.020 1.033 1.015

12,909 13,620 12,835 22,874 15,059 19,018 16,757 12,691 13,919 13,933 18,000

SErb 521.9 529.7 740.6 555.1 486.6 523.7 598.4 1088.9 492.0 591.2 534.6

Other symbols are defined in the text. a b c d

Estimates. Standard error of estimates. Upper limit of . Lower limit of N.

measurement error {KE[Ni ]}2 (variance of the measurement error proportional to the prediction). The weighted leastsquares criterion may be considered as an approximation of the likelihood function when the discrepancies between ni and E[Ni ] follow a Gaussian distribution (Seber and Wild, 1989). In the abundance spectra analysis, J1 or J2 are preferred to the weighted least-squares criterion because they lead without additional hypothesis to the likelihood function for parameter estimation. The weighted least-squares criterion can be used whatever the essence of random variables, continuous or discrete. It could be directly applied for observed abundances, whereas J1 or J2 concern only discrete random variables. Since > 1, i.e. R > , the input rate is slightly greater than the predation rate. This is a necessary condition for a nonempty compartmental configuration. The total number of organisms N implied in the process under steady-state conditions is much larger than the observed number recorded by the OPC instrumentation. For instance and for sample number 2 of observed data, N > 13,620 although the observed number equals 1004. This may reflect the inability of OPC to record the abundances for each particle size measured by integer values in micrometers. From Eq. (2), the growing constant rate indicates size doubling without predation, i.e., ln 2/ˇ. For instance and for sample number 2, one gets 339 compartments or 339 ␮m.

4.

Discussion

Since the beginning of the 1980s, the interest in ecological system dynamics has greatly increased. In fact, these systems are characterised by a high degree of heterogeneity in space and time, by many interacting factors and by feedback mechanisms. Such heterogeneity and complexity result in unpredictability of the process behaviour and in stochastic fluctuations (Jorgensen et al., 1996). In these cases, the variation of body size distribution of different organisms, e.g., the pair zooplankton–phytoplankton, appear to be a key parame-

ter in the dynamical behaviour of ecological systems (Gaedke, 1993). Moreover, it was recently evidenced that the growth rate of zooplankton may drive the system at the edge of chaos (Jorgensen, 1995; Mosekilde, 1996). Therefore, modelling zooplankton dynamics as a part of a holistic ecosystem gained increased interest. Commonly, planktonic ecosystems are described by compartmental models, each of them representing a trophic level or taxonomic group. In such models, the process rates used in the dynamic model vary according to the body size of planktons. Environmental stochasticity was accounted for by expressing the growth rates by probability density functions to reflect the heterogeneity and observed variability in samples (Engen et al., 1998; Chaloupka, 2002). On the other hand, stochastic compartmental configurations have largely been used for modelling biological, medicinal, pharmaceutical, and ecological processes (Matis et al., 1979; Jacquez, 1996; Matis and Kiffe, 2000; Macheras and Iliadis, 2005) and applied to zooplankton models (Carlotti et al., 2000; Carlotti, 2001). However, such models applied to zooplankton have not fully explored the use of stochastic approach to analyse model properties in connection with data. The goal of this study was to develop a compartmental stochastic model able to describe reliably the observed abundance spectra of a zooplankton population in marine stations in the Bay of Marseilles. The present study was greatly inspired by the pioneering works of Platt and Denman (1977, 1978), Heath (1995), and Zhou and Huntley (1997) establishing the first models concerning a pelagic ecosystem, and the contributions of Matis (Matis and Tolley, 1979; Matis and Wehrly, 1979; Matis and Kiffe, 2000) to stochastic compartmental modelling theory. Therefore, the obtained model mixes deterministic and stochastic elements from the existing compartment theories. As in the theory of Platt and Denman (1977, 1978), the model we developed is based on the fundamental assumption of the steady state. To analyse abundance spectra in stochastic context, we proposed a compartmental configuration where,

e c o l o g i c a l m o d e l l i n g 2 0 4 ( 2 0 0 7 ) 372–378

at a given time, marine organisms are located according to their size in the compartments of a compartmental chain. The movement of organisms across compartments is described in probabilistic terms expressing transfer probabilities. We also assumed that predation does not depend on the size of the population, although it is generally accepted that the fewer the number of organisms, the higher the predation (prey deficiency). The model is purely descriptive, aiming to fit the observed data in order to translate them in a reduced set of model parameters. The estimated parameters allow complete and consistent description of the huge mass of observed data. The main innovation in that model is the stochastic background used to introduce the process uncertainty and thereby describe the distribution of data by an inversed cone on the abundance spectra. In fact, the expectation of the model is identical to the prediction from the equivalent deterministic model. However, the fundamental difference lies in the process uncertainty introduced by the stochastic model. The principal technical difficulties were first, the selection of the best model associated with an adequate parameterisation and, second, the derivation of the estimation criteria according to the maximum likelihood principle. The model with a simple boundary condition proved the most reliable, and the Poisson distribution allowing parametric identifiability served as the basis of the estimation criterion. With respect to the literature, two contributions characterise the present study: • The presence of the process uncertainty embodying through the confidence corridors the observed variability in abundance spectra, and • The non-involvement of the subjective spectrum normalisation found in current analyses. On the contrary, abundance measurements must be converted into integer numbers, a necessary condition for using J1 or J2 estimation criteria. The model presented may have several applications. For instance, it may be used for comparison of several collections of data, for analysis of environmental factors influencing the marine populations (temperature, salinity, hydrodynamics, etc.), or for developments involving longitudinal data where fluctuations in observed abundances are recorded over time in pelagic ecosystems (Mosekilde, 1996; Ray et al., 2001).

Acknowledgments This research was supported by the AMPLI (F. Carlotti) and ZOOPNEC (J.P. Durbec and F. Carlotti) projects of the French ˆ programme “Programme National Environnement Cotier”.

references

Carlotti, F., 2001. Plankton: population dynamics models. In: Steele, J.H., Turekian, K.K., Thorpe, S.A. (Eds.), Encyclopedia of Ocean Sciences. Academic Press, London, pp. 2247–2256.

377

Carlotti, F., Giske, J., Werner, F., 2000. Modelling zooplankton dynamics. In: Harris, R.P., Wiebe, P., Lenz, J., Skjoldan, H.R., Huntley, M. (Eds.), Zooplankton Methodology Manual. Academic Press, London. Chaloupka, M., 2002. Stochastic simulation modelling of southern Great Barrier Reef green turtle population dynamics. Ecol. Model. 148, 79–109. Engen, S., Bakke, O., Islam, A., 1998. Demographic and environmental stochasticity: concepts and definitions. Biometrics 54, 840–846. Gaedke, U., 1993. Ecosystem analysis based on biomass size distribution: a case study of a plankton community in a large lake. Limnol. Oceanogr. 38, 112–127. Heath, M.R., 1995. Size spectrum dynamics and the planktonic ecosystem of Loch Linnhe. ICES J. Mar. Sci. 52, 627–642. Herman, A.W., 1988. Simultaneous measurement of zooplankton and light attenuance with a new optical plankton counter. Cont. Self Res. 8, 205–221. Herman, A.W., Beanlands, B., Phillips, E.F., 2004. The next generation of optical plankton counter: the laser-OPC. J. Plankton Res. 26, 1135–1145. Jacquez, J.A., 1996. Compartmental Analysis in Biology and Medicine, 3rd ed. BioMedware, Ann Arbor. Jorgensen, S.E., 1995. The growth rate of zooplankton at the edge of chaos: ecological models. J. Theor. Biol. 175, 13–21. Jorgensen, S.E., Halling-Sorensen, B., Nielsen, S.N., 1996. Handbook of Environmental and Ecological Modelling. CRC Press, New York. Ljung, L., 1999. System Identification: Theory for the User, 2nd ed. Prentice-Hall, New Jersey. Macheras, P., Iliadis, A., 2005. Modeling in Biopharmaceutics, Pharmacokinetics and Pharmacodynamics, Homogeneous and Heterogeneous Approaches. Spinger-Verlag, New York. Matis, J.H., Kiffe, T.R., 2000. Stochastic Population Models. A Compartmental Perspective. Springer-Verlag, New York. Matis, J.H., Patten, B.C., White, G.C., 1979. Compartmental Analysis of Ecosystem Models. International Co-operative Publishing House, Burtonsville. Matis, J.H., Tolley, H.D., 1979. Compartmental models with multiple sources of stochastic variability: the one compartment, time invariant hazard rate case. B. Math. Biol. 41, 491–515. Matis, J.H., Wehrly, T.E., 1979. Stochastic models of compartmental systems. Biometrics 35, 199–220. MATLAB, 2004. High-performance Numeric Computation and Visualization Software. The Math Works, Natick, MA. Mosekilde, E., 1996. Topics in Nonlinear Dynamics: Applications to Physics, Biology and Economic Systems. World Scientific. Norton, J.P., 1986. An Introduction to Identification. Academic Press, London. Platt, T., Denman, K., 1977. Organisation in the pelagic ecosystem. Helgolander Meeresun 30, 575–581. Platt, T., Denman, K., 1978. The structure of pelagic marine ecosystems. Rapp. P-v Reun. Cons. Int. Explor. Mer. 173, 60–65. Ray, S., Berec, L., Strakraba, M., Jorgensen, S.E., 2001. Optimization of energy and implications of body sizes of phytoplankton and zooplankton in an aquatic ecosystem model. Ecol. Model. 140, 219–234. Riandey, V., Champalbert, G., Carlotti, F., Taupier-Letage, I., Thibault-Botha, D., 2005. Zooplankton distribution related to the hydrodynamic features in the Algerian Basin (western Mediterranean Sea) in summer 1997. Deep-Sea Res. 52, 2029–2048. Rodriguez, J., Mullin, M.M., 1986. Relation between biomass and body weight of plankton in a steady-state oceanic ecosystem. Limnol. Oceanogr. 31, 361–370.

378

e c o l o g i c a l m o d e l l i n g 2 0 4 ( 2 0 0 7 ) 372–378

Seber, G.A.F., Wild, C.J., 1989. Nonlinear Regression Analysis. John Wiley, New York. Sheldon, R.W., Prakask, A., Sutcliffe, W.H., 1972. The size distribution of particles in the ocean. Limnol. Oceanogr. 17, 327–340. Sourisseau, M., Carlotti, F., 2006. Spatial distribution of zooplankton size spectra on the French continental shelf of

the Bay of Biscay during spring 2000 and 2001. J. Geophys. Res. 111, 12. Zhou, M., 2006. What determines the slope of a plankton biomass spectrum? J. Plankton Res. 25, 437–448. Zhou, M., Huntley, M.E., 1997. Population dynamics theory of plankton based on biomass spectra. Mar. Ecol. Prog. Ser. 159, 61–73.