Plankton functional type modelling: running before we can walk?

Abstract

Biogeochemical cycling in marine systems is intimately linked to the activity of specific plankton functional types (PFTs) such as diatoms, coccolithophores and nitrogen fixers, thereby providing a focus for contemporary modelling studies. Incorporating extra complexity beyond simple nutrient-phytoplankton-zooplankton-detritus (NPZD) models is, however, fraught with difficulties: poorly understood ecology; lack of data; aggregating diversity within functional groups into meaningful state variables and constants; sensitivity of output to the parameterizations in question and their physical and chemical environment. Although regional models addressing the seasonal succession of plankton types have achieved some degree of success, predicted distributions of PFTs in global biogeochemical models have thus far been less than convincing. While the continued articulation of detail in ecosystem models is surely the way forward, I argue that this can only be so with due care and attention to the formulations employed and a healthy dose of scepticism regarding model outcomes. Future directions should emphasize building up complexity gradually, objective assessment of the resulting parameterizations, and variety in approach such as the use of empirical alternatives to the fully dynamic representation of PFTs in models.

Received July 12, 2005; accepted in principle September 8, 2005; accepted for publication September 28, 2005; published online October 5, 2005
 Communicating editor: K.J. Flynn

INTRODUCTION

Plankton play an integral role in marine biogeochemical cycles, and in particular the so-called ‘biological pump’, whereby photosynthesis by phytoplankton in surface waters, followed by sinking and transport of organic matter to depth, enriches the deep ocean in carbon. Modelling provides an important tool for understanding the dynamics of this pump and its possible response to changing climate. Early models of marine ecosystems were necessarily simple due to the infancy of the discipline and the technological constraints of available computing facilities, generally conforming to the so-called NPZ or NPZD configurations in which nutrients, phytoplankton, zooplankton and detritus form the main model structure (e.g. Steele, 1974; Wroblewski et al., 1988). NPZD-type models have subsequently formed the basis of many 3D general circulation model (GCM) biogeochemical modelling studies (e.g. Sarmiento et al., 1993; Six and Maier-Reimer, 1996; Oschlies and Garcon, 1998; Palmer and Totterdell, 2001). But, when modelling biogeochemical cycles and the flows of energy and material through marine ecosystems, can reliable predictions be expected when aggregating the myriad of plankton species into P and Z? Phytoplankton, for example, can be divided into various plankton functional types (PFTs) including diatoms, coccolithophorids, nitrogen fixers, picophytoplankton, phytoflagellates and dinoflagellates (Totterdell et al., 1993). The current trend in marine ecosystem modelling, and, in particular, models examining the current biogeochemical state of the ocean and its potential future response to anthropogenic perturbations, is to incorporate PFTs, the argument being that realistic predictions can only be expected if key processes associated with system feedbacks are represented (Doney, 1999; Pomeroy, 2001; Dearman et al., 2003). Complexity in nature is mirrored by complexity in models. The advent of modern computers has made possible adoption of the PFT approach even in biogeochemical modelling studies at the global scale (Moore et al., 2002; Gregg et al., 2003; Le Quéré et al., in press).

Simply building up model complexity does not necessarily guarantee improved predictions unless parameterization is sufficiently accurate and robust. Do we understand the ecology of PFTs well enough to do this? And, the proof being in the pudding, so to speak, are models in which PFTs are represented actually performing better than the NPZD models that they are designed to supersede? In this article, I ask whether the gathering momentum towards using ever more complex PFT models, particularly in GCMs, is metaphorically in danger of making us try to run before we can walk. Possible directions for future research are suggested with a view to promoting the representation of PFTs in biogeochemical models.

DEALING WITH COMPLEXITY: STATE OF THE ART

Are NPZD models too simple?

NPZD-type models continue to demonstrate their worth in both regional and global modelling studies (e.g. Kawamiya et al., 2000; Friedrichs and Hofmann, 2001; Palmer and Totterdell, 2001; Spitz et al., 2001; Aumont et al., 2002; Anderson and Pondaven, 2003; Hood et al., 2003; Schartau and Oschlies, 2003). Although predictions are by no means perfect, these models are by and large successful at capturing bulk system properties such as chlorophyll and primary production which are constrained by nutrient availability, light and grazing. For example, simulations of the timing and magnitude of phytoplankton blooms in seasonal regimes such as the Arabian Sea (e.g. Olascoaga et al., 2005), as well as year-round low chlorophyll in high-nutrient low-chlorophyll areas such as the Subarctic Pacific (e.g. Fasham, 1995), generally show good agreement with data. Errors are likely caused as much by deficiencies in model physics as by problems in ecosystem model parameterization.

The precise nature of the food web is, however, important as regards many aspects of biogeochemical cycling, export flux being a typical example. Accurate representation of export in models might require small and large phytoplankton to be distinguished, the former leading primarily to recycling via the microbial loop (Azam et al., 1983), the latter to export via sinking particles. Diatoms, in particular, are commonly believed to play a major role in export (Sancetta et al., 1991; Kemp et al., 2000). Climate feedbacks may be linked to the ecology of individual plankton groups. Enhanced stratification in the ocean may, for example, favour the growth of nitrogen fixers that are most prevalent in the stratified waters at low latitudes (Karl et al., 1997). Or changes in seawater carbonate chemistry as a result of anthropogenic CO₂ invading the ocean may decrease calcification in major planktonic groups such as the coccolithophores and foraminiferans (Riebesell et al., 2000; Barker and Elderfield, 2002). Calcifiers in turn impact on alkalinity and carbonate chemistry, as well as contributing to some of the densest ballasts observed in sinking particles (Klaas and Archer, 2002). The case for superseding NPZD-type models with those that include PFT is clear: biogeochemical cycling in marine systems is intimately linked to particular plankton groups if not in some instances primarily to individual plankton genera or species (e.g. Trichodesmium, Emiliania, Phaeocystis).

How good are existing PFT models?

The seasonal succession of plankton species provided much of the initial focus for PFT modelling studies. An early example is the FYFY model developed by Van den Berg et al. (Van den Berg et al. 1996) which included six PFTs—diatoms, N specialists and P specialists, each in grazed and non-grazed forms—and was used to study the temporal and spatial variability of phytoplankton biomass and succession in the southern North Sea. FYFY was superseded by, the European Regional Seas Ecosystem Model (ERSEM) (Baretta-Bekker et al., 1997), the latest version of which includes picoflagellates, flagellates, dinoflagellates and diatoms as phytoplankton functional groups. Initial predictions for PFTs using this model met with only moderate success (Ebenhöh et al., 1997; Pätsch and Radach, 1997). Recent simulations have, however, provided an improved representation of the seasonal progression from diatoms to flagellates in the North Sea (Allen et al., 2004; Archer et al., 2004), the emphasis put on ensuring a good physical representation of the system likely playing a role in this success. Another model targeted at understanding seasonal succession is the SeaWater Microbial Community model (SWAMCO) which was used successfully to study bloom development at the Polar Front (Lancelot et al., 2000) and blooms of diatoms and Phaeocystis in the North Sea (Lancelot et al., 2005). In similar fashion, Merico et al. (Merico et al., 2004) developed a model of seasonal succession on the Bering Sea shelf that included diatoms, flagellates, dinoflagellates and the coccolithophore Emiliania huxleyi, the timing and duration of E. huxleyi blooms agreeing with remotely sensed data.

Diatoms appear to be reasonably well simulated in most of the models listed above, perhaps, a consequence of relatively straightforward parameterization—they have high growth rates seemingly allowing them to outcompete other phytoplankton when dissolved silicate is readily available. Other PFTs are, however, less easy to constrain. Regional modelling studies have in many instances achieved a reasonable match with data, but it remains to be shown that PFT parameterizations are robust, i.e. that they will perform well when subjected to new scenarios or altered forcing, without adjustment of parameter values. Consider, for example, separate models examining the distribution of nitrogen fixers in the subtropical Atlantic and Caribbean by Hood et al. (Hood et al., 2004) and Lenes et al. (Lenes et al., 2005), respectively. Both models broadly reproduced spatial distributions of Trichodesmium but used entirely different parameterizations, Hood et al. (Hood et al., 2004) focusing solely on the effect of light on growth, Lenes et al. (Lenes et al. 2005) including also the effects of limitation by phosphorus and iron. It would be interesting to see how well each parameterization would perform if supplanted into the other’s geographical domain.

Several recent biogeochemical modelling studies have involved implementation of the PFT approach in global GCMs. The first such model was developed by Moore et al. (Moore et al., 2002, 2004), including small phytoplankton, diatoms and nitrogen fixers as functional types. Predicted patterns of surface nutrients, primary production and export were generally good. Despite a somewhat limited comparison with data, the authors claim that predictions for PFTs captured ‘the known ecological contributions of key phytoplankton functional groups’ (Moore et al., 2004) and, moreover, that results indicated that regional floristic shifts as a result of climate change may be as important as alterations in bulk integrated global productivity (Boyd and Doney, 2002). A second global model incorporating PFTs was published by Gregg et al. (Gregg et al. 2003) with diatoms, chlorophytes, cyanobacteria and coccolithophores as phytoplankton groups. Predicted PFT distributions, other than coccolithophores (see below) conformed with general expectations: diatoms prevailing at high latitudes, coastal and equatorial upwelling regions, cyanobacteria in the central ocean gyres, and chlorophytes inhabiting the transitional regions. Finally, a third global PFT model, the Dynamic Green Ocean Model (DGOM) is being developed by Le Quéré et al. (Le Quéré et al. in press). They identified six key phytoplankton functional types, namely picoautotrophs, nitrogen fixers, DMS producers, calcifiers, diatoms and mixed phytoplankton and modified the PISCES biogeochemistry model (Bopp et al., 2003) to run a prototype version of DGOM that included the last three PFTs in this list.

Comparing the models in terms of predictions for calcifiers, the DGOM model predicts blooms between 40°N and 40°S, at odds with observations showing coccolithophore blooms in the northern North Atlantic and subarctic Pacific (Holligan et al., 1993; Brown and Yoder, 1994). In contrast, the Gregg model correctly predicts blooms in the North Atlantic but also predicts their occurrence south of the Polar Front in the Southern Ocean during the austral summer. Neither evidence from cruises (Hasle, 1960; Findlay and Giraudeau, 2000), nor from the sediment record (McIntyre and Bé, 1967), supports the presence of coccolithophore blooms in this area. Moreover, current wisdom indicates that whereas blooms of the coccolithophore E. huxleyi occur at high latitudes, much of global calcification is, in fact, carried out by larger low-latitude species (Baumann et al., 2004). It appears that neither model captures this phenomenon. Regarding parameterization, both models assigned low maximum growth rates and high affinity for inorganic nutrients to calcifiers, but neither can address the potentially important impact of N:P ratio (see below ‘Poorly understood ecology’). Parameterizations seem at best speculative, predictions as yet unreliable, in these albeit preliminary attempts to incorporate PFTs into global biogeochemical GCMs.

SCALE OF THE PROBLEM

Poorly understood ecology

Confidence in model parameter values is achieved by demonstrating understanding of the ecology of the organisms in question. Focusing again on the calcifiers, there is a tendency to think of this group as being dominated by the coccolithophores, particularly E. huxleyi. A number of hypotheses have been used to explain the conditions which favour the growth of this species including high light, high NO₃:PO₄ ratio, low dissolved CO₂ and high carbonate saturation state (Tyrrell and Merico, 2004). Recent work by Lessard et al. (Lessard et al. 2005) has, however, indicated that E. huxleyi may also bloom when NO₃:PO₄ ratios are low and that attention should be directed towards examining the role of organic substrates as a source of nutrients. Even if the factors limiting the growth of E. huxleyi could be established with confidence, top–down processes also have a bearing on population numbers, selective predation by microzooplankton (Verity, 1991) being of particular relevance. High dimethylsulphoniopropionate-lyase activity by E. huxleyi may act as a chemical deterrent to grazers (Strom et al., 2003) with losses of this species being relatively low compared to other phytoplankton (Fileman et al., 2002). Viruses may also be an important contribution to the demise of E. huxleyi blooms (Wilson et al., 2002). One can equally ask how well we know the ecology of other phytoplankton types such as nitrogen fixers and picophytoplankton.

Aggregation

The aggregation of diversity in plant and animal life into a limited number of state variables poses a perennial problem for ecosystem modellers, the representation of PFTs being no exception. Once again, calcifiers are a case in point. E. huxleyi is by no means the only coccolithophore species in the ocean, with other species such as Florisphaera profunda, Umbellosphaera irregularis and Umbellosphaera tenuis contributing significantly to algal abundance at the Hawaiian Ocean Time series and Bermuda Atlantic Time Series sites (Cortes et al., 2001; Haidar and Thierstein, 2001). These different species have quite different niches to E. huxleyi, F. profunda, for example, contributing to the shade flora of the thermocline. Moreover, organisms other than the coccolithophores also produce calcium carbonate, notably the planktic foraminefera that may contribute 23–56% of the global CaCO₃ flux at 100 m (Schiebel, 2002). Being heterotrophic, these organisms have an entirely different, and even less well known, ecology to the coccolithophores.

Similar arguments can be made for other phytoplankton groups. Synechococcus and Prochlorococcus, for example, both major contributors to the picoplankton, show strong niche segregation with respect to light and other factors (Agusti, 2004). Moreover, whereas I have focused primarily on representation of phytoplankton types in models, equally tough decisions need to be made regarding zooplankton and even bacteria. It may seem reasonable to divide zooplankton between micro- and mesozooplankton functional types (e.g. Le Quéré et al. in press), but the reality is that each of these groups is highly diverse. The ciliates alone exhibit species-specific traits such that they cannot be considered as a functional group any more than can all phytoplankton or mesozooplankton (e.g. Montagnes, 1996). Large herbivorous grazers include copepods, salps and euphausiids, often with complex life histories that differ markedly between the dominant species of different ecosystems (e.g. Parsons and Lalli, 1988). Furthermore, the distinction between phytoplankton and zooplankton is in reality blurred, mixotrophy being widespread in the marine environment (Sanders et al., 2000). The argument can be extended to bacteria. Cottrell and Kirchman (Cottrell and Kirchman, 2000) conducted experiments indicating that degradation of dissolved organic matter in the ocean depends on a diverse assemblage of bacteria, recommending that carbon cycling by bacteria should be most accurately described by using three groups instead of the single compartment currently used in biogeochemical models. Deriving distinct functional types with unique parameterizations poses a serious problem for modellers.

All in the interactions

Modelling embraces the reductionist philosophy that the behaviour at a chosen level in a system can be predicted from rules governing the behaviour of elements at lower hierarchical levels. However, the key to understanding system behaviour lies not only in describing the subunits, but depends crucially on the interconnections between them. This is particularly so in complex systems which are characterized by many interacting parts, consistent with the Latin ‘complexus’ which signifies being entwined or twisted together. The outcome of these interactions is emergent behaviour giving rise to global dynamics that outlast any of the component parts (Bar-Yam, 1997).

Understanding the nature of emergence is a key component of modelling complex systems. Surprising and unforeseen artefacts are always possible when elementary subprocesses are coupled together to form a larger system, such that subtle and poorly understood interactions may restrict us from analysing system behaviour using the reductionist approach (Casti, 1994). Various biological modelling studies have indeed shown remarkable sensitivity to the exact form of the equations used (Wood and Thomas, 1999; Gross et al., 2004; Fussmann and Blasius, 2005). In similar fashion, model sensitivity to PFT parameterizations is likely when functional groups are represented as a series of interacting differential equations in models. The implication is clear: accuracy is required in the PFT parameterizations and, moreover, in the representation of the physico-chemical environment. The latter is worthy of note given the coarse resolution of many physical models and the knowledge that much of the action occurs at the mesoscale or smaller (Levy et al., 2001). The magnitude of these problems should not be underestimated.

Salvation in tuning parameters?

The standard modelling trick if one is left hunting for parameter values is to tune them in order that model output shows agreement with data. The vexatious problem of underdetermination then rears its ugly head deriving from the idea that, for any given set of observations, it is always possible to construct many different and incompatible theories that fit the data equally well (Quine, 1975). Complex models with their many degrees of freedom are in principle most easily fitted to any particular data set, but the number of parameters that must be fitted quickly surpasses our ability to constrain them properly from observations (Denman, 2003). This point is forcefully made by Matear (Matear, 1995) who optimized the parameters of three ecosystem models to fit nitrate, phytoplankton, mesozooplankton and primary production data at Station P in the subarctic Pacific. Error analysis indicated that the data were only sufficient to resolve up to ten independent model parameters.

The current problem is, however, very often not that acceptable parameter values can be conjured at will by tuning methods but rather to identify any single parameter set or model that performs well for a chosen scenario. Finding, for example, any parameterization of coccolithophores that permits their global distribution to be realistically simulated would represent some sort of progress. Further, undertaking such exercises might be useful in stimulating ideas for further research. A certain amount of exploring parameter space may therefore be justified, even when using global GCMs. Nevertheless, one always has to ask whether the resulting parameter values are to be trusted, the charge often being raised that modellers can produce any outcome that they so desire (Aber, 1997). Validation, as always, remains a key issue for modellers, who are frequently criticized for putting most effort into calibrating models and not testing genuine validity by comparing output with independent data (e.g. Arhonditsis and Brett, 2004). It is at least necessary to demonstrate robustness through parameterizations that perform well under a wide range of conditions rather than having parameters merely tuned for individual locations.

POINTERS FOR PROGRESS

The road ahead is tortuous for the aspiring ecosystem modeller; the development of reliable and robust parameterizations of PFTs is a challenge indeed. Unfortunately, there are no simple solutions to the problems of how to implement complex PFT models in GCMs when addressing issues such as global biogeochemical cycles and their possible response to climate change. A few pointers for progress are therefore to indicate possible future directions for those pursuing the PFT modelling approach:

• The development of new and improved ecosystem model formulations remains a priority. Attention to detail matters. Overly simplistic formulations should not be used if demonstrably superior alternatives are available, an example being the erroneous use of Monod models when simulating multi-nutrient limitation of phytoplankton (Flynn, 2003). It would seem logical that complexity is built up gradually, finding out what is important and what is not. New formulations for particular organisms may, for example, be tested as stand-alone subroutines before inclusion in site-specific, regional and finally global models.

• As compared to chlorophyll, availability of PFT data for model validation is very much at a premium. Nevertheless, data sets do exist (e.g. Widdicombe et al., 2002) and, in particular, there is considerable information on plankton seasonal succession that is surely of relevance for PFT modellers. Remotely sensed chlorophyll data have proved useful for validating global models but extending to PFTs poses a significant challenge. One promising approach in this respect is the development of bio-optical methods, various algorithms having been developed for mapping the distributions of coccolithophore blooms (Brown and Yoder, 1994; Brown and Podestá, 1997), Trichodesmium (Subramaniam et al., 2002) and diatoms (Sathyendranath et al., 2004). Empirical relationships relating measurements of in situ pigment to those of remotely sensed ocean colour were worked out by Alvain et al. (Alvain et al., 2005) who proceeded to generate global maps of haptophytes, Prochlorococcus, Synechococcus-like cyanobacteria and diatoms. Further work is required to assess the reliability of these bio-optical approaches.

• Empirical representation of PFTs in models provides a possible alternative to a fully dynamic set of interacting differential equations. One such case is the common method of taking into account the biogeochemical role of calcifiers by assuming a fixed calcification ratio relative to primary production. For example, a base calcification rate equal to 5% of photosynthesis by small phytoplankton was set by Moore et al. (Moore et al. 2002) who concluded that ‘at present it is not possible to model dynamically or predict calcite formation in the ocean’. Further extensions of this type of approach might be considered. The systematization of phytoplankton ‘life-forms’ into a mandala based on intensity of physical forcing, pioneered by Ramon Margalef, may, for example, provide some potential for empirical representation of PFTs in models. Focusing on seasonal succession, Margalef ordinated diatoms and dinoflagellates on the basis of nutrient supply and intensity of turbulence, the former thriving in turbulent, nutrient-rich waters and the latter preferring stratified oligotrophic conditions (Margalef, 1978; Margalef et al., 1979). Coccolithophores may occupy an intermediate position in the mandala (Estrada and Berdalet, 1997). More complex mandalas, incorporating various plankton groups and environmental factors, have followed on from Margalef’s original inception of the idea (Smayda and Reynolds, 2001; Balch, 2004). This ordination of PFTs on the basis of environment is similar in nature to terrestrial biogeography models in which plant functional types such as ‘tropical evergreen’ or ‘cool-temperate conifer’ are predicted as a function of climate and soils (e.g. Prentice et al., 1992).

• Improved methods are required for quantifying the uncertainty associated with model outcomes (uncertainty analysis) and the main factors—model structure, parameter values, data and forcing functions—that contribute to that uncertainty (sensitivity analysis). Model intercomparison is one way of assessing the suitability of different ecosystem model structures and formulations, requiring optimization of parameters and objective evaluation of performance against independent data (Friedrichs et al., in press). Estimation of uncertainty associated with model structure, as opposed simply to parameter values, is a challenging task. Advanced methods are slowly coming to the fore, such as Bayesian hierarchical modelling (Clark, 2005) which accommodates complexity by dissecting a problem into levels, leading to quantitative measures of uncertainty associated with both model structure and parameterization. The development of statistical approximations of large models known as statistical emulators (e.g. Logemann et al., 2004) in order to make multiple simulations computationally feasible is a further development of significance.

CONCLUDING REMARKS

The case for developing PFT models is clear: potentially important system feedbacks are linked to particular plankton groups. The immediate challenge is to demonstrate that these models are performing better than the NPZD models that they are designed to supersede. Progress is certainly being made with regional models looking at seasonal succession of plankton types, although demonstration of the robustness of the parameterizations involved remains a priority if such models are to be used in predictive mode or applied to wider geographical domains. Validation of these parameterizations must focus on the PFTs themselves and not just the successful reproduction of bulk properties such as chlorophyll and nutrients, which is often the case.

There tends to be a tradeoff between complexity and generality in models (Levins, 1966). It, therefore, seems surprising that the most complex models, with multiple PFTs, are often being incorporated directly into global GCMs without thorough testing at the regional scale first. It is here that I believe that we are indeed in danger of trying to run before we can walk. Poorly understood ecology, difficulties of aggregating diverse organisms into unique state variables, lack of PFT data along with the need for objective model data comparison methodologies and the unpredictable nature of the emergent dynamics of interacting equations are problems not to be taken lightly. Whereas bulk properties such as total chlorophyll are constrained by nutrient availability, light and grazing, the relative balance of phytoplankton types within that total, as well as different groups of zooplankton and bacteria, may be sensitive to the parameterizations in question. The formulations and parameter values assigned to PFTs are often less than convincing. When referring to the current generation of prognostic models being used to determine the response of ocean biology to warming, Sarmiento et al. (Sarmiento et al., 2004) comment that ‘the ecosystem models on which such studies are based are immature and that much work needs to be done before their results can be trusted’. Of course, this comment could be applied to simple and complex ecosystem models alike but nevertheless emphasizes the scale of the difficulties currently facing the modelling community. A further aspect is that complexity in biology needs to be matched with an appropriate complexity in the representation of the physical and chemical environment. Many GCMs remain coarse resolution and do not resolve high frequency forcing.

When modelling biogeochemical cycles of nutrients and carbon in the ocean, we are starting from a position of strength provided by NPZD models that in general do a reasonable job at capturing distributions of chlorophyll, primary production and nutrients. Incorporating extra complexity, where achievable and relevant, is undoubtedly the way forward, building new parameterizations in a gradual and progressive manner. But any such model development should be tempered by a healthy dose of scepticism regarding model outcomes, the tendency being that what little validation is done will overstate the case for belief in model results (Mulligan and Wainwright, 2004). Let us remember that there is an optimum level of complexity for ecosystem models beyond which any potential gain as a result of additional articulation is outweighed by the penalty of lowered accuracy due to inadequate parameterization (Constanza and Sklar, 1985; Anderson and Totterdell, 2004). A somewhat sobering thought is that perhaps this optimum may in some instances not reach the level of articulation required to accurately predict system feedbacks of interest, at least not without further experimental research to improve our knowledge base first.

ACKNOWLEDGEMENTS

TRA is funded by the Natural Environment Research Council, UK. I thank my many colleagues and friends in the scientific community for helpful discussions and critique of the manuscript, as well as four referees for their comments.

REFERENCES