Potential increasing dominance of heterotrophy in the global ocean

The metabolic theory of ecology (Brown et al 2004) predicts that the exponential relationship existing between temperature and metabolism (Gillooly et al 2001) of the individual is scalable to ecosystems. Marine planktonic communities are no exception. It has been demonstrated that the differential temperature dependencies of phytoplankton community respiration (a heterotrophic process) and primary production (an autotrophic process), of which the dependency in respiration is stronger (Duarte and Agusti 1998, López-Urrutia et al 2006, Regaudie-de Gioux and Duarte 2012), determines the spatial and temporal metabolic balance of the global surface ocean (Duarte and Agusti 1998, Hoppe et al 2002, López-Urrutia et al 2006). Growth of phytoplankton, however, is also determined by nutrient limitation, and most phytoplankton communities in the modern ocean are resource and not temperature limited (Marañón et al 2014). In the modern surface ocean, regions of net heterotrophy (where community respiration exceeds gross primary production) are compensated for by regions of high primary production (Duarte and Agusti 1998, Hoppe et al 2002). However, whether zones such as the oligotrophic subtropical gyres are truly net heterotrophic remains an ongoing debate (e.g., Ducklow et al 2013, Duarte et al 2013, Williams et al 2013). The attenuation of the vertical flux of particulate organic carbon is likewise strongly influenced by temperature, which regulates remineralization (bacterial respiration) rates below the euphotic zone (Marsay et al 2015).

Regions with high rates of community respiration tend to be classified as heterotrophic, and also act as net sources of carbon to the atmosphere (Duarte and Agusti 1998). There are obvious implications for the global ocean carbon sink in a warmer world, if higher near-surface remineralization rates lead to a reduction in oceanic carbon uptake, producing a positive feedback on greenhouse warming and 'short-circuiting' of the carbon cycle (López-Urrutia et al 2006, Regaudie-de Gioux and Duarte 2012, Marsay et al 2015). Indeed, isotopic reconstructions of remineralization profiles during the warmer Eocene epoch (55.5–33.7 Ma) reveal both shallower remineralization and higher metabolic rates than found in the modern ocean (Olivarez Lyle and Lyle 2006, John et al 2013), and are suspected of contributing to the maintenance of higher global temperatures over the period (Olivarez Lyle and Lyle 2006). Earth system modelling of the Eocene using temperature-dependent remineralization rates by John et al (2014) recently achieved good agreement with the isotopic reconstructions, supporting the hypothesis that enhanced metabolic rates played a critical role in reducing observed carbon fluxes.

Major changes are currently underway in open ocean biogeochemistry owing to anthropogenic additions of greenhouse gases to the atmosphere. These changes are affected through three independent but synergistic phenomena; ocean warming, acidification and deoxygenation (summarized by Gruber 2011). Recent interest in the potential net biogeochemical impacts of these phenomena has resulted in several studies (e.g., Schmittner et al 2008, Steinacher et al 2010, Bopp et al 2013, Moore et al 2013) examining the combined biogeochemical response in coupled Earth system models. The common biogeochemical metrics in these studies are change in net primary production (NPP), carbon export to the deep ocean, and seawater oxygen content. There is general agreement across models and modelling studies that increasing atmospheric CO₂ concentrations over the 21st century will increase sea surface temperatures (SSTs) and stratification, and decrease carbon export production (Schmittner et al 2008, Steinacher et al 2010, Bopp et al 2013, Moore et al 2013, Yool et al 2013). Simulated 21st century NPP also generally decreases, although Schmittner et al (2008) found an increase using the University of Victoria Earth System Climate Model (UVic ESCM). Although model structures vary and regional differences exist, the main reasons for the decline in NPP and export production are similar and related to increasing restrictions on autotrophic primary production. Increased water column stratification in the low and middle latitudes diminishes resupply of nutrients to the euphotic zone from the deep ocean (Bopp et al 2001). Global NPP declines as nutrient limitation increases in these regions (with declining carbon export as a consequence), even though high latitudes might show increases in NPP as light and temperature limitation recedes (e.g., Schmittner et al 2008, Steinacher et al 2010, Yool et al 2013).

Longer-term ocean biogeochemical trends under a high carbon emissions scenario have been speculated to follow those of the 21st century (e.g., Moore et al 2013). However, Schmittner et al (2008) demonstrated the combination of reduced stratification and warming-enhanced biological processes after the 21st century result in large increases in global NPP beyond 2100, though carbon export still declines in their study. Taucher and Oschlies (2011) also demonstrated in the UVic ESCM an increase in temperature-dependent microbial fast recycling by 2100 that increases temperature-dependent NPP, which does not occur when these processes are considered temperature-independent. A recent update of the UVic ESCM biogeochemistry model by Keller et al (2012) refined the phytoplankton growth rate and added new formulations for grazing and iron limitation. These modifications have resolved the discrepancy between the sign of 21st century NPP trends in the Schmittner et al (2008) study and other model studies, though the longer-term increases in NPP remain (described in section 3). Kvale et al (2015) recently implemented a calcifying phytoplankton functional type and prognostic particulate calcium carbonate (CaCO₃) in the model, which increased the ecological dynamism and improved the mechanistic realism of the carbon export processes. These modifications impact the model's behavior under anthropogenic forcing but do not eliminate the long-term trend in NPP (also described in section 3). Here we argue that the apparently robust simulated trend reversal in NPP in the UVic ESCM must be given serious consideration in light of evidence presented by paleoproxy and modern observations that phytoplankton metabolic processes have the capacity to short-circuit the marine carbon cycle.

We compare the transient ocean biogeochemical response of three University of Victoria Earth System Climate Model (UVic ESCM; Weaver et al 2001, Meissner et al 2003, Eby et al 2009) configurations (MIXED, CAL, and NOCACO3TR). The UVic ESCM is a coarse-resolution (1.8° × 3.6° × 19 ocean depth layers) ocean-atmosphere-biosphere-cryosphere-geosphere model. The standard UVic ESCM version 2.9 includes a NPZD submodel with 'mixed phytoplankton' and 'diazotroph' phytoplankton functional types and one 'zooplankton' functional type (described in Schmittner et al 2008). The MIXED model is the standard UVic ESCM NPZD submodel, updated with the changes made by Keller et al (2012) and referred to as 'NOCAL' in Kvale et al (2015). The CAL model builds upon MIXED and additionally includes a 'small and calcifying phytoplankton' functional type (referred to throughout this paper simply as 'calcifiers') and a prognostic CaCO₃ tracer (Kvale et al 2015). The NOCACO3TR model only includes the additional 'small and calcifying phytoplankton' functional type but does not include the prognostic CaCO₃ tracer. Table 1 lists relevant biogeochemical model tracers and functional types for the three model configurations used here. These configurations were chosen for this study because of hypothesized nutrient-driven shifts in future community composition between 'small' (commonly assumed as a proxy for calcifiers) and 'large' (commonly assumed as a proxy for diatoms and other non-calcifiers) phytoplankton under anthropogenic climate change forcing (e.g., Cermeno et al 2008, Marinov et al 2010, 2013). Shifts in relative abundance of phytoplankton functional types have the potential to alter carbon export to the deep ocean (Bopp et al 2005), for which carbonate from calcifiers is a significant vector (Jin et al 2006). Please refer to Kvale et al (2015) for a complete description of the CAL and MIXED ecosystem and biogeochemical models.

The UVic ESCM phytoplankton biomass equations utilize a modified Eppley curve (Eppley 1972) to calculate maximum possible growth rate as a function of seawater temperature and iron availability. Growth limitation by nitrate and phosphate are calculated using half saturation constants and nutrient concentration, where calcifiers in the CAL and NOCACO3TR models are assigned nutrient affinities higher than mixed phytoplankton (Le Quéré et al 2005). Actual growth rates are determined by the most limiting factor, given the maximum possible growth rate and availability of nutrients and light. Biomass loss terms include a temperature-independent linear mortality term, grazing by zooplankton, and a loss term referred to as 'microbial fast recycling'. There are three temperature-dependent heterotrophic loss terms in the UVic ESCM: microbial fast recycling, zooplankton excretion, and detrital remineralization. The microbial fast recycling term accounts for near-surface losses to the microbial loop and dissolved organic matter cycling. Similar processes govern detrital remineralization and both are parameterized using the Eppley curve. Zooplankton excretion is a proportion of zooplankton grazing, which is subject to prey availability and thus phytoplankton growth rate, as well as the zooplankton growth rate which is temperature dependent until capped at $20\;^\circ {\rm{C}}$ (Keller et al 2012). Respiration due to metabolic maintenance is not explicitly accounted for, but is implicitly included in primary production. Redfield stoichiometry is used throughout.

Organic carbon export in all model configurations occurs through a prognostic detrital tracer made up of dead plankton that sinks at a rate that increases linearly with depth and remineralizes at a rate dependent on temperature. Any organic detritus that reaches the sediments is assumed to dissolve instantly back into the water column. In the CAL model, a portion of this organic detritus is 'protected' from remineralization by a fixed fraction of the associated particulate inorganic CaCO₃ as a rough parameterization of ballasting (Klaas and Archer 2002). This protected detritus loses its protection at the rate of CaCO₃ dissolution, which is a function of the degree of carbonate saturation. CaCO₃ is produced at a fixed ratio based on dead plankton (mixed phytoplankton and zooplankton in MIXED, calcifiers and zooplankton in CAL and NOCACO3TR). It is exported and dissolved instantaneously in MIXED and NOCACO3TR (Schmittner et al 2008, Keller et al 2012), but it is traced and assigned prognostic sinking and dissolution rates in CAL (Kvale et al 2015). In all model configurations, any CaCO₃ reaching the sediments is buried.

The three model versions were first integrated for ten thousand years in fully coupled mode using a constant pre-industrial CO₂ atmospheric concentration of 283.8 ppm to establish equilibrium. Historical CO₂ forcing, as well as historical agricultural, volcanic, sulphate aerosol and CFC emissions, and changes to land ice and solar forcing were then applied from year 1800 to 2005 using the PMIP3 (Paleoclimate Modelling Intercomparison Project Phase 3) data compilation (Machida et al 1995, Battle et al 1996, Etheridge et al 1996, 1998, Flückiger et al 1999, 2004, Ferretti et al 2005, Meure et al 2006). From year 2005 to 2400 the models were forced using increasing CO₂ concentrations and radiative forcing from all non-CO₂ greenhouse gases, fractions of the land surface devoted to agricultural uses, and the direct effect of sulphate aerosols as an alteration of the surface albedo following 'business-as-usual' RCP scenario 8.5 (RCP8.5, Riahi et al 2007, Meinshausen et al 2011, figure 1). Weathering fluxes were held constant at the pre-industrial rate (Meissner et al 2012). The model was driven by seasonal variations in solar insolation at the top of the atmosphere and seasonally varying wind stress and wind fields (Kalnay 1996). The wind fields are geostrophically adjusted to air temperature anomalies (Weaver et al 2001).

Zoom In Zoom Out Reset image size

The physical response (i.e., changes in temperature, salinity, ocean circulation, sea ice) to radiative forcing is the same across all three model configurations. Atlantic meridional overturning rapidly weakens from about 22 Sv (10⁶ m³ s⁻¹) before year 2000 to less than 8 Sv by year 2300 (not shown). Antarctic Bottom Water formation increases slightly after year 2000 by 1 Sv (not shown). Zonally averaged SSTs warm up to $11{\;}^{\circ }{\rm{C}}$ relative to the year 1800 value by year 2400, and strong near-surface stratification occurs on a global scale (figure 1). Differences in model response to climate forcing demonstrate the effect of biogeochemical model structure only, such as presence or absence of calcifiers (compare MIXED to NOCACO3TR) and instant CaCO₃ dissolution and sinking versus full carbonate tracer (compare NOCACO3TR to CAL).

3.1. Changes in plankton biomass

3.2. Changes in biogeochemical pathways

Carbon and nutrient recycling fundamentally shifts after the 20th century in all model configurations. The ratio of NPP to near-surface (to 130 m depth) heterotrophic respiration declines rapidly (figure 2). While all three respiration terms increase over time (figure 2), only fast recycling increases for all three models during the 21st century. While the overall contribution of fast recycling to respiration is small (not shown), it is the most sensitive to 21st century increases in SSTs and dominates the shift towards increasing respiration rates until zooplankton excretion and total water column remineralization rapidly increase in the 22nd century. Fast recycling rates of mixed phytoplankton never decline in spite of an overall decline in NPP (figure 3) and detritus production (figure 4) in the first half of the MIXED simulation. The accelerating effect of microbial fast recycling on increasing global NPP and (decreasing the NPP): (heterotrophic respiration) ratio is most apparent in the NOCACO3TR simulation. Calcifiers replace mixed phytoplankton in the lower latitudes, which leads to more fast recycling and higher nutrient availability. More nutrients allow for higher NPP (including more calcifiers), which maintains the zooplankton population in these regions throughout the simulation. This allows zooplankton excretion rates to increase faster in NOCACO3TR than in the other models. A mitigating effect by CaCO₃ ballast is apparent in the CAL simulation when compared to NOCACO3TR, where eventual increases in all heterotrophic rates are not as large because of the protection from remineralization provided to a portion of sinking detritus.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

3.3. Changes in NPP

3.4. Changes in carbon export

3.5. Potential caveats

The experiments described here highlight the importance of biogeochemical model parameterizations for the determination of the behavior of the ocean carbon cycle during climate transitions. Though each model analysed in this study (MIXED, NOCACO3TR, and CAL) has very similar equilibrium states, their respective biogeochemical behaviors under transient forcing are unique both in global and regional trends. Overall, the models share a common response of a transition between two dominant regimes. The first regime is physically-dominated and regulates autotrophy for lower changes in SSTs. Reductions in nutrient resupply from the deep ocean cause mixed phytoplankton populations to decline in the low and middle latitudes, while increased warming and reduced sea ice cover increase the mixed phytoplankton population in high latitudes. In the MIXED model, global NPP, particulate organic carbon and CaCO₃ export production decline over the first half of the simulation as a result of the low latitude trend. Models including an additional plankton functional type representing small phytoplankton and calcifiers (NOCACO3TR and CAL) show a partial replacement of mixed phytoplankton by calcifiers in the stratifying low latitudes due to their higher nutrient affinities. Consequently, global NPP, particulate organic carbon and CaCO₃ export production change less. All models start to transition into the second regime around the start of the 21st century, when temperature-dependent heterotrophic respiration rates increase relative to NPP.

During the second regime a biologically-driven feedback on NPP dominates the physically-driven one that defined the 21st century. It is driven by warming-enhanced heterotrophic respiration rates, which compensate for reductions in nutrient resupply by rapidly recycling carbon and nutrients in the near-surface. This feedback increases global NPP and CaCO₃ export production in all model configurations and represents a fundamental ecological shift towards increasing heterotrophy. Microbial fast recycling in the near-surface is the most sensitive to warming and accelerates the transition to this second regime. The simulated transition to this second regime in all model configurations raises the possibility trends observed in NPP over recent decades cannot be extrapolated far into the future. The net carbon export into the deep ocean during this second regime is dependent on whether the model accounts for detrital ballasting by calcite. The model that contains a ballasting parameterization (CAL) shows a mitigated increase in near-surface particulate organic carbon production, and a mitigated decline in deep ocean particulate organic carbon export because the calcite protects particulate organic carbon from near-surface remineralization and continues to effectively remove it to the deep ocean throughout the simulations. No CMIP5 model currently accounts for the ballasting effect, which suggests they are missing feedbacks that are important to long-term trends in particulate organic carbon export.

Including small phytoplankton and calcifiers as a distinct phytoplankton functional type causes the model to become more sensitive to warming-enhanced remineralization because small phytoplankton maintain NPP in the warmest (low and middle latitude) regions. A $2{\;}^{\circ }{\rm{C}}$ change is a well-established threshold for dangerous climate change (e.g., Graßl et al 2003), and this analysis adds another reason to support this guardrail. This analysis also suggests that the low to middle latitude NPP response is going to be one of the more important indicators of system sensitivity to observe as ocean change continues, and that accurately reproducing phytoplankton competition in models is of utmost importance for understanding tipping points and dangerous thresholds in the global carbon cycle. Given the impact this feedback has on global NPP, further investigation of temperature-dependent heterotrophy in models is warranted. In particular, to what degree a transition between 'physical' and 'biological' regimes represents an ecological or climatological tipping point with a demonstrable hysteresis (Duarte et al 2012) is an open question, better answered with a model that resolves autotrophic respiration and dissolved organic matter explicitly, or multiple size classes of zooplankton.

Modellers tend to shy away from discussion of longer-term and larger changes in their models but we argue that increasing heterotrophy might be possible because (1) the parameterizations used here are based on measurements of phytoplankton for a range of temperatures from 0–40${\;}^{\circ }{\rm{C}}$ (Gillooly et al 2001), the upper limit of which is never exceeded in our simulations, (2) two recent publications suggest a temperature threshold into net heterotrophy for Arctic ecosystems within 5–6${}^{\circ }{\rm{C}}$ (Duarte et al 2012, Holding et al 2013), which roughly agrees with our simulated global response, and (3) Eocene proxy reconstructions off Tanzania and eastern Mexico reveal a transition to high rates of production and low rates of export at SSTs within $4\;^\circ {\rm{C}}$ warmer than present day (John et al 2013), though abyssal temperatures were likely $10\;^\circ {\rm{C}}$ warmer (Olivarez Lyle and Lyle 2006, John et al 2013).

The potential for a rapid global transition from a physically dominated to a biologically dominated biogeochemical regime in the not-too-distant future urges an improved quantification of the associated uncertainties. The sensitivity of the transition to physical model response and biogeochemical model formulation is better addressed with a suite of models including all of the major carbon cycle components and phytoplankton functional types included here. Optimality-based models, e.g., Pahlow et al (2013), are not limited by fixed stoichiometries and could offer particular insight into systemic sensitivities. Broadly similar marine ecosystem models have been shown to produce a range of physical and biological drivers of NPP trends over the 21st century (Laufkötter 2015), and it is therefore logical to assume the timing and severity of a global shift in ecosystem regimes would likewise be influenced by a range of drivers. The authors would welcome simulations of other models extending beyond 2100 and a systematic examination of the associated drivers of NPP.

This work was supported by an award under the Merit Allocation Scheme on the NCI National Facility at the ANU. KJM is grateful for support under the ARC Future Fellowship programme (FT100100443). KFK is grateful for support from UNSW through a University International Postgraduate Award, and the ARC Centre of Excellence for Climate System Science. The authors would also like to thank Andreas Schmittner for his thoughtful commentary on an early version of the manuscript.