Examination of the role of the microbial loop in regulating lake nutrient stoichiometry and phytoplankton dynamics

The recycling of organic material through bacteria and microzooplankton to higher trophic levels, known as the “microbial loop”, is an important process in aquatic ecosystems. Here the significance of the microbial loop in influencing nutrient supply to phytoplankton has been investigated in Lake Kinneret (Israel) using a coupled hydrodynamic– ecosystem model. The model was designed to simulate the dynamic cycling of carbon, nitrogen and phosphorus through bacteria, phytoplankton and zooplankton functional groups, with each pool having unique C : N : P dynamics. Three microbial loop sub-model configurations were used to isolate mechanisms by which the microbial loop could influence phytoplankton biomass, considering (i) the role of bacterial mineralisation, (ii) the effect of micrograzer excretion, and (iii) bacterial ability to compete for dissolved inorganic nutrients. The nutrient flux pathways between the abiotic pools and biotic groups and the patterns of biomass and nutrient limitation of the different phytoplankton groups were quantified for the different model configurations. Considerable variation in phytoplankton biomass and dissolved organic matter demonstrated the sensitivity of predictions to assumptions about microbial loop operation and the specific mechanisms by which phytoplankton growth was affected. Comparison of the simulations identified that the microbial loop most significantly altered phytoplankton growth by periodically amplifying internal phosphorus limitation due to bacterial competition for phosphate to satisfy their own stoichiometric requirements. Importantly, each configuration led to a unique prediction of the overall community composition, and we conclude that the microbial loop plays an important role in nutrient recycling by regulating not only the quantity, but also the stoichiometry of available N and P that is available to primary producers. The results demonstrate how commonly employed simplifying assumptions about model structure can lead to large uncertainty in phytoplankton community predictions and highlight the need for aquatic ecosystem models to carefully resolve the variable stoichiometry dynamics of microbial interactions.


Introduction
One of the principal objectives for water quality management of freshwater bodies is to reduce the magnitude and frequency of nuisance algal blooms.Excess nutrients are generally implicated in the production of nuisance blooms since they fuel primary production and organic matter accumulation (Elser, 1999).In trying to understand these processes much work in limnology is based on the classic "N-P-Z-D" (nutrients-phytoplankton-zooplankton-detritus) Published by Copernicus Publications on behalf of the European Geosciences Union.

Y. Li et al.: Microbial loop effects on lake stoichiometry
paradigm, which assumes a relatively simple flow of nutrients to autotrophic and then heterotrophic pools.However, it is now well-documented both in oceanographic and, to a lesser extent, in limnological applications, that higher order predators such as crustacean zooplankton or fish can be supported by two paths: the so-called "green" (algal-based) and "brown" (detrital-based) food web components (Moore et al., 2004).The latter refers to the dynamics of the heterotrophic bacteria and the microzooplankton grazers (defined here as size less than 125 µm that account for rotifers, ciliates and juvenile macrograzers; Thatcher et al., 1993) -often termed the "microbial loop".This has been shown to play an important role in shaping carbon fluxes in lakes and in enhancing nutrient cycling at the base of food webs (Gaedke et al., 2002), including in Lake Kinneret which is the focus in this study (Stone et al., 1993;Hart et al., 2000;Hambright et al., 2007;Berman et al., 2010).
Less well understood is how the microbial loop affects phytoplankton growth and thus potentially shape patterns of phytoplankton succession.There are several main mechanisms by which microbial loop processes are thought to influence phytoplankton dynamics: (i) the provision of bacterially mineralised nutrients for phytoplankton growth; (ii) the excretion of readily available nutrients by micrograzers that support primary production (Johannes, 1965;Wang et al., 2009); (iii) the competition of bacteria with phytoplankton for inorganic nutrients when organic detritus becomes nutrient depleted (Barsdate et al., 1974;Bratbak and Thingstad, 1985;Stone, 1990;Kirchman, 1994;Caron, 1994;Joint et al., 2002;Danger et al., 2007).Additionally, a potential fourth indirect mechanism is that bacteria provide an alternative food source for micrograzers, thus alleviating some grazing pressure from small primary producers.The relative significance of each of these mechanisms, and in particular how they interact in a dynamic environment to shape microbial community composition and influence net productivity remains unclear.
Models of lake ecosystems are increasingly common to support management and analysis of water quality problems, acting as 'virtual' laboratories for exploring ecosystem processes particularly for questions where empirical studies would be difficult to undertake (Van Nes and Scheffer, 2005;Mooij et al., 2010).In most models published to date it is generally assumed that the biomass of heterotrophic bacteria is fairly stable and that the majority of bacterial production is lost to respiration (Cole, 1999).As a result, most quantitative models of carbon and nutrient fluxes in freshwater ecosystems essentially simplify microbial loop processes by assuming a relatively static mineralisation rate of organic material and simulating direct zooplankton consumption of detritus as a proxy for microzooplankton consumption of bacteria (e.g.Janse et al., 1992;Saito et al., 2001;Bruce et al., 2006;Mooij et al., 2010).These simplifications do not capture the range of nutrient 'adjustments' that occur during microbial loop processes, since stoichiometric composition of organisms and the fluxes between them are in reality neither uniform nor static (Elser and Urabe, 1999;Sterner and Elser, 2002).Whilst representation of microbial loop processes has been developed in marine ecosystem models (e.g.Faure et al., 2010), their uptake in freshwater ecosystem models has been limited and none to our knowledge simultaneously resolve the microbial loop and the dynamic stoichiometry of carbon (C), nitrogen (N) and phosphorus (P).
As a background to this study, there have been several attempts to incorporate the microbial loop into Lake Kinneret ecosystem models.Initially, a steady-state C flux model was developed to examine C cycling through the planktonic biota, including consideration of the microbial loop (Stone et al., 1993;Hart et al., 2000).A one-dimensional (1-D) coupled hydrodynamic-ecosystem model (DYRESM-CAEDYM) was presented by Bruce et al. (2006), which focused specifically on the zooplankton dynamics and their contribution to nutrient recycling.However, the model presented by Bruce et al. (2006) had a simplistic representation of the microbial loop dynamics, like many contemporary lake ecosystem models, and also did not individually simulate two important cyanobacterial species, Microcystis sp. and Aphanizomenon sp., which are important to the health of the ecosystem and sensitive to stoichiometric constraints within the food web (Zohary, 2004).Gal et al. (2009) expanded this model to include a dynamic microbial loop parameterisation and accounted for the two cyanobacterial species listed above and validated the model approach against a comprehensive data set.The relationship between phytoplankton stoichiometry and patterns in the stoichiometry of available nutrients was further analysed by Li et al. (2013), who noted that the microbial loop parameterisation approach could adjust both the quantity and stoichiometry of nutrient transfers.
This research builds on these studies and applies the validated model with the general aim of isolating the significance of the microbial loop on the phytoplankton patterns within the lake.Specifically, three different microbial loop model structural configurations were designed and analysed to unravel how the microbial loop processes identified above combine to affect (a) the stoichiometry of nutrient transfers through the planktonic food web, and (b) phytoplankton growth and community composition.The results highlight the importance of resolving the variable stoichiometry of microbial interactions in aquatic models and suggest that commonly used simplifying assumptions may compromise model function.

Site description
Lake Kinneret (Sea of Galilee) is a large monomictic lake located in the Syrian-African Rift Valley in northeastern Israel.It covers an area of 170 km 2 , is 21 km long and 16 km wide, has a maximum depth of 43 m, and has been the focus of considerable limnological research over the past few decades.Major phytoplankton groups present in the lake include Peridinium sp., Aulacoseira sp., Aphanizomenon sp., Microcystis sp. and nanophytoplankton.A number of zooplankton species occur in the lake and can be grouped as rotifers, ciliates, and herbivorous (cladocerans and copepodites) and predatory zooplankton (adult copepods).The maximum ciliate abundance is observed in autumn, generally preceding a metazooplankton peak.Heterotrophic nanoflagellates are most abundant in winter and spring, and least abundant in autumn.Bacteria numbers are highest during the decline of the Peridinium gatunense (hereafter referred to as Peridinium) bloom and are the lowest during the winter (Hadas et al., 1998).Lake Kinneret was once well known for seasonal blooms of Peridinium that regularly occurred until the late 1990s (Zohary et al., 1998;Zohary, 2004;Roelke et al., 2007).However, observations over the last decade have seen a major decline in Peridinium and a disruption in the historically stable patterns of phytoplankton succession (Zohary, 2004).In response, the biomass of Aulacoseira blooms has changed and the contribution of cyanobacteria and nanophytoplankton to the total phytoplankton biomass has increased in summer.Due to reduced water quality, the occurrence of nuisance cyanobacterial blooms is an increasing concern (Ballot et al., 2011).

Model overview and simulation approach
In this study, the 1-D hydrodynamic-ecological model DYRESM-CAEDYM was applied to the lake for the period from 1997 to 2001.The configuration adopted here is a continuation from Gal et al. (2009), with three microbial loop sub-model configurations applied (Fig. 1), as described below.The model simulates phytoplankton dynamics, bacterial production, carbon and nutrient recycling, sedimentwater interactions, and relevant inflow, outflow and mixing processes.In each simulation conducted, five phytoplankton groups, A, are included, each with three state variables (internal C, N, and P, denoted as A, AIN, and AIP, respectively): Peridinium (A 1 ), Microcystis (A 2 ), Aphanizomenon (A 3 ), nanophytoplankton (A 4 ) and Aulacoseira (A 5 ).Three zooplankton functional groups, Z, each with fixed internal nutrient ratios, were also simulated: predatory copepods (Z 1 ), macrograzers (Z 2 ) and microzooplankton (Z 3 ).Bacteria (B) were modelled as a separate state variable for two of the microbial loop configurations.An additional 10 nutrient variables (FRP, NO 3 , NH 4 , DIC, DOC, DON, DOP, POC, PON, POP), and dissolved oxygen (DO) were also modelled, giving a total of 30 biogeochemical state variables (Table 1).Field data was used to initialise the vertical profiles of all major state variables.
Earlier versions of the model (Gal et al., 2009, Makler-Pick et al., 2011a, b;Li et al., 2013) have been thoroughly validated against field data and process measurements (Ap-pendix A).Here we use the best-calibrated model version from these studies to explore the impact of changes in bacterial dynamics on patterns of phytoplankton growth.Within a well-documented set of core ecological process parameters determined elsewhere, we vary the structure and function of the microbial loop to assess how these changes would impact broader ecosystem biogeochemistry.Therefore, while this is essentially a theoretical study, it remains nested in a robust modelling framework built on a strong process understanding of the Lake Kinneret ecosystem.

Bacteria and microbial loop sub-models
Three alternative microbial loop sub-model configurations were compared to evaluate the relative importance of the three key mechanisms by which the microbial loop can affect phytoplankton dynamics (Table 2).Note that in this study we are not further considering the role of the fourth indirect mechanism listed in the introduction, since in Lake Kinneret the micrograzer food source is thought to be predominantly bacteria.The three simulations are differentiated by having (1) an assumed constant bacteria biomass state variable using static organic matter mineralisation rates and microzooplankton grazing directly on POM (NOBAC hereafter); (2) bacteria simulated with dynamic biomass and hence mineralisation rates, but unable to take up dissolved inorganic N and P (BAC − DIM hereafter); (3) dynamic bacteria (as per 2) with an additional ability for supplementing their internal nutrient requirement with dissolved inorganic N and P (PO 4 and NO 3 / NH 4 ) if the available organic matter becomes nutrient deplete (BAC + DIM hereafter).
The general mathematical description of the mass balance for each of the variables and associated notations are in Table 3.For each configuration, parameterisation of the common microbial loop process pathways are described in detail next and summarised in Table 4.For other CAEDYM variable descriptions, process equations and parameter values and justifications, readers are referred to Gal et al. (2009).

POM hydrolysis
This process considers the enzymatic hydrolysis and decomposition (D POM ) of particulate detrital material, limited by dissolved oxygen concentration (DO) and bacterial biomass (B) if bacteria are simulated: where µ POM max is the maximum transfer of POM to DOM, and refers to one of µ POC max , µ PON max , or µ POP max (Table 4).

DOM mineralisation
Whilst the mineralisation of DOM to DIM is common to all configurations, when the bacteria state variable is included the process adopts a two-stage breakdown pathway as shown in the subsequent details of the BAC − DIM and BAC + DIM configurations.The general rate of DOM breakdown/uptake (U DOC ) is simulated as

Micrograzer grazing
All simulations include microzooplankton (Z 3 ), which graze either on a lumped detrital pool (NOBAC) or directly on bacteria (BAC − DIM and BAC + DIM).For simplicity, microzooplankton are considered to graze on either bacteria or detritus, since the rate of grazing on small size phytoplankton (A 4 ) has been reported to be relatively low compared to the rate of bacterial grazing (∼ 10 % of total microzooplankton diet in Lake Kinneret; Hambright et al., 2007).

Micrograzer excretion and respiration
In all configurations micrograzers respire (R) and excrete (E) labile organic matter: where k Zr is the respiration rate and k Ze is the DOC excretion rate.Since micrograzers are configured to have a stable C:N:P requirement, their excretion of N and P is variable in order to balance the other output nutrient fluxes.This is numerically achieved by performing the excretion at the end of the time step after other terms have been accounted for: where where where k ZIN is the internal ratio of N to C and k ZIP is the internal ratio of P to C of the particular zooplankton class.

Configuration 1 -NOBAC
This configuration assumes organic matter is mineralised at a rate that is not dependent on the bacterial biomass (i.e. the bacterial biomass is assumed non-limiting and f B (B) in Eq. ( 1) is fixed at 1).This approach moves C, N and P fluxes between DOM and DIM proportionally.Since there are no bacteria simulated for micrograzers to graze upon, the grazing preferences were adjusted to consume POM in place of bacteria, thereby assuming bacterial biomass is lumped within the detrital pool.The grazing rate of microzooplankton simplifies to where POC is used to determine the grazing rate and PON and POP are consumed at a rate commensurate with their local stoichiometry at the time of grazing.The grazing rate parameter (g MAX ) was adjusted to make G Z3 (POC) in NOBAC approximately equal to G Z3 (B) in BAC + DIM (Table 4), to keep the general C flow and biomass patterns comparable between these simulations.

Configuration 2 -BAC − DIM
This configuration includes the heterotrophic bacteria state variable, B, however, they are restricted to DOM uptake during the mineralisation process.Under this scenario, the bacterial biomass and their mineralisation rate increase and decrease depending on temperature and organic matter availability, but their nutrient requirement must be satisfied from the DOM pool.The basic equations for BAC − DIM are similar to NOBAC, except the inclusion of the bacterial equation and their associated growth and loss processes (Table 3).Bacterial uptake of DOC is similarly defined using Eq. ( 2) with f B (B) defined as Bacterial uptake of DON and DOP is based on the C mineralisation rate, converted according to the stoichiometric requirement of N and P (k BIN and k BIP ), but limited to the available pool to enforce mass conservation: Note that if they cannot support the stoichiometric requirement in line with the U DOC from the DON and DOP pool, then they take what is available and U DOC will be reduced accordingly.In this configuration, POM decomposition is also dependent on the changing bacterial biomass through f B (B) and micrograzers graze on bacteria (B) rather than POM.Therefore G Z3 (B) is set as

Configuration 3 -BAC + DIM
This configuration is an extension of BAC − DIM where bacteria compete with phytoplankton by supplementing their internal nutrient requirements through the uptake of inorganic nutrients when there is insufficient N and P in the DOM pool to support growth.The bacterial uptake of N and P requires the following additional terms (Table 3): If there is insufficient organic and inorganic N or P to support the carbon uptake rate, U DOC , the growth is limited to enforce mass balance as in configuration 2.

Structural sensitivity
The averages of a number of variables from the upper 10 m of the water column were computed over the simulated period (1997)(1998)(1999)(2000)(2001)) to be consistent with Gal et al. (2009).
The physical (T , DO), chemical (TN, TP, NO 3 , NH 4 , PO 4 ) and biological variables (A 1-5 , Z 1-3 ) of the NOBAC, BAC − DIM, and BAC + DIM were statistically compared by one-way ANOVA (5 % significance level, SPSS software version 18.0) and Multiple Comparisons (POST HOC, SPSS software version 18.0) to determine significant differences between the outputs of the alternative microbial loop submodels.Since model time series are not suited to ANOVA, our approach was to conduct monthly averages of surface layer model output to match the frequency of observational data that was used for validation.Given that the timescale for many processes is the order of days to weeks (e.g.Recknagel et al., 2013), this was done to reduce the degree of temporal auto-correlation between consecutive model points.

Parameter sensitivity
In addition, a sensitivity analysis of the impact of the microbial loop parameters on the simulated C, N and P cycles for BAC + DIM was conducted, since this configuration was considered to be the most similar to the actual dynamics of the lake.The limited selection of parameters were chosen based on the detailed global analysis of the complete set of ecological parameters by Makler-Pick et al. (2011a), and relevance to the microbial loop processes investigated here.A simple "one-at-a-time" sensitivity analysis was undertaken by scaling the parameters individually by +20 % and −20 %, similar to Bruce et al. (2006), and the degree of sensitivity of state variable concentrations and major process pathways for C, N and P cycles were compared.

Quantification of pools, fluxes and limitation
To determine the influence of the microbial loop on the food web, the numerous pools and fluxes of C, N and P were Table 3. Equations for C, N and P within nutrients, organic matter, bacteria and zooplankton pools.Note that the pools and processes related to phytoplankton are not included here for brevity since they are not different between the three configurations.For the complete balance equations including the dynamics of phytoplankton readers are referred to Gal et al. (2009).
is particulate organic matter sedimentation, S B is bacterial sedimentation), G Z3 is grazing by microzooplankton, Mz is zooplankton mortality and messy feeding, Ma is mortality of phytoplankton, R B is bacterial respiration, R Z3 is respiration of microzooplankton, P Z1 is predation by Z 1 , E A is phytoplankton excretion of DOM, E PO4 and E NH4 refer to bacterial mineralisation of nutrients, E DOM is DOM excretion from zooplankton, DSF is dissolved sediment flux, NIT is nitrification, DEN is denitrification, U DOM is dissolved organic matter uptake, either independent or linked to B biomass in the case of NOBAC and the other simulations, respectively.U NH4 , U NO3 and U PO4 refer to inorganic nutrient uptake, and the functions are designed to account for phytoplankton uptake only in the case of NOBAC and BAC − DIM, U (A), and phytoplankton and bacteria in the case of BAC + DIM, U (A, B).  2009) values adopted.
averaged over the simulation period, with both nutrient and biological state variables and fluxes being vertically integrated to provide lake-wide averages.
For each of the phytoplankton groups, the nutrient limitation functions, f a (N) and f a (P), at a depth of 1 m below the water surface were assessed to explore the impact of the microbial loop on phytoplankton nutrient limitation.The functions were calculated by the model based on the internal nutrient concentrations (Li et al., 2013): which range from 0 (extreme limitation) to 1 (no limitation).

Comparison of model simulations
As expected, the simulated temperature and dissolved oxygen patterns were similar in the three models and matched the field data equally well (Appendix A).The simulated major nutrient results (TN, TP, NO 3 , NH 4 and PO 4 ) for the three configurations were noticeably different in the surface waters, but were similar in the hypolimnion where sediment fluxes dominate (Fig. 2a).Most noticeable was the reduced surface water concentrations of NH 4 and NO 3 and seasonal PO 4 spikes in the simulated output of BAC − DIM, relative to BAC + DIM.Also, relative to BAC + DIM, increased levels of TN were simulated in NOBAC and BAC − DIM.In terms of the microbes, all three configurations followed similar seasonal trends, however, noticeable differences were reduced bacteria and Peridinium and increased Aphanizomenon concentrations in BAC − DIM (Fig. 2b).The impact of the three alternative microbial loop configurations on the 15 physical, chemical and biological variables was statistically analysed by one-way ANOVA and Multiple Comparisons (Table 5).Although the simulated results for T and DO were not significantly different in the three configurations, the simulated results for nutrients were significantly different (p < 0.05): NH 4 , TN and TP of BAC + DIM were significantly different from BAC − DIM and NOBAC; NO 3 and PO 4 of BAC + DIM were significantly different from BAC − DIM, but similar to NOBAC.Biological variables were also significantly different between these microbial loop configurations: Peridinium, Aphanizomenon, and microzooplankton of BAC + DIM were significantly different from NOBAC and BAC − DIM; predatory zooplankton within BAC + DIM were significantly different from NOBAC but similar to BAC − DIM; macrograzers within BAC + DIM and BAC − DIM were significantly different from NOBAC but similar to each other; Microcystis of BAC + DIM was also significantly different from BAC − DIM.
Validation metrics for the three simulations are presented in Appendix A. Using these typical measures of performance, they show that the nutrient variables are significantly better predicted by BAC + DIM, however, the microbe predictions are generally comparable across the three configurations despite the significant differences reported above.

Model parameter sensitivity analysis
Several phytoplankton state variables, microzooplankton, and the various process pathways that connected them, were particularly sensitive to a number of key microbial loop parameters (above the 20 % sensitivity level) (Fig. 3).In particular, Peridinium and Microcystis were sensitive to the diameter of POM particles (d POM ) and the bacterial optimum temperature (T OPTB ).In addition to d POM and T OPTB , Microcystis was sensitive to the zooplankton internal N : C ratio (k ZIN ), and Aphanizomenon was also highly sensitive to T OPTB (> 50 %).Microzooplankton biomass, bacterial grazing rates and zooplankton excretion rates were strongly sensitive to K Ze (> 30 %), with mild sensitivity to d POM , T OPTB , K ZIN and the half-saturation constant for bacterial function (K B ).The DOM concentration was sensitive to d POM , particularly for N (> 50 %), and the maximum bacterial DOC uptake rate (µ DECDOC ), and K B and k ZIN (> 30 %).Looking specifically at the process pathways, rates of algal excretion and algal uptake were sensitive to T OPTB , particularly in the P cycle (> 30 %).To summarise, the model output was most sensitive to changes in the microbial loop parameters d POM , T OPTB and K Ze , which had a significant effect on DOM, the biomass of Peridinium, cyanobacteria, heterotrophic bacteria and microzooplankton.

Nutrient pools
The multi-annual and lake-wide nutrient pools were compared between the three microbial loop configurations to understand how the microbial loop shifts the partitioning of nutrients between different ecosystem compartments (Table 6).In each configuration, the stoichiometry of the POM, DOM and DIM pools was free to change, whereas the stoichiometry of zooplankton and bacteria were fixed, and the stoichiometry of phytoplankton was allowed to vary only within the range prescribed by the minimum and maximum parameters of internal nutrient ratios.In each configuration, the DIC pools were similar, but the DOC pool in BAC + DIM was significantly lower (1.79 mg C L −1 ) than in NOBAC (9.56 mg C L −1 ) and BAC − DIM (7.81 mg C L −1 ).Similarly the DON and DOP pools in BAC + DIM were also lower than the corresponding pools in NOBAC and BAC − DIM, even though bacteria were able to take up DIN and DIP to meet their nutrient needs in this configuration.The N:P ratio of DOM in NOBAC was 307 : 1, and with bacteria included (both BAC − DIM and BAC + DIM), the N : P ratios increased significantly to 28 475 : 1 and 3543 : 1 respectively.For configurations with dynamically simulated bacteria, the DIP pools in BAC − DIM (6.4 × 10 −3 mg P L −1 ) and BAC + DIM (5.2 × 10 −3 mg P L −1 ) were higher than that in NOBAC (3.6 × 10 −3 mg P L −1 ), suggesting enhanced P availability for phytoplankton uptake when bacteria were present.The POM pools in BAC − DIM and BAC + DIM were also higher than those in NOBAC.
The biomass of bacteria and zooplankton varied in the different microbial loop configurations, although the N : P stoichiometry of zooplankton and bacteria were fixed at 5 : 1 (bacteria), 27 : 1 (Z 1 , predatory zooplankton), 20 : 1 (Z 2 , macrograzers), and 28 : 1 (Z 3 , microzooplankton).When bacteria were able to uptake dissolved inorganic nutrients in BAC + DIM, the total bacterial biomass was 2.7 times larger than for BAC − DIM.For zooplankton, biomass of microzooplankton (Z 3 ) was similar in NOBAC and BAC + DIM and significantly lower in BAC − DIM.For predatory zooplankton (Z 1 ) simulated biomass was greatest in NOBAC and lowest in BAC + DIM and for macrograzers (Z 2 ) it was greatest in NOBAC and lowest in BAC − DIM.
The biomass and N : P stoichiometry of the five simulated phytoplankton groups each varied in response to the presence of bacteria in BAC − DIM and then with the addition of bacterial uptake of inorganic nutrients in BAC + DIM.For Microcystis, Peridinium, nanophytoplankton and Aulacoseira, the total C, N, and P content in BAC + DIM were higher than those in BAC − DIM.Similarly, the molar N : P ratios of phytoplankton in BAC + DIM (   8 : 1; nanophytoplankton 47 : 1; Aulacoseira 16 : 1) were also higher than their N : P ratios in BAC − DIM (Peridinium 59 : 1; Microcystis 4 : 1; nanophytoplankton 18 : 1; Aulacoseira 10 : 1).Conversely, for Aphanizomenon, simulated biomass in BAC − DIM was higher than in BAC + DIM, but no change was observed in their molar N : P ratios (4 : 1).Overall, the total phytoplankton biomass in BAC + DIM was higher than that in BAC − DIM despite this simulation, including competition for inorganic nutrients by bacteria.

Nutrient fluxes
Simulated fluxes of C, N and P from the three microbial loop configurations representing the dominant C, N and P recycling pathways demonstrate significant differences in the relative magnitude of bacterial mineralisation, zooplankton excretion, zooplankton grazing, and bacterial competition with phytoplankton for inorganic nutrients (Fig. 4).The simu-lated rate of algal primary productivity in the BAC + DIM and NOBAC configurations was higher than that simulated in BAC − DIM.Relative to the algal CO 2 fixation rate (defined as 100 % for each simulation), in NOBAC, bacterial respiration returned 32.7 % of the total DIC assimilated by phytoplankton, which was fuelled by DOC from microzooplankton excretion (29.7 %), hydrolysis of particulate detritus (26.7 %) and algal exudation (12.0 %); in BAC − DIM, bacterial respiration returned 43.3 %, fuelled less by DOC from microzooplankton excretion (10.8 %), and more from hydrolysis of particulate detritus (62.8 %) and a similar amount from algal exudation (10.9 %); in BAC + DIM, the magnitude of bacterial respiration was 77.3 % of the total algal fixed carbon, with similar proportions as in BAC − DIM; DOC from microzooplankton excretion (10.2 %), hydrolysis of particulate detritus (69.5 %) and algal exudation (10.7 %).
For nitrogen, the algal DIN uptake rates in BAC + DIM and NOBAC were similar, though greater than that in BAC − DIM.The algal DIN uptake percentage in BAC + DIM was 97.0 % relative to the 3.0 % uptake by bacteria (i.e.rates are normalised by the total inorganic N uptake rate).In NOBAC, static mineralisation recycled 77.4 % of the total DIN taken up by phytoplankton, with microzooplankton excretion being the primary source of organic N with a similar relative magnitude (68.4 %).In BAC − DIM, bacterial mineralisation recycled 47.2 % of N, however only 17.6 % was supplied through microzooplankton excretion due to the lower Z 3 biomass overall.In BAC + DIM, the bacterial mineralisation returned 74.3 %, with microzooplankton excretion supplying 21.5 %.As for carbon, in this simulation hydrolysis of particulate detritus was the major source of labile organic nitrogen (> 50 %) relative to that from zooplankton and phytoplankton excretion.
For phosphorus, the algal DIP uptake rate in BAC − DIM was higher than in BAC + DIM and NOBAC.In NOBAC, bacterial mineralisation replaced 84.2 % of total DIP assimilated by phytoplankton, and zooplankton excretion provided 29.3 % of this P to bacteria, with 41.5 % coming from POM hydrolysis and 10.5 % from algal exudation.In BAC − DIM, however, bacterial mineralisation recycled 94.0 %, with zooplankton excretion contributing just 12.8 % and the remainder coming from POM hydrolysis (67.1 %) and algal exudation (14.1 %).When uptake of DIM by bacteria was simulated (BAC + DIM), DIP uptake shifted significantly to 27.8 % by algae and 72.2 % by bacteria, when normalised relative to the total PO 4 uptake rate.Of this total consumed PO 4 , bacterial mineralisation was responsible for replacing 95.9 %, and DOM supplied by zooplankton excretion contributed 10.9 %, algal exudation contributed 3.4 % and POM hydrolysis 28.3 %.These fractions were less than that for C and N due to the large rate of supplementation of PO 4 .
Note that in all cases the amount of dissolved inorganic N and P that came from recycling compared to the inflows and sediment fluxes was very high.For example, the BAC + DIM model predicted that 95.9 % of dissolved inorganic P was sourced from recycling within the water column, only 4.4 % from the sediments, and less from the inflows.For N, the model predicted a reduced dependence on recycling (approximately 67 %), higher sediment flux (22.3 %) and a similar low contribution (0.7 %) from the inflows.

Patterns of phytoplankton biomass
In conjunction with variability in temperature, light and vertical mixing, changes in nutrient availability resulting from the dynamic nutrient recycling processes led to variation in phytoplankton nutrient uptake and their nutrient limitation functions, f a (N) and f a (P).The different patterns of seasonal variation in the nutrient limitation of the five simulated phytoplankton groups within the three model configurations highlight the potential for microbial loop sub-model structures to influence phytoplankton growth response (Fig. 5).For example, Peridinium in NOBAC and BAC + DIM was predicted to have periodic N and P co-limitation, however, in BAC − DIM, N limitation was predicted to dominate most of the year.For Aulacoseira, in BAC − DIM, N and P colimitation was experienced most of year, but in NOBAC and BAC + DIM, it displayed more P limitation.For Microcystis, in NOBAC, P was the limiting factor for algal growth, however, in BAC − DIM, it was predicted to switch from P limitation to significant N limitation, and in BAC + DIM it experienced significant P limitation with an annual occurrence of N and P co-limitation in spring.For the nanophytoplankton, in NOBAC and BAC + DIM, its growth was P limited with annual N and P co-limitation, but in BAC − DIM the growth switched between N limitation and P limitation annually.For Aphanizomenon, in all three configurations, P limitation dominated growth, since it is a N 2 -fixing species.

Model performance and suitability
Given the complexity of interactions affecting phytoplankton succession and bloom dynamics, our ability to predict complex microbial food webs accurately remains a challenge.To date, there are limited modelling examples for a complete lake ecosystem that confidently simulate the successional dynamics of phytoplankton and zooplankton at the level of multiple trophic complexity as described here.This is due to the nonlinearity of these complex models and a large number of uncertain processes and parameters with limited validation data (Arhonditsis and Brett, 2004;Rigosi et al., 2010;Mooij et al., 2010).Nonetheless, the simulations were successful in capturing the seasonal dynamics and some of the interannual variation of the key plankton functional groups in Lake Kinneret, though their absolute concentrations tended to be underpredicted.This is not unexpected given we have adopted a laterally averaged one-dimensional approach which is being compared to inherently patchy field data, known to be particularly relevant during Peridinium blooms (e.g.see Hillmer et al., 2008;Ng et al., 2011).However, the models were able to match the annual sequence and timing of the predicted peaks of these blooms, particularly the BAC + DIM configuration, which we consider to have the most biologically meaningful configuration.Within this simulation, the timescales of growth or decay of the biomass of biological variables generally matched the observed data, and seasonal trends were accurately captured for physical and chemical variables since the model responds significantly to the strong seasonal forcing of the lake (Makler-Pick et al., 2011a).While we acknowledge further improvements could be made, the focus of our study was not specifically to find the model with the best fit to the data but to use the dynamic model to help us gain insights into the significance of microbial loop processes on phytoplankton growth in accordance with the approach suggested by van Nes and Scheffer (2005) for application of ecological models to explore theory.For this purpose, the model captures the variability of key physical, chemical and biological processes to a suitable level to allow us to investigate the mechanisms governing the microbial interactions between the configurations.
Accordingly, different microbial loop configurations were found to have a significant impact on the sensitivity of most state variables based on the ANOVA and Multiple Comparison analysis.The predicted surface-water nutrient concentrations appeared to be the most sensitive variables to microbial configuration, with particular sensitivity noted in the concentrations of inorganic nutrients available for phytoplankton growth.Generally, it was noted that in BAC − DIM inorganic nutrients were lower on average even though the total nutrients were higher, due to a larger accumulation of organic matter over time, indicating stoichiometric constraints reduce the efficiency of mineralisation (discussed below).
The differences in predicted surface nutrient concentrations between the simulations led to the differences in predicted plankton biomass and growth rates.The structure of the microbial loop model had a significant impact on the total phytoplankton biomass, similar to the results of Faure et al. (2010) who explored the effect of microbial loop on DIN and phytoplankton for a coastal ecosystem.In this study, the analysis also includes phosphorus and several different functional groups of phytoplankton and zooplankton, and it was found that nutrients, Peridinium, Aphanizomenon and zooplankton were the main variables that showed sensitivity to assumptions related to microbial loop configuration.Indeed, the evolution of the three model structures reported here was the product of adding process complexity based on perceived deficiencies in the calibration and this is reflected with BAC + DIM having the closest representation to the real lake ecosystem (see Appendix A).
We note numerous recent developments related to modelling the impact of food quality on grazing rates, related to both prey stoichiometry (Mitra et al., 2007) and prey fatty acid (HUFA) content (Perhar et al., 2013).Given the specific focus of this study, our investigation has centred around the grazing rate of micrograzers who have grazed exclusively upon bacteria and have relatively stable biochemical composition compared to phytoplankton (Sterner and Elser, 2002).Our grazing rates also become limited as zooplankton become unable to meet their stoichiometric requirements from prey, however, we acknowledge incorporating food quality dynamics in the general grazing formulation of the zooplankton as an important area of further model development.

Role of the microbial loop in regulating nutrient flows
Specifically it was our aim to understand the mechanisms by which bacterial and microbial loop processes influence phytoplankton via changes in carbon and nutrient cycles: (i) bacterial mineralisation of organic nutrients, (ii) zooplankton excretion of labile material, and (iii) bacterial competition with phytoplankton for inorganic nutrients when organic matter quality is poor (i.e.nutrient deplete).By comparing fluxes between pools of C, N and P we were able to gain insights into the role of the microbial loop in the recycling of nutrients.
Bacterial mineralisation had a strong regulatory effect on nutrient recycling, and the BAC + DIM model predicted that 70 % of N and around 95 % of P available for phytoplankton growth was from bacterial mineralisation of organic mat-ter.These figures are based on a long-term simulation average and relative contributions were found to vary seasonally in response to temperature and organic matter availability.However, in terms of carbon biomass, the bacterial population was found to be relatively stable.A key result emerging from the simulations is that the lowest concentration of DOC occurred in BAC + DIM, suggesting bacterial metabolism is enhanced when nutrient supplementation is accounted for.Although reports of bacterial growth being C limited have been made in many lakes (Coveney et al., 1992), bacteria in our simulations were mainly limited by P and also occasionally co-limited by N and P, as indicated by the relative use of inorganic nutrients.In the model, the DOM is assumed to be relatively labile, however, in reality different bioavailability of the various organic matter constituents may become limited due to a lack of suitably bioavailable carbon.There is therefore scope for further extension of the model to understand how processes of mineralisation compare when multiple lability fractions of organic matter are considered.
It is well established that microzooplankton can transfer energy and nutrients via bacterial grazing to higher trophic levels due to their small size and high mass-specific grazing rates (Hart et al., 2000;Loladze et al., 2000), thereby playing an important role in carbon and nutrient recycling (Stone et al., 1993;Dolan, 1997;Hambright et al., 2007).In turn, larger zooplankton grazing on microzooplankton further provide organic matter for bacterial growth through excretion of nutrient rich organic compounds (DOM) and fecal pellet production (POM) (Peduzzi and Herndl, 1992).From this point of view, the recycling of organic nutrients is facilitated by bacterial consumers rather than bacteria themselves, known as consumer-driven nutrient recycling (CNR) (Elser and Urabe, 1999).The model simulations presented in this study have allowed us to estimate the significance of this pathway and characterise the relative contributions of upward and downward nutrient fluxes of these pathways.For example, in BAC + DIM, as a fraction of algal uptake, microzooplankton excretion was predicted to account for 10 % of C, 22 % of N and 11 % of P returned for mineralisation, which was significantly larger than that supplied from algal excretion for N and P (but not for C), and therefore different from the relative proportions consumed through bacterial grazing (18 % for C, 15 % for N and 20 % for P).This highlights the dissimilarity and decoupling of the C, N and P cycles and the importance of nutrient adjustments that occur during these microbial interactions, and is consistent with empirical work (Gaedke et al., 2002).The parameter sensitivity analysis indicated that both microzooplankton biomass and bacterial growth were particularly sensitive to the excretion fraction of the ingested material (K Ze ) grazed by microzooplankton.Adjusting this excretion fraction parameter not only impacted their own biomass and grazing rates, but also impacted the biomass of other zooplankton groups and the phytoplankton community more broadly, including Peridinium.Although Peridinium is not grazed directly by zooplankton, any reduction in nutrient supply from micrograzers leads to reduced P availability and ultimately reduced growth.These results reinforce that the interaction between phytoplankton and zooplankton is nonlinear and that there is both top-down (i.e.grazing-mediated) and bottom-up (i.e.microbial loop nutrient supply) processes shaping the phytoplankton community.Interestingly, the smaller microzooplankton have a significant overall impact shaping the food web structure in the model simulations despite having the lowest biomass.These findings are in line with the empirical studies of Hart et al. (2000) and Hambright et al. (2007), who highlighted the critical role of small micrograzers in the microbial loop processes.Since there exists a range of uncertainty surrounding the parameterisation of microzooplankton excretion, with large ranges being reported (Fasham et al., 1999;Faure et al., 2010), correct model parameterisation remains an important challenge for modellers.
In freshwater ecosystems, the concentration of DON can often be higher than that of DIN, and the DON pool plays an important role in providing N to both bacteria and algae (Berman, 1997(Berman, , 2001;;Berman and Bronk, 2003), though the latter is not considered in our model conceptualisation.In the present study, concentrations of DON were higher than those of DIN in NOBAC and BAC − DIM, which fits with observations by Berman and Bronk (2003).However, DON was lower than DIN in BAC + DIM where bacteria biomass and mineralisation rates were higher.As a result of increased DIN, DOP became the limiting factor when competition by bacteria for inorganic nutrients was included in the model configuration.Therefore, the variable stoichiometry of organic matter, and different stoichiometric requirements of various process pathways, leads to a complex interplay between the groups, and future studies should further consider the significance of organic matter stoichiometry, microzooplankton excretion rates and rates of nutrient immobilisation by bacteria when modelling planktonic food webs.We note the potential for N-fixation by heterotrophs, which is not considered in our model, but expect the results would be relatively insensitive to this since bacteria are mostly P limited (refer to the high N:P ratios for DOM in Table 6).

Impact of the microbial loop on the phytoplankton community
Bacterial competition for inorganic nutrients has a two-fold effect on phytoplankton growth by limiting nutrient supply and regulating the N : P ratio of available nutrients.The time series of nutrient limitation functions for the five simulated phytoplankton groups for each of the three alternative microbial loop configurations were used to decipher the effect of bacterial competition on phytoplankton growth.Whilst most freshwaters are considered to be P limited (Schindler et al., 2008), Elser et al. (2007) discusses that N and P colimitation is commonly also prevalent.During the simulation period in this study, Lake Kinneret had an average TN : TP ratio ∼ 50 : 1, indicating strong P limitation as suggested by other authors.However, Gophen (2011) argues N limitation also occurs, potentially due to large fractions of unavailable organic nitrogen distorting nutrient ratios (Ptanick et al., 2010).The model predictions of the five simulated phytoplankton groups were predominantly P limited, with potential for periodic N and P co-limitation depending on the microbial loop configuration.When organic matter became P depleted, it could not support bacterial growth and therefore PO 4 supplementation of bacteria was evident in the increased uptake rates (Fig. 4c).Bacteria generally have faster P uptake rates relative to phytoplankton (Berman, 1985), which in our model was captured by not limiting the rate of uptake of PO 4 by bacteria; as a result they were able to effectively out-compete the phytoplankton.The model results indicate that this stoichiometric regulation of bacteria through DIM supplementation shifted patterns of phytoplankton nutrient limitation.Several phytoplankton groups experienced differences in the degree of N and P limitation when bacteria were configured not to take up inorganic nutrients, as opposed to when bacteria were also consuming inorganic nutrients.The effect of competition on Peridinium manifested as N limitation when bacteria could not compete, but with competition, periods of phosphorus limitation emerged following periods of accelerated growth.For Aulacoseira, the lack of bacterial competition for nutrients led to severe N limitation relative to predominant P limitation, coinciding with the period of the Peridinium bloom.Similarly, for Microcystis and the mixed nanophytoplankton community, stronger N limitation simulated in BAC − DIM switched to predominant P limitation in BAC + DIM.It follows that bacteria-induced shifts in nutrient limitation can ultimately influence the overall biomass and composition of the phytoplankton community (Andersen et al., 2004).Indeed, here we noted relative differences in the simulated phytoplankton biomass, and in particular, when competition with bacteria for inorganic nutrients was simulated, Peridinium dominated and Aulacoseira also occurred in significant numbers.When this competition was not included, reduced Peridinium and Aulacoseira biomass were evident, with a corresponding significant increase in Aphanizomenon.Microcystis were also slightly reduced and the nanophytoplankton appeared to exhibit greater seasonality.Interestingly, focusing on the total phytoplankton biomass, the increased competition for nutrients by bacteria somewhat paradoxically led to higher total phytoplankton biomass overall.Therefore the effect of the competition has a role in shaping the community structure and timing of blooms (see also Li et al., 2013), but overall the micrograzer driven recycling of N and P positively promotes phytoplankton productivity.
The results highlight the importance of understanding the role that the microbial loop plays in nutrient recycling as a critical model component that must be considered when simulating phytoplankton dynamics in freshwater ecosystems.the role of the microbial loop in nutrient recycling (Mooij et al., 2010), with many studies based on extensions of the "N-P-Z-D" approach.In the present study, we identified that the phosphorus content of organic matter is a critical factor driving microbial loop processes, yet this is rarely parameterised in detail within lake ecosystem models, which generally maintain structures equivalent to our NOBAC simulation.We conclude therefore that not only should microbial loop processes and stoichiometric constraints between groups be further considered in future model studies but also that the parameterisation of these processes be supported with targeted empirical studies able to resolve the complexities of these interactions.

Appendix A: Model validation and suitability assessment
The performance of DYRESM-CAEDYM simulations has been reported in Gal et al. (2009) and Li et al. (2013) for the equivalent of our BAC + DIM simulation.The parameters for this simulation were derived from a manual calibration over the period 1997-1998, and error metrics were then compared over the full period 1997-2001.Parameter choice was based on a detailed synthesis of the large amount of laboratory and field process investigations conducted by the Kinneret Limnological Laboratory covering all aspects of model function (Table 2 in Gal et al., 2009).Fine-tuning of model parameters was conducted within acceptable ranges for numerous parameters, however no formal calibration and error minimisation routine was used.The performance of BAC + DIM and the NOBAC and BAC − DIM simulations is summarised in Table A1 for 16 observed state variables.Data were available from weekly sampling at several locations and depths and used to generate epilimnion and hypolimnion monthly averages.Note that only surface water data are compared here as no major differences were noted between the three simulations in the lake hypolimnion (except for NO 3 as seen in Fig. 2b).
Additional Spearman rank correlation assessments were also reported in the previous publications to indicate the performance of the model in capturing the general timing and interannual variation in phytoplankton biomass.Model performance assessments by Li et al. (2013) also included comparison of indirect metrics, including suitability for simulating water column and phytoplankton N : P ratios.
Macrograzers: Cladocerans) internal P concentration mg P L −1 ZIP 3 Zooplankton #3 (Micrograzers: Rotifers/Ciliates) internal P concentration mg P L −1 B BAC Heterotrophic bacterial C biomass concentration mg C L −1 BIN Heterotrophic bacterial internal nitrogen concentration mg N L −1 BIP Heterotrophic bacterial internal phosphorus concentration mg P L −1

Figure 2 .
Figure 2. Comparison of model simulations for (a) nutrient variables in the surface 10 m (left) and the bottom 10 m (right) of the water column (mg L −1 ), and (b) for the nine simulated biotic groups (mg C L −1 for A 1-5 and B, and mg C m −2 for Z 1-3 ).Data represents the monthly mean of samples collected over these depths.

Figure 3 .
Figure 3. Local sensitivity analysis of simulated state variables and process rates for the C, N and P cycles presented as the lake average absolute change after a ±20 % parameter shift.

Figure 4 .
Figure 4. Summary of simulated annual C (black), N (red) and P (blue) pathways for the three configurations: (a) NOBAC, (b) BAC − DIM, and (c) BAC + DIM.Note the dotted, dashed, and dash-dot lines emphasise configuration specific pathways.Selected fluxes relevant to the analysis are displayed as the lake-wide average % relative to the total DIM taken up by phytoplankton and bacteria (where relevant), with the flux rate in brackets (× 10 −5 mg L −1 d −1 ).

Figure 5 .
Figure 5.Comparison of nutrient limitation functions f a (N) and f a (P) respectively for the five simulated phytoplankton groups in a) NOBAC, b) BAC − DIM and c) BAC + DIM.
Few lake biogeochemical model studies, directly simulate www

Table 2 .
Summary of three microbial loop simulations configured.
Used in model studies where bacteria are simulated but stoichiometry is not specifically a constraint on bacterial production Most likely the closest representation to reality with bacteria biomass variable and inorganic nutrient uptake used to support bacterial growth requirement

Table 4 .
Gal et al. (2009)lated parameters used in the three model simulations (refer toGal et al. (2009)for other parameter values).

Table 5 .
Statistical analysis of water quality variables in the three microbial loop configurations by ANOVA and Multiple Comparisons.

Table 6 .
Summary of average values (1997Summary of average values ( -2001) )for C, N and P contents (mg L −1 ) and N : P molar ratios of the various food web components in different microbial loop configurations.

Table A1 .
Comparison of R 2 and RMSE model validation metrics for the surface 10 m of Lake Kinneret for the main simulated variables.The best performing model configuration is highlighted bold.DIM) (BAC − DIM) (NOBAC) (BAC + DIM) (BAC − DIM) (NOBAC)