## Abstract

A dynamic energy budget (DEB) model for microalgae is proposed. This model deviates from the standard DEB model as it needs more reserves to cope with the variation of assimilation pathways, requiring a different approach to growth based on the synthesizing unit (SU) theory for multiple substrates. It is shown that the model is able to accurately predict experimental data in constant and light-varying conditions with most of the parameter values taken directly from the literature. Also, model simulations are shown to be consistent with stylized facts (SFs) concerning N∶C ratio. These SFs are reinterpreted and the general conclusion is that all forcing variables (dilution rate, temperature and irradiance) impose changes in the nitrogen or carbon limitation status of the population, and consequently on reserve densities. Model predictions are also evaluated in comparison with SFs on chlorophyll concentration. It is proposed that an extra structure, more dependent on the nitrogen reserve, is required to accurately model chlorophyll dynamics. Finally, SFs concerning extracellular polymeric substances (EPSs) production by benthic diatoms are collected and interpreted and a formulation based on product synthesis and rejection flux is proposed for the EPSs production rate.

## 1. Introduction

A common approach in microalgae research is to obtain and generalize empirical patterns through experimental work. These patterns, which in this paper will be referred as stylized facts (SFs), relate some independent variable, such as irradiance or dilution rate and some output, typically carbon or chlorophyll. However, despite their importance in microalgae research, few authors give insight on the mechanics behind these SFs. In this paper, a model based in dynamic energy budget (DEB) theory is proposed to demonstrate and give an explanation based on physical and chemical principles on several SFs found in the literature (table 1).

Dynamic energy budget theory provides a theoretical framework (Sousa *et al.* 2008; Kooijman 2009; reference Editorial of this issue) to model metabolism based in shared biological characteristics. The existence of such a theoretical framework is suggested by some empirical laws, such as von Bertalanffy's growth curve and Kleiber's law on metabolic rate, which apply to a great number of organisms (Sousa *et al.* 2008). It is based on a set of simple physical and chemical rules to describe the uptake and the use of energy and nutrients and how these are related to the organism's physiology through its life cycle (Kooijman 2009). The proposed DEB model is based on the multivariate model presented by Kooijman (2009), which already allows multiple and simultaneous substrate limitation, but extends it by also including photosynthesis, as the carbon assimilation pathway and providing time-dependent results. The model is parametrized and validated using data and parameter values from the literature.

The SFs in table 1 are also used to give insights into some open questions. In this paper, it is shown how SFs can be used to extract conclusions on the existence of a chlorophyll-specific variable. Chlorophyll is typically used as a proxy for microalgae biomass, but this approximation is known to be rather rough since there is not a linear relationship between these variables. Furthermore, the proposed model is used to interpret SFs concerning extracellular polymeric substances (EPSs) production by benthic diatoms. EPSs synthesis rate has been found to be dependent on growth conditions (Underwood *et al.* 2004)—a fact that can be further understood using DEB theory.

The present paper is structured as follows. In §2, we detail the developed model, with focus on photosynthesis, reserve dynamics, growth and temperature effects. In §3, parameter values are presented and model simulations are compared with experimental data from dynamic settings. Furthermore, it is presented an explanation of the SFs based on the proposed model, some insight into the behaviour of chlorophyll as a model variable is given, and EPSs production is related to microalgal metabolism. Finally, §4 finishes the paper with a summary of the conclusions of the presented work.

## 2. Model development

There are several features that differentiate the microalgae model from the standard DEB model. These are the need to consider multiple reserves and multiple substrates, the modelling of microalgae as V1-morphs, i.e. organisms, which have a surface area proportional to volume, and, precisely for this reason, the modelling of the whole population as a single organism with state variables that are the sum of the state variables of all the individuals.

The need to consider multiple reserves arises from the strong homeostasis assumption—the chemical composition of a reserve or a structure does not change during the organism's life cycle. When the assimilation pathways are essentially independent, such as in microalgae and other autotrophs, the constant chemical composition of reserve can only be ensured if there are different reserves for each independent assimilation pathway. In other words, strong homeostasis implies a reserve for each independent assimilation pathway. We considered only two nutrients, carbon and nitrogen, each with a corresponding assimilation pathway and reserve. Further reserves could be added to model other nutrients, such as phosphorous, and their dynamic would be similar to that of nitrogen. For simplicity, only nitrogen is considered, which is equivalent to assuming that all other nutrients are non-limiting.

The nitrogen case is modelled following the standard assimilation formulation, as presented in §2*a*. It is well documented that microalgae can feed from various forms of the same nutrient (e.g. NH_{3} and NO_{3}), but for the sake of simplicity, we considered only one general nutrient source for nitrogen. The assimilation of carbon through the photosynthesis pathway, common to most autotrophs, requires a different formulation, one based in multiple substrates. This mechanism—described by Saito *et al.* (2008) as type I colimitation—requires the presence of two substrates, photosynthetic radiation and CO_{2} (photosynthetic radiation is not treated as other substrates, since its availability is not determined by mass balance equations).

Physiological processes are either proportional to surface area, e.g. uptake of nutrients or proportional to volume, e.g. maintenance costs. Thus, the relationship between surface area and volume controls the metabolism of an organism (Hein *et al.* 1995). Organisms that propagate through division, such as microalgae, are treated by DEB theory as V1-morphs (Kooijman 2009), given that their small size range makes the differences between an isomorphic organism and a V1-morph negligible, from the population point of view. In the case of V1-morphs, the population is modelled as a single individual, which has volume (surface) equal to the sum of the volume (surface) of all individuals. The remaining section describes processes for an individual cell as function of structural mass *M*_{V} (which is related to volume and therefore to surface area). These descriptions are then directly used for the population.

Figure 1 represents the overall model with its main compartments, synthesizing units (SUs) (Kooijman 2009) and fluxes.

### (a) Assimilation and photosynthesis

In standard DEB theory, the assimilation of nutrients follows a Monod-type function. Thus, the formulation for the specific assimilation rate of nitrogen is:
2.1
where *j*_{ENAm} is the maximum specific nitrogen assimilation rate. Throughout the paper we use the standard DEB notation where *j* stands for a specific rate, i.e. a rate per C-mol of structure. Moreover, [N] stands for the nitrogen source concentration and *K*_{N} is the half-saturation concentration, i.e. the concentration that supports half of the maximum flux.

However, for an assimilation of carbon to exist one needs both a carbon source (CO_{2}) and photons. Therefore, we used the concepts of complementary assimilation of multiple substrates to build a photosynthesis model. The photosynthetic radiation is caught by photosynthetic units (PSUs) by excitation of electrons that fuel the production of NADPH^{+} and ATP, the fuel batteries that supply energy to the Calvin–Benson cycle. We model these reactions as a single-step reaction, enabling the computation of the relation between photon flux density (PFD) and photosynthesis rate, which is usually described in a photosynthesis–irradiance curve (PI-curve). To model photosynthesis with photoinhibition, the approach presented in Wu & Merchuk (2001)—see figure 2—is followed.

The specific relaxation rate *j*_{ph} in steady-state is given by equation (2.2) (see appendix for derivation).
2.2

Where *α*, *β*, *γ* and *δ* are photosynthetic coefficients that are described in figure 2, and *ρ*_{PSU} is the PSU density, i.e. the number of PSUs per unit of structural mass *M*_{V}, which is assumed to be a constant parameter. This formulation is consistent to those presented by Eilers & Peeters (1993), Megard *et al.* (1984), Zonneveld (1998*b*), although they take equation (2.2) to be the photosynthesis rate. In fact, the underlying assumption of these authors is that carbon dioxide is not limiting the photosynthesis rate. We combine the effects of CO_{2} concentration and irradiance using SUs theory for complementary and parallel substrates Poggiale *et al.* (2010). The resulting specific carbon assimilation flux *j*_{ECA} is given by
2.3
where *j*_{ECAm} is the maximum specific carbon assimilation rate. The carbon assimilation flux is more often called the photosynthesis rate, and we will also use this denomination. Finally, the specific CO_{2} uptake flux *j*_{CO2} follows also a Monod-type formulation,
2.4

### (b) Reserve dynamics

Microalgae use the energy and mass stored in reserves to fuel growth and maintenance processes, including costs with maturity maintenance. The specific catabolic flux *j*_{EiC} (with *i* = *C*, *N* for either carbon or nitrogen) represents the mobilized reserve per unit of structural mass that will be used for the aforementioned processes, with maintenance taking priority over growth, i.e. growth can only occur after maintenance requirements are fulfilled. From several assumptions of DEB theory (Kooijman 2009), it follows that reserves have first-order dynamics for the reserve density and that the structure specific catabolic flux is given by
2.5
where *ṙ* is the specific growth rate
2.6
*m*_{Ei} is the reserve density of reserve *i*, in C-mole of reserve per C-mole of structure, and *k̇*_{E} is the reserve turnover rate, which is constant. Moreover *k̇*_{E} is taken to be the same for the various reserves, since data on reserve dynamics in algae (Kooijman 2009) have shown that reserve turnover rates are equal even when the assimilatory pathways are not coupled (Kooijman *et al.* 2003). From equation (2.5) one can conclude that the mobilization of reserves is (i) independent of food availability, which providing the organism with an increased resistance against food availability fluctuations and an increased control over metabolism (Sousa *et al.* 2008), (ii) and supply driven, i.e. an increase (decrease) in the size of the reserve leads to an increase (decrease) in catabolic flux.

From figure 1, the dynamics of the reserves can be calculated from a mass balance, as follows:
2.7
where the last term corresponds to dilution by structural growth (with structural growth the density decreases, even if the reserve mass remains constant). The term *κ*_{Ei}*j*_{EiR} corresponds to a reabsorption of a fraction of the rejected flux that will be explained in §2*c*.

As stated before, the catabolic flux will be used first to pay maintenance requirements. Maintenance can be understood as the group of processes necessary for the survival of an organism, including structure turnover, maintaining constant gradient concentrations, among others. As in the standard DEB model, maintenance is taken to be proportional to structural volume and *j*_{EiM} is the constant specific maintenance cost paid by reserve *i*. If maintenance requirements cannot be fulfilled with the catabolic flux, a great number of organisms use structure as a mass and energy source. In the proposed model, we assume that the organism has the capacity to use structure, which may result in shrinking (reduction in structural mass), as explained in §2*c*.

### (c) Growth and rejection

The fraction of catabolic flux that is not used for maintenance is the specific growth flux *j*_{EiG},
2.8

The growth-SU merges all specific growth fluxes *j*_{EiG} stoichiometrically to synthesize a C-mole of structure, but these stoichiometric constraints imply that the growth-SU rejects some of the arriving molecules at a rate *j*_{EiR}. The resulting flux of synthesized structure is represented by the specific growth flux *j*_{G}. Following SU theory for parallel and complementary substrates (Kooijman 2009), we can derive an equation for the specific growth flux *j*_{G}:
2.9
where *y*_{EiG} is the yield coefficient that represents the number of moles of the aggregated compound of reserve *i* required to synthesize one C-mole of structure. Equation (2.9) neglects the synthesizing time assuming that this is irrelevant when compared with the binding time (Kooijman 2009). Moreover, when *j*_{EiC} is not sufficient to pay the maintenance costs, the model allows the remainder to be paid from structure at a rate *j*_{VM}. From the mass balance, the equation for the specific growth rate *ṙ* is given by
2.10

The inclusion of the *j*_{VM} term allows the specific growth rate *ṙ* to be negative. This term is composed by two fluxes *j*_{VM} = ∑_{i}^{C,N} *j*_{VMi}. Following the Switch Model proposed by Kooijman (2009), each of the fluxes *j*_{VMi} is given by
2.11

In population dynamics this term is extremely important, since negative specific growth rates are as common as positive. Still, this will only occur when the maintenance flux paid by the reserves is not sufficient. Using equations (2.5), (2.8) and (2.9) rewrite equation (2.10) as 2.12

If one or more growth fluxes *j*_{EiG} are limiting the growth rate, the non-limiting growth fluxes will occupy the respective SU's binding sites without being promptly processed. As the available binding sites decrease, the probability of these non-limiting growth fluxes *j*_{EiG} being rejected increases. The specific rejected flux *j*_{EiR} is thus given by:
2.13
where the first term represents the flux available for growth, while the last represents the structural growth, i.e. the flux that was actually used for growth. Substituting *j*_{EiG} in equation (2.13) with equations (2.5) and (2.8), we obtain:
2.14

Each rejected reserve molecule is fed back to the reserves with probability *κ*_{Ei} or excreted to the environment with probability (1 − *κ*_{Ei}), hence contributing to the dynamics of reserves given by equation (2.7). If one would not consider excretion, i.e. *κ*_{Ei} = 1, an organism limited by one of the nutrients would accumulate the non-limiting nutrients without any restriction, leading to impossibly large amounts of these reserves.

### (d) Product synthesis

In DEB theory all mass and energy fluxes result from a linear combination of three types of fluxes: assimilation, dissipation and growth. For the proposed model, with two reserves and one structure, and assuming that the catabolic flux is sufficient to pay the maintenance costs (i.e. *j*_{VM} = 0), the specific synthesis rate of product *i* is given by
2.15
where *y*_{Pi*} are coefficients that couples product i formation rate to mass flux *j*_{*}.

### (e) Temperature

Temperature largely affects metabolism by changing the rates of physiological processes. Zonneveld (1998*a*) points out that differences in biochemical cell composition arise as various processes are affected by temperature in different ways. The temperature dependence can be described with the Arrhenius relationship, within given temperature ranges, as follows
2.16
where *T*_{A} is the Arrhenius temperature, *T*_{0} a chosen reference temperature and *k̇*_{0} a general physiological rate at temperature *T*_{0}. Outside a given temperature range, the Arrhenius relationship cannot be applied as physiological rates start to behave too differently and the organism is no longer able to coordinate processes (Kooijman 2009).

It is important to stress that the Arrhenius relationship by definition can only be applied to parameters describing rates. In this case, this includes maximum specific assimilation rates, the specific maintenance costs and the reserve turnover rate. These ultimately produce a nonlinear effect in more complex rates (such as the specific growth rate), since not all processes are affected equally by temperature change. For example, photon capture by the photosystems is not affected by temperature, contrary to CO_{2} binding by RuBisCO, thus the overall carbon assimilation should not be directly corrected by the Arrhenius formulation. Thus, it was assumed that parameters *α*, *β*, *γ* and *δ* are temperature independent and that the Arrhenius relationship is not applicable.

### (f) Mass balance equations

To finalize the description of the proposed model, we present the mass model equations for substrate *X*_{i} and excreted *X*_{i}^{e}, reserve density *m*_{Ei} and structural mass *M*_{V}. A more detailed explanation for these equations is given in appendix A.
where the variable *ḣ* is the dilution rate of the medium and *X*_{Ii} is the concentration of substrate *i* in the input medium.

## 3. Results and discussion

### (a) Comparison with experimental data

Model simulations were compared with data from chemostat experiments for *Thalassiosira weissflogii* growing in nitrogen and light-limiting conditions (Pawlowski 2004). The experimental setting is as follows.

The *T. weissflogii* population grew during 4 days in batch conditions with constant irradiance (*I* = 250 µEm^{−2}s^{−1}). After the fourth day, dilution (*h* = 0.4 d^{−1}) and a light/dark circadian cycle (*I* = 180exp[3(cos(2*π*(*t* + 0.4)) − 1)] µEm^{−2} s^{−1}) were activated. Twenty-four days after the beginning of the experiment, irradiance was again forced constant (*I* = 100 µEm^{−2}s^{−1}). Nitrate concentration in the input was *X*_{Nr} = 70 µM and CO_{2} concentration was assumed to be *X*_{Cr} = 1 mM (corresponding to *p*CO_{2} ∼ 350 µatm). This last value had to be assumed since no information was given. Also, no information about initial biomass concentrations was found, thus these had to be assumed (*M*_{V} = 0.1 mM; *m*_{EN} = 0.01; *m*_{EC} = 0.1). The assumptions for the initial conditions had no impact on end results, given that the first 4 days in batch allowed the system to forget initial conditions. After 20 days from the beginning of the experiments, measurements of particulate carbon and particulate nitrogen were made. These correspond to the sum of carbon or nitrogen in the structure and reserves.

The presented model simulations were made using parameter values obtained from the literature (table 2), none of which was DEB specific. Also, most parameters found in the literature are given per C-moles of biomass, different from the presented approach in which these parameters are treated as per C-mole of structure. When possible, parameter values specific for *T. weissflogii* were used.

Comparisons of simulation results and experimental data are presented in figure 3. These show a good fit between empirical and predicted outputs, even more so if one considers that parameter values were directly obtained from the literature.

Particulate carbon, which in the proposed model corresponds to the sum of carbon in structure and C-reserve, shows a negative trend and decreasing amplitude that is not captured by the model. This may be caused by insufficient waiting time for the population to reach a periodic steady state or by the agglomeration of biofilm in the container's walls that would reduce outgoing biomass. Nevertheless, model predictions are in phase with the experimental data and show approximately the same amplitude. This cyclic behaviour is due not only to structure net growth but also mainly to carbon accumulation in C-reserve. The existence of this carbon reserve allows a lag between the carbon concentration peaks and irradiance. To understand it we refer to equation (2.7). During irradiance increase, *j*_{ECA} is bolstered causing the increase of the C-reserve density; after the irradiance peak, *j*_{ECA} steadily decreases with irradiance but is still higher than the increasing sum of *j*_{ECC} and *ṙ**m*_{EC} (which are proportional to reserve density) thus reserve density still increases; finally, when *j*_{ECA} is no longer capable of supporting the high values of *j*_{ECC} plus *ṙ**m*_{EC}, the C-reserve density reaches its peak and then starts to decrease. The lag corresponds to the time between maximum *j*_{ECA} and *j*_{ECA} = *j*_{ECC} + *ṙ**m*_{EC} − *κ*_{EC}*j*_{ECR} (i.e. (d/d*t*)*m*_{Ei} = 0). The last term *κ*_{EC}*j*_{ECR} was neglected from discussion analysis since it is nearly constant j.

After day 24, when irradiance is forced constant, model predictions continue to behave in accordance with experimental data. An exponential growth phase is observed approximately until day 26, after which a balanced growth phase is reached.

Concerning particulate nitrogen concentration, the proposed model predicts a near constant and slightly overestimated value and does not fully capture the negative trend, also visible in the carbon data. This near constancy is somewhat expected since N-reserve density is not directly affected by irradiance but indirectly through C-reserve density (which decreases the rejection flux *j*_{ENR} and affects the specific growth rate *ṙ*). It was observed that structure closely follows C-reserve density, hinting that nitrogen was hardly limiting. This hypothesis was verified by comparing the normalized growth flux *j*_{EiG}/*y*_{EiV} for the two reserves. The conclusion was that the C-reserve was more limiting even in the peak of irradiance.

Overall, and considering that parameter values did not result from any statistical optimization procedure, it is commendable that model simulations are so similar to experimental results.

### (b) Comparison with stylized facts

Several SFs, concerning chemical composition and growth rate, were collected from the literature (table 1) and compared with model outputs.

#### (i) N∶C ratio

The most common observed variable when assessing variation in the chemical composition of microalgae is its N∶C ratio. It was observed that the N∶C ratio increases with dilution rate when nutrient limited (Laws & Bannister 1980; Chalup & Laws 1990; Geider *et al.* 1998) and decreases when light limited (Laws & Bannister 1980). Furthermore, the N∶C ratio decreases with temperature when nutrient or light limited and with irradiance. These four SFs were simulated with our model (figure 4). The N∶C ratio can be computed using the following equation
3.1
where *n*_{NV} is the elemental coefficient of nitrogen in structure and *M*_{Ei} stands for the total mass of reserve *i*. Equation (3.1) shows that to assess N∶C ratio one should look to the dependence of reserve densities on forcing variables.

The first curve (figure 4*a*) shows that the N∶C ratio increases with dilution rate in nutrient limiting (in this case nitrogen limiting) conditions, as predicted by SF1. This happens not only because of the increase of N-reserve density (which obviously is needed to support higher specific growth rates) but also because C-reserve density decreases. The first point is simply owing to the increase of nitrogen source concentration and the consequent increase of *j*_{ENA}. Furthermore, higher dilution rates reduce the rejected flux for the non-limiting reserve—in this case C-reserve (equation (2.14))—which then leads to a decrease in reserve density (equation (2.7)) since carbon assimilation is constant. The same pattern is not observed in the N-reserve because *j*_{ENA} increases with dilution rate, contrary to the *j*_{ECA}, which is practically constant (since CO_{2} is non-limiting and light is forced constant).

SF2—N∶C ratio decreases with dilution rate when light limited (Laws & Bannister 1980)—is harder to simulate and understand, mainly because of the definition of light limited. If a population is grown with dilution rate *ḣ* in light-limiting conditions there are only two possible outcomes: (i) irradiance is enough to support *ṙ* > *ḣ* and, since *j*_{ph} is constant per unit of *M*_{V} (equation (2.2)), biomass will grow and steady state will only be attained when one mineral (carbon or nitrogen) limits growth; (ii) irradiance is not enough to support *ṙ* > *ḣ* and the population will only reach steady state at *M*_{V} = 0. This means that in a chemostat light can only co-limit. The alternative, not contemplated in the proposed model, is that cell-shading effects can decrease irradiance per unit of *M*_{V} in a similar way that carbon and nitrogen sources decrease with increasing biomass.

We reformulated the SF so that it includes carbon limitation—N∶C ratio decreases with growth when light or carbon limited. In this case, it can be shown that its rationale is similar to curve in figure 4*a*. When dilution rate is increased, CO_{2} available to microalgae also increases, which partially compensates the constant low-light conditions (equation (2.3)) thus increasing C-reserve density. On the other hand, higher dilution rates (which in steady state implies higher specific growth rates) reduce the rejection flux corresponding to N-reserve and, since *j*_{ENA} is already at its maximum (nitrogen is non-limiting), the N-reserve density decreases with dilution rate. Owing to the increase in C-reserve density and the decrease in N-reserve density, it is clear that the N∶C ratio decreases with dilution rate in light or carbon-limiting conditions. However, it must be noted that the N∶C values for this simulation are outside the typical range of values owing to the excessive accumulation of nitrogen. This suggests that either the fraction of rejection flux incorporated in N-reserve *κ*_{EN} is too high or that there are other relevant processes not considered in the proposed model (e.g. energy costs for nitrogen assimilation).

SF3—N∶C ratio decreases with temperature in nutrient and light-limiting conditions—is observed in figure 4*c* and can be related to an increasing limitation of the N-reserve. Since higher temperatures lead to higher rates, including assimilation, maintenance and catabolic rates, it is expected that the nitrogen limitation becomes more important with temperature. While all rates increase with temperature, the assimilation fluxes are still restricted to environmental concentrations—thus, in nitrogen-limiting conditions, *j*_{ENA} increases much less than other rates and the N-reserve density decreases with temperature. Our model also predicts that total nitrogen (i.e. nitrogen in reserve and in structure) remains near constant (not shown), a result observed experimentally by Berges *et al.* (2002). The C-reserve density increases with temperature owing to the increase of *j*_{ECAm} and *j*_{CO2} and of the rejection flux *j*_{ECR} until approximately *T* = 295 K, after which C-reserve density steadily declines since *j*_{ECA} increases less than other relevant rates, particularly *j*_{ECC}, since temperature does not affect photon uptake (see §2*e*). This is not in agreement with Berges *et al.* (2002) statements, but does not contradict his experimental data, which had a more limited temperature range (approximately from 280 to 300 K). The N∶C ratio that results from these dynamics has a minimum at *T* = 300 K.

Finally, SF4—N∶C ratio decreases with irradiance in nutrient replete conditions—is represented in figure 4*d*. The rationale is similar to the one for curve 4*b*. The microalgal population will grow if irradiance is such that *ṙ* > *ḣ*. If such condition is fulfilled then the population will grow until another factor is limiting, and only then it will reach a steady state. Figure 4*d* presents the N∶C ratio for a simulation in which nitrogen becomes limiting for all values of irradiance considered (typical environmental setting). Irradiance bolsters *j*_{ECA} and thus the C-reserve density. The consequent increase of *M*_{V} and decrease of *j*_{ENR} reduce N-reserve density. The overall result is a decrease of the N∶C ratio.

#### (ii) Chlorophyll concentration

SFs state that chlorophyll concentration decreases with irradiance but increases with the dilution rate (SF7 and SF8). Until this point no explicit considerations were made concerning chlorophyll; it was implicitly assumed that it was a part of structure (since photon capture *j*_{ph} is assumed to be proportional to structural mass *M*_{V}). Keeping with that assumption, it is possible to analyse chlorophyll concentration if one considers it as proportional to structural mass *M*_{V}.

Values for structural mass and reserve density are presented in figure 5 for a simulation in which the effect of increasing irradiance is tested. After day 20 the irradiance is forced to be twice the previous value. The expected effect in C-reserve is confirmed—it increases until reaching a new steady-state. Structural mass follows the same behaviour since it is also supported by a greater C-reserve. Finally, N-reserve density does not get bolstered by increased irradiance but needs to support a larger structure, thus stabilizing at a lower value.

Since it is observed that chlorophyll concentration decreases with irradiance (SF8), our model results suggest that chlorophyll can not be solely in structure or in C-reserve (both increase with increasing irradiance). If this was verified then structure or C-reserve would have varying chemical composition thus breaking the strong homeostasis assumption. This leaves the hypothesis that a relevant fraction of chlorophyll be a part of the N-reserve, which is indirectly supported by the SF9—chl∶N ratio decreases with irradiance for constant dilution rate in nutrient replete conditions (which is equivalent to say that N-reserve density decreases with irradiance—see discussion on SF4). Similar conclusion can be drawn from SF5 and SF6. These two SFs are similar to SF1 and SF2, which suggest that the ratio chl∶C behaves similarly to the ratio N∶C in those conditions (figure 4*a*,*b*), i.e. chlorophyll increases (decreases) when N-reserve density increases (decreases). However, chlorophyll being part of N-reserve would imply that microalgae were too susceptible to environmental nitrogen variations.

A better hypothesis is to consider an additional structure, which is more dependent on the N-reserve than on the C-reserve. This hypothesis still verifies the mentioned SFs and is coherent from an evolutionary perspective. In fact, many phytoplankton models follow this idea of attributing a state variable to chlorophyll (Laws & Chalup 1990; Geider *et al.* 1997, 1998; Baklouti *et al.* 2006). Typically the assimilation of carbon is made to a carbon compartment while nitrogen is assimilated to chlorophyll, thus chemical variations arise as the assimilation fluxes vary.

If one considered this chlorophyll structure, represented by *M*_{chl}, the specific deactivation rate *j*_{ph} would be linked to this structure. Then it would be possible to model photoacclimation—the process by which cells adjust chlorophyll content in accordance to irradiance (Zonneveld 1998*b*; Macintyre *et al.* 2002; Papadakis *et al.* 2005)—by establishing a growth-SU for the chlorophyll structure in which light directly or indirectly inhibits growth. Alternatively, (Papadakis *et al.* 2005) proposed varying parameter values to cope with the variation of chlorophyll content by photoacclimation. In the proposed model, this could be accomplished by making PSU density *ρ*_{PSU} and the activation coefficient *α* functions of irradiance at which cells were acclimated.

The addition of a chlorophyll structure would obviously increase model precision at the cost of additional parameters. The variation of chlorophyll content is rather important in natural setting and in large timescales, but for relatively constant or cyclic irradiances the chlorophyll concentration varies little from an average value (e.g. Owens *et al.* (1980), Geider *et al.* (1998) and Pawlowski (2004) show variations of approx. ±15% in controlled settings, while for natural settings works by Fuhrman *et al.* (1985) and Hitchcock (1980) report variations of approx. ±30%), and thus can be assumed to be proportional to the general structure, as was done in the proposed model.

#### (iii) EPSs production

The major issue when considering EPSs-related SFs is the lack of agreement between studies (e.g. Staats *et al.* (2000) found that EPSs production rate is null in the darkness while Smith & Graham Underwood (1998) state that it may even increase). Underwood *et al.* (2004) made an effort to clarify the effect of environmental aspects in production rates, production pathways and chemical composition of EPSs. Notwithstanding its merits, the work lacks an unified approach, presented by those such as Laspidou (2002) for microbial EPSs production. The formulation of an unified theory for microalgal EPSs production is outside the scope of this work; our aim is to discuss how SFs can be understood in the context of DEB theory. These SFs are based on the works of Smith & Graham Underwood (1998) and Underwood *et al.* (2004) for EPSs production by benthic diatoms.

To understand how EPSs production is considered in DEB theory, one must consider the two possible contributors for a production rate: product synthesis *j*_{P} (equation (2.15)) and rejection flux *j*_{EiR} (equation (2.14)). Thus *j*_{EPS} can be given by
3.2
where *θ* represents the fraction of rejected flux *j*_{ECR} that is considered to be EPSs. If one takes into account that microalgal EPSs are mostly polysaccharides (Underwood *et al.* 2004), its dynamics should be closely related to the C-reserve dynamics. It is then reasonable to assume that product synthesis rate *j*_{P} is given by equation (2.15):
3.3
where *y*_{P*} couples EPSs production to the flux represented by asterisk, i.e. assimilation, maintenance or growth.

We use equations (3.2) and (3.3) to discuss SF10, SF11 and SF12. SF10—no EPSs production results directly from light-dependent carbon assimilation—clearly relates the EPSs production rate *j*_{EPS} to the carbon assimilation rate *j*_{ECA}. The lack of correlation immediately forces the coefficient *y*_{PECA} to be null. However, photosynthesis indirectly affects *j*_{EPS} through reserve dynamics—if *j*_{ECA} is low or null (e.g. dark conditions) then C-reserve density will decrease implying that all other possible sources of EPSs also decrease.

SF11—cellular production rates of EPSs and other extracellular carbohydrates increase with specific growth rate—can be interpreted in the context of DEB theory as the increase of *j*_{EPS} with specific growth rate. This fits perfectly with equation (3.3), in which *j*_{P} has a component proportional to the specific growth rate represented by *j*_{G}. The proportionality is given by *y*_{PG}, which has to be positive.

Finally, SF12—a different type of EPSs (i.e. different chemical composition) is produced under nutrient-limiting conditions—identifies a new type of EPSs related to nutrient-limiting conditions. Underwood *et al.* (2004) differentiate this EPSs as type 2, in opposition to the EPSs type 1, which is produced under all conditions. From the several possibilities—growth, maintenance and rejection—only one is directly related to nutrient-limiting conditions, and that is the rejection flux *j*_{ECR}.

Until this point, no relation between maintenance and EPSs production was considered. The fact that EPSs type 1 is produced under a wide range of conditions, even when growth is almost null, suggests that this type of EPSs is also related to maintenance. Furthermore, Underwood *et al.* (2004) suggest that EPSs production is also related to motility in some diatoms, a process associated with maintenance in the context of DEB theory.

Relating all these statements to equations (3.4) and (3.3), one obtains the following formulation for *j*_{EPS}
3.4

Rewriting equation (3.4) so that structural mass *M*_{V} becomes explicit, we have
3.5
where *j*_{G} was replaced by *ṙ*. This formulation is similar to the ones presented by Characklis & Marshall (1990) and Hsieh *et al.* (1994) for microbial EPSs production. Both formulations have a term for growth-related EPSs production with a general form *ṙ**X* with *ṙ* dependent on substrate concentration, which in equation (3.5) is *y*_{PG}*ṙ**M*_{V}. The difference lies in the consideration of structural mass *M*_{V} as biomass *X* and in the formulation of *ṙ*, which in the proposed model is given by equation (2.12). Also, both articles consider a term with the general form *f*_{XP}*k*_{d}*X* in which *f*_{XP} is the fraction of biomass converted to EPSs and *k*_{d} is the biomass decay coefficient. Again, a similar term exists in equation (3.5): *y*_{PECM}*j*_{ECM} *M*_{V}. The constants *y*_{PECM} and *j*_{ECM} correspond, respectively, to *f*_{XP} and *k*_{d}, and, once again, structural mass *M*_{V} stands for biomass *X*. The last term in equation (3.5) is related to rejection and has no direct parallel in the formulations proposed by Characklis & Marshall (1990) and Hsieh *et al.* (1994), which was to be expected since microbial metabolism is more similar to the standard DEB model (which has a single reserve), which does not need to consider rejection fluxes. For microalgae this term is necessary to fit SF12—EPSs type 2 production increases in nutrient-limiting conditions. Hsieh *et al.* (1994) comply with this fact using an alternative approach to the decay term: this term is taken to be proportional to nutrient-limitation, i.e. when in nutrient-limiting conditions the decay-related fraction of EPSs production rate will increase.

## 4. Conclusions

A DEB-based model for microalgae was proposed. This model deviates from the standard DEB model as it needs more reserves to cope with the variation of assimilation pathways, implying a different approach to growth. It was shown that the model is able to accurately predict experimental data in light-varying conditions, especially if one considers that most parameter values were taken directly from the literature.

In this paper, we presented a DEB-based model for microalgae and used it to explain several patterns collected throughout the literature—the SFs. The proposed model deviates from the standard DEB model since it considers two reserves and therefore needs a different formulation for growth. This framework was already proposed by Kooijman (2009), but we took it a step further by including a carbon assimilation pathway dependent on light and by proposing a set of parameter values that can be used as a baseline for future parametrization works.

The parameter values taken from the literature (mostly for *T. weissfloggi*), some of which required a translation procedure to comply with the definitions of DEB theory, resulted in simulations which fit the data for total carbon and total nitrogen concentration. The fits were surprisingly good, given that there was no optimization procedure.

The model was also qualitatively validated with SF1 to SF4—model simulations were compared with the behaviour of the N∶C ratio in four cases. It was found that the model complies to all cases, with the partial exception of SF2. This SF is qualitatively verified, but model predictions are not consistent with actual values of the N∶C ratio. The model clearly overestimates the N∶C value, which may be owing to the fact that we did not consider some important process in the N-reserve dynamics, such as the energy costs for nitrogen assimilation (Flynn 1990).

The fact that reserve dynamics determine the N∶C ratio also allowed us to analyse and explain how these chemical variations are imposed by the forcing variables. The general conclusion from the analysis of SF1 to SF4 is that all independent variables change the limitation status of the population. The increase of dilution rate makes the limiting factor (N or C) less limiting, and consequently the respective reserve density increases (SF1 and SF2). On the other hand, temperature increase makes nutrient-limitation more severe that obviously leads to the decrease of N-reserve density (SF3). Finally, it was found that for constant irradiance per unit of structural mass, steady state can only be achieved when N or C become limiting (SF4). This forced us to reformulate SF2 to include C limitation and SF4 to include C or N limitation.

Model results were also compared with SFs on chlorophyll concentration. These comparisons allowed us to conclude that the proposed DEB model is not able to simulate the decrease in chlorophyll concentration with irradiance (SF8). Based on these findings, we conclude that a specific structure for chlorophyll, which would be strongly dependent on the N-reserve, should be included in a more detailed DEB model. This approach is consistent with other authors' proposals (Laws & Chalup 1990; Geider *et al.* 1997; Geider *et al.* 1998; Baklouti *et al.* 2006).

Finally, the proposed model was used to investigate SFs on EPSs production by benthic diatoms. It was found that all SFs can be interpreted considering that EPSs production is a combination of products and rejection fluxes, as defined by DEB theory. This means that EPSs production is a linear combination of maintenance, growth and rejection owing to stoichiometric constraints on growth and differences between species can be captured with species-specific coefficient values. These differences can explain seemingly contradictory SFs found in the literature, such as the findings by Staats *et al.* (2000) that found no EPSs production in the dark and thus concluded that oxygenic photosynthesis as a prerequisite for secretion of polysaccharides, while others such as Smith & Graham Underwood (1998) and Underwood *et al.* (2004) argue that there is not a direct connection between photosynthesis and EPSs production.

This paper introduced several interesting issues that have not yet been fully answered both in the context of DEB theory and microalgae modelling. One of them is the analysis of the importance of modelling chlorophyll as an independent structure in different spatial and time scales. This aspect should be linked with a more detailed study of the photosynthesis process, particularly of acclimation and photorespiration. An important result would be a rule of thumb that would relate chlorophyll concentration in natural environments (easy to measure) with phytoplankton biomass. Another important issue is the development of a coherent and inclusive framework for the production of not only EPSs but also all extracellular carbohydrates. This will require an exhaustive collection of empirical patterns and data to test existing and new models for the production extracellular carbohydrates. A third issue that should be completed is the collection and comparison of parameter values for microalgae species. This is simplified in the context of DEB theory since it gives a framework in which parameter values can be systematically categorized and compared for most microalgae species.

## Acknowledgements

The authors are grateful to Prof. J. C. Poggiale for his helpful comments and M. Baklouti for providing the experimental data presented. The authors would also like to thank three anonymous reviewers for their useful comments. This research was financially supported by FCT through Grant No. PPCDT/AMB/55701/2004. The work of G.M.M. was financed by FCT through grant SFRH/BPD/27174/2006.

## Appendix A. Mass balance equations

Besides the organism-related state variables, the model also includes environmental concentrations variables. Each reserve has a corresponding substrate and a rejected compound, the latter with the same chemical composition of the aggregrated compound of the reserve. The dynamics of the substrate of reserve *i* can be given by a mass balance as follows:
A 1
where *ḣ* is the dilution rate, *X*_{ir} is the concentration of substrate in the input, the yield factor *y*_{EiV} represents the moles of *i*-reserve required to synthesize 1 mole of structure and *n*_{iV} is the chemical index of element *i* in the structure compound. The first term represents the ingoing and the outgoing fluxes in the chemostat, the second term corresponds to the mass of substrate assimilated and, finally, the last term represents the subproducts that arise from the transformation of reserve compound into structure.

The quantity of structural mass *M*_{V} in the chemostat is simply given by the following mass balance
A 2
where the specific growth rate *ṙ* represents the generation of structural mass, while the dilution rate *ḣ* stands for the fraction of structural mass leaving the chemostat.

## Appendix B. Derivation of the equation for specific relaxation rate

Equation (2.2) can be obtained using the following system of differential equations: B 1

The solution for steady state (d*χ*_{1}/d*t* = 0) defines the rate at which a PSU transits from *χ*_{2} to *χ*_{1}, i.e. the quantity *χ*_{2}*γ*:
B 2

Finally to obtain the structure-specific deactivation rate *j*_{ph} we multiply (B2) by the PSU density *ρ*_{PSU}:
B 3

## Footnotes

One contribution of 14 to a Theme Issue ‘Developments in dynamic energy budget theory and its applications’.

- © 2010 The Royal Society