Articles | Volume 9, issue 2
Research article
 | Highlight paper
13 Jun 2018
Research article | Highlight paper |  | 13 Jun 2018

Hazards of decreasing marine oxygen: the near-term and millennial-scale benefits of meeting the Paris climate targets

Gianna Battaglia and Fortunat Joos

Ocean deoxygenation is recognized as key ecosystem stressor of the future ocean and associated climate-related ocean risks are relevant for current policy decisions. In particular, benefits of reaching the ambitious 1.5 C warming target mentioned by the Paris Agreement compared to higher temperature targets are of high interest. Here, we model oceanic oxygen, warming and their compound hazard in terms of metabolic conditions on multi-millennial timescales for a range of equilibrium temperature targets. Scenarios where radiative forcing is stabilized by 2300 are used in ensemble simulations with the Bern3D Earth System Model of Intermediate Complexity. Transiently, the global mean ocean oxygen concentration decreases by a few percent under low forcing and by 40 % under high forcing. Deoxygenation peaks about a thousand years after stabilization of radiative forcing and new steady-state conditions are established after AD 8000 in our model. Hypoxic waters expand over the next millennium and recovery is slow and remains incomplete under high forcing. Largest transient decreases in oxygen are projected for the deep sea. Distinct and near-linear relationships between the equilibrium temperature response and marine O2 loss emerge. These point to the effectiveness of the Paris climate target in reducing marine hazards and risks. Mitigation measures are projected to reduce peak decreases in oceanic oxygen inventory by 4.4 % C−1 of avoided equilibrium warming. In the upper ocean, the decline of a metabolic index, quantified by the ratio of O2 supply to an organism's O2 demand, is reduced by 6.2 % C−1 of avoided equilibrium warming. Definitions of peak hypoxia demonstrate strong sensitivity to additional warming. Volumes of water with less than 50 mmol O2 m−3, for instance, increase between 36 % and 76 % C−1 of equilibrium temperature response. Our results show that millennial-scale responses should be considered in assessments of ocean deoxygenation and associated climate-related ocean risks. Peak hazards occur long after stabilization of radiative forcing and new steady-state conditions establish after AD 8000.

1 Introduction

Oxygen (O2) is a sparingly soluble gas and its abundance in the ocean is decreasing under ongoing global warming (IPCC2013). Decreasing O2 concentrations, warming and changes in other environmental parameters forced by anthropogenic greenhouse gas (GHG) emissions pose high risks for marine ecosystems (Gattuso et al.2015). The parties to the United Nations Framework Convention on Climate Change (UNFCCC2018) and to the Paris Agreement note “the importance of ensuring the integrity of all ecosystems, including Oceans” (Paris Agreement2018) and “the threats of irreversible damage” (UNFCCC2018). Marine changes are projected to evolve over multi-century and millennial timescales with peak impacts occurring potentially long after stabilization of atmospheric GHG concentrations and peak temperatures. However, only a few studies have assessed millennial-scale impacts of anthropogenic GHG emissions on the ocean and the reversibility of marine changes in oxygen. Explicit quantification of the benefits of meeting the 2 or 1.5 C climate targets mentioned by the Paris Agreement with respect to the reversibility and avoidance of implied impacts on marine oxygen and related environmental parameters, including ocean circulation, ocean warming, metabolic viability and biological productivity on multi-millennial timescales is missing.

Typical thresholds leading to O2 stress for many macro-organisms (hypoxia) are around 50 mmol O2 m−3. Water with lower O2 concentrations are effectively dead zones for many higher animals (Breitburg et al.2018; Keeling et al.2010; Storch et al.2014). Species are also sensitive to thermal stress (Gattuso et al.2015) and their sensitivity to hypoxia increases with higher temperatures (Pörtner2010). In the modern ocean, oxygen-poor zones with O2< 50 mmol m−3 occupy about 5 % of its volume (Bianchi et al.2012; Garcia et al.2014). The expansion of oxygen-poor waters leads to habitat compression, mortality and major changes in community structure where energy preferentially flows into microbial pathways to the detriment of higher trophic levels. Suboxic (< 5 mmol O2 m−3) or anaerobic conditions can also lead to production of poisonous H2S within sediments (Breitburg et al.2018; Diaz and Rosenberg2008), and decreasing O2 concentrations potentially lead to higher production and emissions of the greenhouse gas nitrous oxide.

Observational (Ito et al.2017; Schmidtko et al.2017) and modeling studies (Oschlies et al.2017) indicate an overall decline in the oceanic oxygen content over past decades. Systematic discrepancies exist for the typically low-oxygen tropical thermocline, where observations suggest O2 has decreased and most models simulate increased O2 levels over the past decades. Model projections to the end of the 21st century consistently project the global ocean oxygen inventory to further decline with anthropogenic climate change (Bopp et al.2002, 2013; Cocco et al.2013; Frölicher et al.2009; Matear2000; Plattner et al.2001). The most recent generation of Earth system models simulate global deoxygenation by the end of the 21st century of around 1.81 % (RCP2.6) to 3.45 % (RCP8.5; IPCC2013). Impact studies have highlighted potential habitat compression (Deutsch et al.2015; Mislan et al.2017) and reduced catch potential (Cheung et al.2016) associated with climate change at the end of the century. Large model–model differences remain in projections of oxygen minimum zones (OMZs; Bopp et al.2013; Cocco et al.2013).

Given the long residence time of anthropogenic CO2 in the atmosphere and long equilibration timescales of the ocean overturning circulation, anthropogenic climate change will grow and persist beyond the end of the 21st century, the typical near-term assessment timescale of climate change (Clark et al.2016). Only a few studies have assessed ocean biogeochemistry and the oceanic oxygen content beyond this near-term timescale. Available studies employ a range of physical and biogeochemical complexity levels from box models to general circulation models (GCMs). Oxygen concentrations are simulated to decline beyond the 21st century on multi-centennial timescales (Hofmann and Schellnhuber2009; Matear and Hirst2003; Mathesius et al.2015). Simulations covering 2 millennia show a recovery phase thereafter (Schmittner et al.2008; Yamamoto et al.2015). In most studies, simulated oxygen concentrations have not reached new steady-state conditions at the end of the simulation. Low-order Earth system models and those of intermediate complexity (EMIC) integrated by up to 100 000 years have demonstrated the potential for long-term ocean oxygen depletion in response to carbon dioxide emissions and the long equilibration timescales of ocean biogeochemical variables in response to carbon emissions (Ridgwell and Schmidt2010; Shaffer et al.2009). Multi-millennial simulations are therefore required to assess the full amplitude of ocean biogeochemical changes and new steady-state conditions due to anthropogenic climate change.

The distribution of O2 in the ocean results from the sum of its solubility component set through air–sea exchange, the effect of O2 production by phytoplankton in the euphotic zone and O2 consumption during organic matter remineralization at depth. In modeling studies, it is possible to identify the drivers of O2 changes by considering changes due to solubility and changes due to oxygen consumption. When assessing the near-term timescale at the end of 21st century, studies have shown that different depths in the water column tend to be associated with different dominant mechanisms of change. In the surface ocean, O2 changes tend to be determined by changes in O2 solubility. In the subsurface, both changes in solubility and utilization may reinforce (mid- and high latitudes) or compensate for each other (tropics; e.g., Bopp et al.2017; Cabre et al.2015). In the deep ocean, simulated O2 changes are dominated by changes in O2 utilization, which is in turn controlled by ocean ventilation (Matear and Hirst2003; Yamamoto et al.2015). Changes in the oceanic heat content and in ocean circulation are therefore crucial for O2 changes.

Well-defined metrics that summarize the Earth system response are useful in many aspects and may facilitate communication in the mitigation policy context of the Paris Agreement. The Transient Climate Response to Cumulative Carbon Emissions (Allen et al.2009) or the Transient Earth System Response to Cumulative Carbon Emissions (Steinacher and Joos2016) are such metrics. These link changes in global surface air temperature and environmental parameters to cumulative CO2 emissions relying on near-linear relationships. Similar metrics may be developed for oceanic oxygen.

Deoxygenation is one of several marine ecosystem stressors including warming, acidification, hypocapnia, changes in food supply and sea-level rise (Bopp et al.2013; Cocco et al.2013; Gattuso et al.2015; Gruber2011; Sweetman et al.2017). Several key marine and coastal ecosystems may face severe risks due to climate change even if the suggestions for low emission pathways are followed through the end of the century (Breitburg et al.2018; Gattuso et al.2015; Magnan et al.2016). This growing body of concern helped drive the 21st Conference of the Parties (COP21) to reach the Paris Agreement. A goal of the Paris Agreement is to “hold the increase in the global average temperature to well below 2 C above pre-industrial levels and pursuing efforts to limit the temperature increase to 1.5 C above pre-industrial levels” (Paris Agreement2018). Article 4 of the Paris Agreement further states that “[in] order to achieve the long-term temperature goal set out in Article 2, Parties aim to reach global peaking of greenhouse gas emissions as soon as possible … and to undertake rapid reductions thereafter in accordance with best available science, so as to achieve a balance between anthropogenic emissions by sources and removals by sinks of greenhouse gases in the second half of this century” (Paris Agreement2018).

In this study, we assess the effectiveness of the Paris climate targets in reducing hazards of decreasing oceanic oxygen, ocean warming and marine export productivity as simulated by the Bern3D Earth System Model of Intermediate Complexity. We prescribe in the model four different, idealized scenarios where anthropogenic GHG forcing is stabilized by AD 2300. The four scenarios are designed to reach an equilibrium warming of 1.5, 1.9, 3.3 and 9.2 C above preindustrial. Simulations are run to year AD 10 000 by which time the ocean has reached new steady-state conditions. This allows us to assess reversibility and the full amplitude of changes, acknowledging the long equilibration timescale of biogeochemical variables with peak hazards potentially occurring long after stabilization of radiative forcing in the atmosphere. We summarize the outcomes developing global metrics which quantify avoided marine hazards per avoided global warming on three different time horizons. The first time horizon is the end of the 21st century, the typical assessment timescale of climate change hazards. Here, changes in a variable are related to changes in surface air temperature (SAT) in the year 2100. Those are contrasted to the millennial-scale perspective where peak changes in the variable in the course of the simulation and equilibrium changes at the end of the simulation are related to the corresponding equilibrium warming.

In Sect. 2, we briefly describe the Bern3D model and the experimental setup. Four different radiative forcing stabilization scenarios to meet four temperature targets (1.5, 1.9, 3.3 and 9.2 C above preindustrial) are considered. The observation-constrained 100-member ensembles used to explore parameter uncertainties for each scenario is introduced. In Sect. 3, physical changes, including changes in overturning, water mass age, sea ice, temperature, salinity and density as well as biogeochemical changes, including changes in global oxygen inventory, the extent of oxygen minimum zones and productivity, are presented. The compound effects of warming and oxygen changes are assessed in the form of a metabolic index (Deutsch et al.2015). Underlying physical and biogeochemical processes and mechanisms are discussed. Following earlier studies, we attribute the contributions of O2 changes from changes in solubility, and the interplay of ocean biology and ventilation by carrying four explicit O2 tracers and an ideal age tracer. The graphical illustration of spatial changes is focused on the 1.5 C equilibrium warming target at the point of peak O2 decline. Additional supporting figures are given in the appendix. In Sect. 4, the relationship between change in global mean surface air temperature (ΔSAT) and selected impact-relevant parameters is quantified. Relationships are established for the near-term (AD 2100), the time of the peak decline in oxygen around AD 3000 to 4000 and at year AD 10 000 when a new equilibrium has been reached in the model. Often a near-linear relationship is found between the change in a variable of interest and the change in SAT as simulated across the range of scenarios and ensemble members at a distinct time. This allows us to develop new metrics to quantify avoided marine hazards per unit change in ΔSAT at different points in time. These quantitatively illustrate the benefits of meeting the Paris target in terms of marine hazards. Each modeling exercise is associated with uncertainties, and in Sect. 5, we discuss relevant uncertainties, mention neglected processes and compare our findings to other studies. Finally, in Sect. 6 we present implications and conclusions, and we summarize our findings graphically for a “1.5 C world” and contrast peak changes across the range of temperature targets.

2 Model and simulations

2.1 Bern3D

The Bern3D Earth System Model of Intermediate Complexity is a three-dimensional frictional geostrophic balance ocean model (Müller et al.2006), which includes a sea-ice component coupled to a single-layer energy and moisture balance model of the atmosphere (Ritz et al.2011) and a prognostic marine biogeochemistry module (Parekh et al.2008; Tschumi et al.2011). A version with 41 × 40 horizontal grid cells and 32 vertical layers is used (see also Battaglia et al.2016; Roth et al.2014 for model evaluation). The NCEP/NCAR monthly wind stress climatology (Kalnay et al.1996) is prescribed at the surface. Air–sea gas exchange, carbonate chemistry and natural Δ14C of dissolved inorganic carbon is modeled according to OCMIP-2 protocols (Najjar et al.1999; Orr and Epitalon2015; Orr and Najjar1999). The global mean air–sea transfer rate is reduced by 19 % compared to OCMIP-2 to match observation-based estimates of natural and bomb-produced radiocarbon (Müller et al.2008).

The biogeochemical module is based on phosphorus and simulates production and remineralization/dissolution of organic matter, calcium carbonate and opal. Production of particulate organic phosphorus (POP) within the euphotic zone (top 75 m) depends on temperature, light availability, phosphate and iron following Doney et al. (2006). POP remineralization within the water column follows a power law profile (Martin et al.1987). Organic matter falling on to the sea floor is remineralized in the deepest box. Two-thirds of organic phosphorus production form dissolved organic phosphorus (DOP), which decays with an e-folding lifetime of 1.5 years. An updated remineralization scheme assigns remineralization of POP and DOP to aerobic and anaerobic pathways depending on the mean grid cell dissolved O2 concentration (see Battaglia and Joos2018). We introduce two power law profiles with two distinct remineralization length scales for aerobic and anaerobic remineralization (αaerob and αdenit). Constant stoichiometric ratios are used for both aerobic and anaerobic remineralization to convert biological P fluxes into carbon and alkalinity fluxes (P : Alk : C = 1 : 17 : 117). The O2 demand for complete aerobic remineralization is 170 molO2molPO4 and no oxygen is consumed for anaerobic remineralization. Accordingly, aerobic remineralization in the ocean is smaller than O2 production in the euphotic zone leading to an O2 outgassing for steady-state conditions. The atmospheric oxygen inventory is constant. This is justified as 99.5 % of the ocean–atmosphere inventory is in the atmosphere and potential net fluxes of O2 from the ocean and land to the atmosphere and fossil fuel burning have a small impact on atmospheric O2. O2 components from O2 production, consumption and solubility are carried as explicit model tracers to attribute changes. Tracers add up to within 10−14 Pmol with mean inventories of 23.2, 239, 430 yielding a total of 214 Pmol, respectively (median values given). O2 components inferred from O2 saturation can result in systematic errors from surface disequilibrium (Ito et al.2004). The use of explicit tracers avoids such systematic errors in the O2 components. As changes in the O2 production term are small, we combine the O2 production and consumption tracers to a O2 biology tracer when displaying sections.

We include evaluation of a metabolic index, Φ, which was proposed by Deutsch et al. (2015). It combines temperature and pO2 as indicators of metabolically viable environments and is defined as the ratio of O2 supply to an organism's resting O2 demand. We consider only relative changes in Φ(t) relative to a reference time, t0 (average over 1870–1899):

(1) Δ Φ ( t ) Φ t 0 = p O 2 ( t ) p O 2 t 0 × exp E 0 k B 1 T ( t ) - 1 T t 0 - 1 ,

where T is the absolute temperature, kB is Boltzmann's constant and the exponential function and the parameter E0 characterize the temperature dependence of the baseline metabolic rate. E0 only weakly affects the relative influence of temperature and O2 gradients, and relative changes in Φ are therefore independent of species (Deutsch et al.2015). Here, we consider E0= 0.87 eV (representative of Atlantic cod). For the calculation of pO2 we pressure-correct the equilibrium constant following Eq. (5) in Weiss (1974). The metabolic index Φ, as proposed by Deutsch et al. (2015), is linear in pO2 (representing the rate of O2 supply) and decreases non-linearly with temperature (indicative of the resting metabolic demand). One may note that the exponential curve varies approximately linearly for typical global warming associated temperature changes as E0kb ( 10 000 K) is large.

The current set up does not include sediment interactions, temperature-dependent remineralization, variable stoichiometry, nitrogen cycle feedbacks, atmospheric nutrient deposition, dynamic wind or freshwater input/albedo changes from melting of continental ice sheets.

2.2 Ensemble and scenarios

To explore potential oxygen changes we set up four 100-member ensembles each targeting a different equilibrium temperature response ( 1.5, 1.9, 3.3 and 9.2 C above preindustrial). A feedback parameter λ (W m−2 K−1; Ritz et al.2011), accounting for climate feedbacks that are not explicitly treated in the Bern3D model, is chosen in combination with radiative forcing from the Representative Concentration Pathways (Meinshausen et al.2011) to achieve these stabilization targets. RCP2.6, stabilizing by 2300, is run with λ values of 0.71 and 1 W m−2 K−1, reaching the 1.5 and 1.9 C targets, respectively. RCP4.5, stabilizing after 2100, is run with 1 W m−2 K−1 yielding a 3.3 C temperature response and RCP8.5, stabilizing in the 23rd century, with 0.71 W m−2 K−1 yielding a 9.2 C response (median values given for temperature targets). Each member is spun up over 5000 years to AD 1765 boundary conditions. The radiative forcing follows RCP scenarios (Meinshausen et al.2011). The RCP scenarios are extended to year AD 10 000 by which time the ocean has equilibrated to new steady-state conditions. Radiative forcing includes an 11-year solar cycle up to year AD 3000. After that, all forcings are kept constant. We employ a single-model setup, and assess uncertainties arising from organic matter remineralization (αaerob and αdenit) and vertical mixing (kdiff−dia). The three parameters are sampled using the Latin hypercube sampling technique (McKay et al.1979). The parameter ranges are chosen such that all members achieve similar skill scores with respect to observation-derived fields of natural radiocarbon (Key et al.2004) and dissolved O2 (Bianchi et al.2012; Garcia et al.2014) and correspond to the values chosen in Battaglia and Joos (2018). A normal distribution is used to sample αaerob with a standard value of 0.83 and a standard deviation of 0.0625. αdenit is sampled uniformly between 0.1 and 0.01. And a lognormal distribution is used to sample kdiff−dia (standard value = 2.25 × 10−5 m2 s−1, shape parameter = 0.2, location parameter = 0). We choose a single ensemble member with parameter values close to the standard values as representative ensemble member to illustrate spatial anomalies (αaerob=0.85, αdenit=0.037, kdiff−dia= 2.05 × 10−5 m2 s−1).

Figure 1Temporal evolution of physical variables relative to 1870–1899 for model ensembles aiming at 1.5, 1.9, 3.3 and 9.2 C warming targets. Lines mark the median and shading marks the 90 % range of the ensemble. The shading reflects uncertainties due to variations in the diapycnal mixing coefficient. (a) Displays the changes in surface air and (b) displays the changes in ocean mean temperatures. (c,) and (d) are annual mean sea-ice areas in the respective hemisphere. (e) Atlantic meridional overturning is the maximum of the Atlantic meridional overturning stream function below 400 m depth and (f) Indo-Pacific meridional overturning is the minimum of the Indo-Pacific meridional overturning stream function below 400 m depth.


2.3 Pre-industrial characteristics

The ensemble produces a range in overturning strengths, remineralization fluxes and O2 distributions. The following numbers represent the 90 % confidence ranges of important model characteristics across the ensemble. The maximum of the Atlantic meridional overturning circulation (AMOC) stream function below 400 m depth ranges from 16.5 to 19.7 Sv. The minimum of the Indo-Pacific meridional overturning stream function below 400 m depth (Indo-Pacific MOC) ranges between 13.6 and 15.6 Sv. Export of particulate organic matter at 75 m ranges from 9.0 to 11.4 Gt C yr−1. The simulated oxygen inventory ranges between 195 and 230 Pmol given the three parameters and the simulated oxygen distribution covers the observational range and spatial pattern well (Battaglia and Joos2018). Biases exist in the simulated extent of OMZs. The volume of suboxic conditions (O2< 5 mmol m−3) is overestimated by a factor of 5 but water column denitrification fluxes are well within current estimates (Battaglia and Joos2018). This is a common model bias in EMICs and GCMs (Bopp et al.2013; Cabre et al.2015; Cocco et al.2013). Vastly enhanced spatial resolution may be required to simulate equatorial physics and ecosystems in better agreement with observations (Bopp et al.2013).

3 Marine changes in temperature, circulation and biogeochemistry

We first describe the evolution of important physical quantities that impact O2 concentrations. Figure 1 displays the temporal changes in global mean surface air and ocean temperature, the evolution of annual mean sea-ice area in the Northern and Southern hemispheres and the Atlantic and Indo-Pacific meridional overturning circulations.

Figure 2Temporal evolution of critical variables relative to 1870–1899 for model ensembles aiming at 1.5, 1.9, 3.3 and 9.2 C warming targets. Lines mark the median and shading marks the 90 % range of the ensemble. The shading reflects uncertainties due to variations in the diapycnal mixing coefficient and the aerobic and anaerobic remineralization length scales of particulate organic matter (αaerob and αdenit). (a) is the change in the total oceanic O2 inventory. (b) O2 solubility is the explicitly traced solubility component of oceanic oxygen, (c) is the explicit O2 production tracer, and (e) is the explicit O2 consumption tracer. (h) Oxygen-poor waters are taken as the volume of water with O2< 50 mmol m−3.


Figure 3Changes in O2 and its components at time of peak O2 decline (AD 3150) relative to preindustrial steady state for a single representative ensemble member reaching the 1.5 C warming target. (a) Change in total O2, (b) change in O2 due to biology, (c) change in O2 due to solubility, (d) change in ideal age.


In response to the RCP scenarios, atmospheric temperatures rise and stabilize after  1000 years (Fig. 1a). The four ensembles reaching 1.5, 1.9, 3.3 and 9.2 C above preindustrial surface air temperature show an equilibrium ocean warming of 1.1, 1.3, 2.0 and 5.5 C, respectively (median values given). Sea ice retreats in both hemispheres (Fig. 1c and d). The retreat is more pronounced for higher forcing. In the Southern Hemisphere, even the lower forcing levels show strong decline in the annual mean sea-ice area and sea ice vanishes for higher forcing. The warming perturbation causes the AMOC and Indo-Pacific MOC to decline transiently (Figs. 1e, f and A1). The larger the forcing and implied changes in stratification, the larger the peak decline in overturning (Fig. 1e and f). The decline is likely driven by upper ocean warming, leading to increasing surface-to-deep density gradients as further modulated by salinity changes (Fig. A2). The deep ocean water mass age increases in response to the slowed overturning (Figs. 2d and 3d). As retreating sea-ice increases wind stress over these newly exposed areas, younger water masses form in the upper ocean of the Southern Ocean (Fig. 3d). Reduced convection may contribute to a younger upper ocean through decreased entrainment of old deep water (not quantified within the scope of this paper). As the model tends to equilibrate under the sustained radiative forcing, the surface-to-deep gradients in the density anomalies diminish (Fig. A2), the meridional overturning circulation recovers (Fig. 1e and f) and anomalies in water mass age become again smaller (Fig. 3d versus Fig. A3). The final circulation state is close but not identical to the preindustrial steady-state circulation. Maximum overturning strength in AMOC and the Indo-Pacific MOC varies by less than ±1 Sv around the initial value. At the new steady state, the maximum in AMOC below 400 m tends to be lower under higher forcing, whereas the maximum in the Indo-Pacific MOC below 400 m tends to be higher under higher forcing (Fig. 1e and f). It is difficult and beyond the scope of this paper to conclusively explain such subtle changes in ocean dynamics and overturning (Fig. A1), likely linked to the complex changes in density (Fig. A2) and sea-ice retreat (Fig. 1c and d). However, these differences have direct consequences for the projected global water mass age and by that for oceanic oxygen (Fig. 2) at the new equilibrium as further discussed below.

The response in oceanic oxygen is complex and characterized by an initial decline followed by a recovery phase (Fig. 2a). In line with earlier studies (Matear and Hirst2003; Mathesius et al.2015; Ridgwell and Schmidt2010; Schmittner et al.2008; Shaffer et al.2009), our results demonstrate the potential for large changes in marine oxygen under anthropogenic forcing, a large inertia in the response and a slow and partially incomplete recovery of the perturbation. Transiently, the whole ocean oxygen inventory decreases by a few percent (6 %) under low forcing and by as much as 40 % under high forcing (median values given). The minimum in oxygen occurs about a thousand years after stabilization of radiative forcing, and it takes several millennia to approach a new equilibrium. Then, the global ocean O2 inventory is a few percent higher than at preindustrial conditions under low and intermediate forcing and remains depleted by around 8 % in the high forcing case.

Figure 2 further explains the temporal evolution and interplay of the underlying drivers. In all cases, the changes in global oxygen inventory (Fig. 2a) strongly correlate with water mass age (Fig. 2d) and are also impacted by gradual oxygen loss due to warming as evidenced by the evolution of the O2 solubility tracer (Fig. 2b). Inventory changes based on the O2 production tracer (Fig. 2c) are negligible; changes equilibrate with the atmosphere and only a small fraction remains in the ocean. The O2 consumption tracer (Fig. 2e) determines the shape of the global O2 signal (Fig. 2a). Its decline and recovery phase is strongly correlated with the evolution of ideal age (Fig. 2e and d, see also Fig. 3b and d). As has been shown previously (Bopp et al.2017; Yamamoto et al.2015), the high correlation between changes in O2 and ideal age and the absence of a direct relationship between changes in remineralization fluxes and O2, indicate that circulation changes are the major contributors to changes in O2. Changes from remineralization fluxes include both changes in absolute aerobic remineralization fluxes and changes in the relative share of denitrification (Fig. 2i). An increased share of denitrification at organic matter remineralization, for instance, effectively constitutes an implicit O2 gain. Denitrification fluxes correlate with the volumetric expansion of OMZs and are also impacted by changes in remineralization fluxes within them (Fig. 2i). The recovery level of the O2 consumption tracer (Fig. 2e) reflects the global recovery level of ideal age (Fig. 2d), where younger water masses are associated with less O2 consumption and therefore higher O2 concentrations. The total O2 recovery level (Fig. 2a), on the other hand, is diminished due to O2 loss from solubility (Fig. 2b). As such, 1.5 to 3.3 C warming targets reach similar global O2 equilibrium levels for different reasons. The 1.9 and 3.3 C warming targets tend to result in younger water masses, which would increase O2 due to less O2 consumption compared to 1.5 C warming targets. As those scenarios are also associated with higher warming, they lose more O2 due to less solubility compared to 1.5 C warming targets and yield similar global anomalies despite more pronounced spatial patterns. The 9.2 C warming target reaches a lower equilibrium O2 inventory compared to preindustrial due to high O2 loss from solubility (44.1 Pmol).

Figure 4Changes in potential ecosystem stressors at peak O2 decline (AD 3150) relative to preindustrial steady state for a single representative ensemble member reaching the 1.5 C warming target. Results are displayed for a cross section through the Atlantic (25 W), across the Southern Ocean (58 S) and into the Pacific (175 W). Changes in POM export at 75 m (c) and in surface PO4 concentrations (d) are displayed along the same section.


We illustrate spatial changes in critical variables for a single representative ensemble member (see Sect. 2.2) at its peak O2 decline, which occurs at year AD 3150 and amounts to 5 % (Fig. 3). The member eventually reaches the 1.5 C warming target. Figure 3 displays anomalies in total O2 (Fig. 3a), and the contributions from biologically mediated changes (termed “biology” below, Fig. 3b) combining the changes in the O2 production and consumption tracer and from changes in solubility (Fig. 3c). In the upper ocean O2 concentrations tend to increase due to biology and decrease due to solubility. Such compensating mechanisms have been documented elsewhere (Bopp et al.2017; Cabre et al.2015). The resulting changes in O2 are less pronounced than the changes in each component. The increase in O2 due to biology stems from younger water masses and less export in the low and mid-latitudes (see next paragraph and Fig. 4c). O2 changes show strong spatial correlation with changes in water mass age (Fig. 3a and d). Largest decreases in O2 are simulated in bottom waters in line with older water mass age. The equilibrium response in O2 for this 1.5 C warming case is characterized by slight O2 decreases in the Atlantic, caused mainly by less solubility, and increases in the Southern Ocean and deep Pacific, caused by higher overturning and less sea-ice coverage in the Southern Ocean compared to preindustrial (Fig. A3).

Global export production is simulated to decline over the first few centuries, and reach higher values under new steady-state conditions (Fig. 2g). The decline is stronger for higher forcing, while the recovery level of global export production is similar across the scenarios. Bern3D transiently simulates decreased export in the mid- and low latitudes (Fig. 4c; see also Battaglia and Joos2018; Steinacher et al.2009) as a result of increased stratification (Fig. A2c, f and i) and reduced nutrient concentrations in the surface ocean (Fig. 4b). In the high latitudes, the model simulates increased export production, as a result of less temperature and light limitation as surface waters warm and sea ice retreats. This pattern of decreased export in mid- and low latitudes and increased export in high latitudes is similar across the scenarios. Export production in the low latitudes fully recovers for lower forcing and partially recovers for higher forcing. The lower recovery level in the low latitudes is compensated by higher increases in the high latitudes for high forcing. The magnitude of positive and negative changes increases with forcing, but the global anomalies remain comparable at the end of the simulation.

Next to changes in export, we consider the evolution of a metabolic index in the upper ocean which integrates effects of changes in O2 and temperature at the organism level (Figs. 2f and 4e). The globally averaged, upper ocean (depth < 400 m) metabolic index declines throughout the simulation dominated by increased temperatures (Fig. 2f). The metabolic index, Φ (Deutsch et al.2015), decreases in most places in line with warming and lower pO2 (Fig. 4a, d and e). The O2 gain in upper ocean waters is able to counteract the adverse effect of warming in some high-latitude environments. In other places with higher pO2, the temperature increase dominates the response in Φ. Near-bottom waters in the Pacific are prone to largest reductions in Φ, driven by large decreases in pO2 (Fig. 4e).

Oxygen-poor waters (O2< 50 mmol m−3, Fig. 2h) are simulated to transiently increase across all scenarios. The response is characterized by high uncertainty as introduced by the sampled parameters. Under new equilibrium conditions, the volume of low O2 waters is reduced for low and intermediate forcing and remains higher than pre-industrial levels in the high forcing case.

Figure 5Changes in marine ecosystem stressors versus changes in global mean surface air temperature at three distinct points in time relative to 1870–1899. The colored dots indicate results of the four 100-member ensembles targeted to reach 1.5 (blue dots), 1.9 (light blue dots), 3.3 (orange dots) and 9.2 C (red dots) warming targets. The lines connect results of individual ensemble members at AD 2100 (gray), at time of peak decline of each variable (light blue) and by the end of the simulation (light green) when a new equilibrium state has been reached. Peak and equilibrium changes in variables of interest are related to the corresponding equilibrium temperature response, while changes at the end of the 21st century are related to the transiently realized warming. Peak and equilibrium responses are indistinguishable in the figure for the metabolic index (d) and the ocean temperature change (e). Speak is the peak sensitivity of each variable per C equilibrium warming.


Turning to uncertainties in our perturbed parameter ensemble, we find that variations in the vertical diffusion parameter (kdiff−dia) dominate the uncertainty in the globally averaged evolution of ideal age, sea-ice cover, temperature and O2. The modeled uncertainty in the volume of low O2 waters is dominated by different values of the αaerob parameter. Whether a threshold in O2 concentration is met depends on the pre-industrial tracer distribution. Longer remineralization length scales bring more remineralization to depth, leading to higher O2 consumption.

4 Metrics linking global warming to marine hazards

The purpose of this section is to quantify the relationship between changes in global mean surface air temperature (SAT), the target variable of the Paris Agreement, with selected aggregated metrics for marine ecosystem stressors. In this way, we link marine hazards to the temperature target of the Paris Agreement and quantify avoided marine change per unit of avoided global warming. Specifically, we investigate the relationship of SAT with changes in the marine O2 inventory, ocean temperature and the metabolic index of Deutsch et al. (2015), the volume occupied by hypoxic water and in low-latitude export production (30 S–30 N) across the range of warming scenarios in our ensemble (Fig. 5). Distinct and often near-linear relationships emerge. Near-linearity allows us to characterize the benefits of avoided warming by single sensitivities, corresponding to the slopes of the relationships displayed in Fig. 5.

The relationships between SAT and marine hazard metrics critically depend on the time horizon considered (Fig. 5). We compare and contrast changes at the end of the 21st century, the typical assessment timescale of climate change, to peak and equilibrium changes at the millennial timescale. Peak and equilibrium changes are related to the corresponding equilibrium temperature response, while changes at the end of the 21st century are related to the transiently realized warming at the end of the 21st century. Larger magnitudes are simulated on millennial timescales compared to the near-term end of the 21st century. Assessment of ocean deoxygenation by the end of the 21st century, therefore, underestimates the full amplitude of change.

Transient (end of 21st century), peak (AD  3000) and equilibrium (AD  8000) oxygen changes exhibit distinct relationships with their corresponding warming (Fig. 5a). At the end of the 21st century, simulated oxygen decreases by 0.68 % C−1 of realized warming (median values). At peak oxygen decline, this sensitivity increases and oxygen decreases by 4.4 % C−1 of equilibrium temperature response. In other words, an avoided warming of 1 C, avoids a peak decline in marine O2 inventory of 4.4 %. The linear relationship breaks down for the equilibrium response. While 1.5 to 3.3 C warming targets lead to similar and higher oxygen levels, the 9.2 C warming target results in lower oxygen levels compared to preindustrial as discussed in the previous section. The relationships generally hold across the sampled parameter space.

The volume of low-oxygen waters is particularly sensitive to warming and parameter uncertainty (Fig. 5b). We illustrate the sensitivities at the example of the volume of waters with O2< 50 mmol m−3. At the end of the 21st century, there is a 1.7 % increase in this volume per C of realized warming. Peak increases scale with 63 % C−1 of equilibrium temperature response. Uncertainties in remineralization cause a spread in this response ranging from 36 to 76 % C−1 of equilibrium temperature response (90 % confidence range): the longer the remineralization length scale, the higher this sensitivity. Pre-existing low O2 waters expand and new low O2 waters may develop in near-bottom environments for higher forcing levels. While the lower temperature targets yield lower volumes of low-oxygen waters, the 9.2 C target yields higher low O2 volumes under new steady-state conditions. In brief, hypoxic waters expand over the next millennium across the scenario range and recovery towards modern conditions is slow and, in the case of high forcing, incomplete. Acknowledging millennial timescales, the hazard of expanding low O2 waters is much larger than when assessed on the near-term timescale.

The sensitivity of changes in low-latitude export production (30 S–30 N) is similar at the end of the 21st century and at time of peak decline. Changes scale with 2.2 % C−1. Export production in the low latitudes recovers for lower forcing, but remains reduced in the high forcing case.

Decreases in the metabolic index of the upper ocean scale linearly with forcing: 5.1 % C−1 of realized warming at the end of the 21st century and 6.2 % C−1 of equilibrium temperature response. Likewise, global mean oceanic temperatures increase by 0.099 C per degree of realized warming at the end of the 21st century and by 0.56 C per degree of warming at equilibrium. In conclusion, the compound hazards related to deoxygenation and warming, as indicated by the metabolic stress index, evolve over millennia and increase with increasing anthropogenic forcing and with time.

The magnitude and intensity, as well as the duration of oxygen-related transient hazards, and thus the severity of the hazards, increase with increasing temperature targets. The severity combines magnitude and duration of a hazard into one quantity. It may be defined as the time integral of a hazard. The severity of the hazard of expanding hypoxic waters, for example, corresponds to the area under the scenario curve shown in Fig. 2h (the area enclosed by the null line and the modeled evolution, here until the end of the simulation). Figure 2a, h, and g illustrate that the severity of the three hazards of decreasing mean oxygen concentration, expanding hypoxic waters and reduced export of particulate organic matter providing food for deep sea organisms increases strongly from low to high temperature targets.

5 Uncertainties in O2 projections

The pattern and magnitude of simulated global O2 changes are determined by the response of the overturning circulation. O2 loss due to less O2 solubility at higher temperatures also gradually decreases oceanic O2. Only a few multi-millennial simulations with GCMs currently exist. The response of the overturning circulation on long timescales differs among available model simulations (including EMICs and GCMs). Uncertainties in the equilibrium climate sensitivity additionally impact projections of O2 loss due to solubility. These uncertainties directly impact projections of oceanic oxygen.

Similar circulation dynamics as simulated here (Fig. 1e and f) were found by Rugenstein et al. (2016) based on EMIC simulations over 10 000 years with ECBILT-CLIO, which features a dynamic, quasi-geostrophic atmosphere. Schmittner et al. (2008) also found similar AMOC and Indo-Pacific MOC characteristics for their EMIC UVic 2.7, which includes an atmospheric energy balance model with fixed wind fields similar to the Bern3D model, over a 2000-year simulation. Yamamoto et al. (2015), on the other hand, found different overturning characteristics in a simulation with a state-of-the art Earth system model (MIROC 3.2 for a 4 ×CO2) over 2000 years. There AMOC slowed down with no recovery, while Antarctic Bottom Water formation decreased only slightly and gradually increased thereafter. Predictions of AMOC have received more attention so far, and AMOC slowdown and partial or full recovery emerges in other multi-millennial simulations (Li et al.2013; Weaver et al.2012; Zickfeld et al.2013). AMOC and Southern Ocean overturning in CMIP5 Earth system models were analyzed by Heuzé et al. (2015). They found AMOC and Southern Ocean overturning is positively correlated in most CMIP5 models by the end of the 21st century. Generally, preindustrial circulation states, magnitudes and timing of changes are highly model and scenario dependent such that the long-term evolution of meridional overturning is uncertain. As oxygen changes are dominated by circulation changes, this makes the oxygen prediction highly model and scenario dependent as well. The simulated timing and strength of peak O2 decrease in Bern3D is similar to what Schmittner et al. (2008) found (AD 3000, 30 % for SRES A2 high emission scenario/SAT  10 C in UVic 2.7). Other comparable simulations show earlier peaks and smaller magnitudes such as AD 2600, 16 % decrease for RCP8.5/ΔSAT  7 C in CLIMBER-3α in Mathesius et al. (2015), and after 800 model years, 10 % for 4 ×CO2/ΔSAT  8.5 C in MIROC 3.2 in Yamamoto et al. (2015).

Major physical limitations of our simulations concern prescribed winds and ice sheets. Future model studies may include sensitivity simulations with prescribed changes in the wind stress over the ocean (Tschumi et al.2008) and prescribed meltwater fluxes or apply Earth system models with interactive atmospheric dynamics and ice sheets. Our study, as is the case for most climate change simulations, does not include melting of continental ice sheets, which would tend to further (transiently) reduce circulation (Bakker et al.2016) and increase the equilibrium climate sensitivity.

Current generation GCMs, such as is the case for Bern3D, have difficulty simulating the current distribution of OMZs due to missing physical processes operating at small spatial scales, such as eddies and zonal jets (Bopp et al.2013; Cocco et al.2013) or missing biogeochemical characteristics. Large model–data and model–model discrepancies exist (Bopp et al.2013). Laufkötter et al. (2017) recently achieved improved representation of OMZs introducing temperature and oxygen dependence of the remineralization profile within a GCM (GFDL ESM2M). In our ensemble, the magnitude of peak increases in low O2 waters depend strongly on the rate of organic matter remineralization. Temperature-dependent feedback mechanisms, neglected here, may be addressed in future studies. Both particulate sinking speed and local remineralization rates, which control the remineralization profile, have been shown to be sensitive to temperature. While higher temperatures increase bacterial activity and therefore remineralization (Bendtsen et al.2014) they decrease viscosity and therefore increase sinking speed (Taucher et al.2014). The net effect on the remineralization profile is correspondingly uncertain. In addition, ecosystem structure influences the size and density of organic particles available for export (Armstrong et al.2001, 2009). Given these existing uncertainties and the coarse-resolution physical models, the projections of OMZs has to be viewed with caution. Despite simulated lower background concentrations of O2 in the subsurface ocean, the volumes of low O2 waters decrease for steady-state conditions in the model.

We neglect a number of biogeochemical feedback mechanisms that could alter biological productivity in the surface ocean and thus remineralization fluxes in the water column. Any mechanisms that would increase remineralization would tend to decrease oceanic oxygen, and mechanisms that decrease remineralization would increase the oceanic oxygen content. Future studies may address feedbacks from sediment interactions and imbalances from riverine input and burial (Niemeyer et al.2017; Roth et al.2014), temperature-dependent remineralization and variable stoichiometry. Further investigations may also address nitrogen cycle dynamics and assess the interplay of denitrification and N fixation and of external atmospheric and terrestrial nitrogen sources.

6 Implications and conclusion

In Bern3D, strong deoxygenation in all basins is projected to peak long after the end of the 21st century, and new steady-state conditions establish after AD 8000 in scenarios where radiative forcing is stabilized in the next century. The equilibration timescale of oceanic oxygen is therefore longer than the thermal equilibration timescale of both the atmosphere ( 1000 years) and the ocean ( 4000 years). Based on CMIP5 models, Sweetman et al. (2017) discuss the deep sea ecosystem implications of climate change by 2100. Deep sea ecosystems provide a range of services including habitat provision, nursery grounds, trophic support and refugia to biodiversity (Sweetman et al.2017). Biogeochemical changes such as deoxygenation, warming, acidification and less food availability will likely be accompanied by exploitation of mineral resources, overfishing and dumping of pollutants and microplastics. We project the largest biogeochemical changes to occur after 2100 and to aggravate over millennia. How these changes will affect deep sea ecosystems is poorly understood, but the adaptation to stress may be limited by their slow growth rates and long generation times (Sweetman et al.2017).

Figure 6Contrasting hazards of ecosystem impacts expressed in definitions of hypoxia, metabolic viability of the upper ocean and food availability. (a) Changes for a 1.5 C warmer world at present, at the end of the century and compared to peak changes. (b) Peak changes across 1.5, 1.9 and 3.3 C temperature targets. Lines correspond to the median response across each ensemble relative to 1870–1899.


Figure 6a contrasts near-term (AD 2100) and peak changes (relative to 1870–1899) in definitions of metabolically viable habitats in the upper ocean, hypoxia and food availability as projected by Bern3D for a 1.5 C warmer world. Export in low latitudes (30 S–30 N) as an indicator of food availability is reduced by maximally 4 % over the course of the simulation in this scenario. Median decreases in the metabolic index, representing viable habitat reductions of the upper ocean, amount to 11 % for a 1.5 warmer world. The volume of low-oxygen waters is particularly sensitive to anthropogenic warming and peak changes occur after the end of the 21st century. The volume of water with O2< 50 mmol m−3 changes by 6.6 % by the end of the century and by 14 % at its peak. Meeting the 1.5 C climate target of the Paris Agreement requires very fast and very stringent emission reductions (Millar et al.2017; Sanderson et al.2016; Steinacher et al.2013). Estimates by Steinacher et al. (2013) for a range of scenarios show that post-2017 allowable carbon emissions from fossil fuel need to be lower than 320 Gt C to meet the 1.5 C target with 66 % probability (320 Gt C is derived by adjusting the emission limit as displayed in Fig. 4 of Steinacher et al.2013, for year 2000 using fossil emissions published by the Global Carbon Project). This corresponds in the most optimistic scenario to only slightly more than 3 decades of current fossil fuel use. The nationally determined contributions, outlining emission mitigation actions by the parties of the Paris Agreement, need to be more ambitious and broader in scope to meet the 1.5 or 2 C targets (Joeri et al.2016). Such efforts would lead not only to lower warming compared to the current emission trajectory, but also have the benefit of reduced marine hazards as investigated here (Fig. 6b).

Higher temperature targets increase the hazard of ecosystem impacts as expressed in the investigated variables. In particular, definitions of peak hypoxia exhibit strong sensitivity to additional warming (Fig. 6b). Definitions of deoxygenation, marine food scarcity and marine aerobic habitat reduction are aggravated for the 2 C compared to the 1.5 C temperature target and investigated hazards are strongly amplified in a world where surface air temperature is stabilized at 3.3 C (Fig. 6b). Unbounded use of carbon emissions from existing fossil resources is projected not only to lead to a global warming on the order of 10 C (Randerson et al.2015; Zickfeld et al.2013), but also to a peak reduction in global mean oxygen inventory by almost a factor of 2 (Figs. 2a and 5a).

We find a near-linear relationship between impact-relevant marine hazards and global mean surface air temperature. This allows us to quantify avoided hazards per unit of avoided global warming. For example, emission mitigation measures would help to reduce peak O2 loss by 4.4 % C−1 of avoided equilibrium warming.

The Earth system response timescale to climate change spans several millennia such that anthropogenic perturbations to greenhouse gas concentrations commit the Earth system to long-term, irreversible climate change (Clark et al.2016). Our simulations show that the long-term fate of oceanic oxygen is characterized by an initial decline followed by a recovery phase. Peak decline and associated potential adverse ecosystem impacts are projected long after stabilization of radiative forcing in the atmosphere. This adds to the list of long-term Earth system commitments including warming, acidification and sea-level rise assessed elsewhere (Clark et al.2016; Eby et al.2009; Lord et al.2016; Pfister and Stocker2016; Ridgwell and Schmidt2010; Steinacher and Joos2016). Long-term, multi-millennial perspectives are thus required for a full account of climate-related ocean risks.

Data availability

Model output is available upon request to the corresponding author (

Spatial properties of a representative ensemble member reaching the 1.5 C warming target

In this appendix we document additional spatial properties of the representative ensemble member reaching the 1.5 C warming target. Figure A1 illustrates the meridional overturning stream function for pre-industrial (PI), year AD 2050, where the AMOC is at its lowest value, and for new steady-state conditions. Figure A2 illustrates the evolution of temperature, salinity and density anomalies across a transect from the Atlantic Ocean, through the Southern Ocean and into the Pacific at AD 2100 and AD 3150, when the O2 inventory values are at their lowest, and for new steady-state conditions. Figure A3 shows the anomalies in total O2 and the contributions from biology and solubility components for new steady conditions relative to PI. In addition, anomalies in ideal age are shown.

Figure A1Meridional overturning stream function in (a–c) the world ocean, (d–f) Atlantic Ocean and (g–i) Indo-Pacific for PI, year 2050, and for new steady-state conditions for a representative ensemble member reaching the 1.5 C warming target (columns). Circulation is clockwise along positive (red) contours and anticlockwise along the negative (blue) contours.


Figure A2Changes in temperature, salinity and density at AD 2100 (a–c), at AD 3150 (d–f) and for new steady-state conditions (g–i) compared to pre-industrial conditions.


Figure A3Changes in O2 and its components for new steady-state conditions relative to preindustrial steady state for a single representative ensemble member reaching the 1.5 C warming target. (a) Change in total O2, (b) change in O2 due to biology, (c) change in O2 due to solubility, (d) changes in ideal age.


Competing interests

The authors declare that they have no conflict of interest.

Special issue statement

This article is part of the special issue “The Earth system at a global warming of 1.5 C and 2.0 C”. It is not associated with a conference.


We thank Andreas Oschlies and Andreas Schmittner for their thorough reviews of our manuscript, which allowed us to better communicate our results. We also thank Patrik Pfister and Thomas Frölicher for fruitful discussions. This work was supported by the Swiss National Science Foundation (200020_172476).

Edited by: Axel Kleidon
Reviewed by: Andreas Oschlies and Andreas Schmittner


Allen, M. R., Frame, D. J., Huntingford, C., Jones, C. D., Lowe, J. A. Meinshausen, M., and Meinshausen, N.: Warming caused by cumulative carbon emissions towards the trillionth tonne, Nature, 458, 1163,, 2009. a

Armstrong, R. A., Lee, C., Hedges, J. I., Honjo, S., and Wakeham, S. G.: A new, mechanistic model for organic carbon fluxes in the ocean based on the quantitative association of POC with ballast minerals, Deep-Sea Res. Pt. II, 49, 219–236,, 2001. a

Armstrong, R. A., Peterson, M. L., Lee, C., and Wakeham, S. G.: Settling velocity spectra and the ballast ratio hypothesis, Deep-Sea Res. Pt. II, 56, 1470–1478,, 2009. a

Bakker, P., Schmittner, A., Lenaerts, J. T. M., Abe-Ouchi, A., Bi, D., van den Broeke, M. R., Chan, W.-L., Hu, A., Beadling, R. L., Marsland, S. J., Mernild, S. H., Saenko, O. A., Swingedouw, D., Sullivan, A., and Yin, J.: Fate of the Atlantic Meridional Overturning Circulation: Strong decline under continued warming and Greenland melting, Geophys. Res. Lett., 43, 12,252–12,260, 2016. a

Battaglia, G. and Joos, F.: Marine N2O Emissions From Nitrification and Denitrification Constrained by Modern Observations and Projected in Multimillennial Global Warming Simulations, Global Biogeochem. Cy., 32, 92–121,, 2018. a, b, c, d, e

Battaglia, G., Steinacher, M., and Joos, F.: A probabilistic assessment of calcium carbonate export and dissolution in the modern ocean, Biogeosciences, 13, 2823–2848,, 2016. a

Bendtsen, J., Hilligsøe, K. M., Hansen, J. L., and Richardson, K.: Analysis of remineralisation, lability, temperature sensitivity and structural composition of organic matter from the upper ocean, Prog. Oceanogr., 130, 125–145,, 2014. a

Bianchi, D., Dunne, J. P., Sarmiento, J. L., and Galbraith, E. D.: Data-based estimates of suboxia, denitrification, and N2O production in the ocean and their sensitivities to dissolved O2, Global Biogeochem. Cy., 26, GB2009,, 2012. a, b

Bopp, L., Le Quéré, C., Heimann, M., Manning, A. C., and Monfray, P.: Climate-induced oceanic oxygen fluxes: Implications for the contemporary carbon budget, Global Biogeochem. Cy., 16, 6-1–6-13, 2002. a

Bopp, L., Resplandy, L., Orr, J. C., Doney, S. C., Dunne, J. P., Gehlen, M., Halloran, P., Heinze, C., Ilyina, T., Séférian, R., Tjiputra, J., and Vichi, M.: Multiple stressors of ocean ecosystems in the 21st century: projections with CMIP5 models, Biogeosciences, 10, 6225–6245,, 2013. a, b, c, d, e, f, g

Bopp, L., Resplandy, L., Untersee, A., Le Mezo, P., and Kageyama, M.: Ocean (de)oxygenation from the Last Glacial Maximum to the twenty-first century: insights from Earth System models, Philos. T. Roy. Soc. Lond. A, 375, 2102,, 2017. a, b, c

Breitburg, D., Levin, L. A., Oschlies, A., Grégoire, M., Chavez, F. P., Conley, D. J., Garçon, V., Gilbert, D., Gutiérrez, D., Isensee, K., Jacinto, G. S., Limburg, K. E., Montes, I., Naqvi, S. W. A., Pitcher, G. C., Rabalais, N. N., Roman, M. R., Rose, K. A., Seibel, B. A., Telszewski, M., Yasuhara, M., and Zhang, J.: Declining oxygen in the global ocean and coastal waters, Science, 359, 6371,, 2018. a, b, c

Cabre, A., Marinov, I., Bernardello, R., and Bianchi, D.: Oxygen minimum zones in the tropical Pacific across CMIP5 models: Mean state differences and climate change trends, Biogeosciences, 12, 5429–5454,, 2015. a, b, c

Cheung, W. W. L., Reygondeau, G., and Frölicher, T. L.: Large benefits to marine fisheries of meeting the 1.5 C global warming target, Science, 354, 1591–1594,, 2016. a

Clark, P. U., Shakun, J. D., Marcott, S. A., Mix, A. C., Eby, M., Kulp, S., Levermann, A., Milne, G. A., Pfister, P. L., Santer, B. D., Schrag, D. P., Solomon, S., Stocker, T. F., Strauss, B. H., Weaver, A. J., Winkelmann, R., Archer, D., Bard, E., Goldner, A., Lambeck, K., Pierrehumbert, R. T., and Plattner, G.-K.: Consequences of twenty-first-century policy for multi-millennial climate and sea-level change, Nat. Clim. Change, 6, 360–369,, 2016. a, b, c

Cocco, V., Joos, F., Steinacher, M., Frölicher, T. L., Bopp, L., Dunne, J., Gehlen, M., Heinze, C., Orr, J., Oschlies, A., Schneider, B., Segschneider, J., and Tjiputra, J.: Oxygen and indicators of stress for marine life in multi-model global warming projections, Biogeosciences, 10, 1849–1868,, 2013. a, b, c, d, e

Deutsch, C., Ferrel, A., Seibel, B., Pörtner, H.-O., and Huey, R. B.: Climate change tightens a metabolic constraint on marine habitats, Science, 348, 1132–1135,, 2015. a, b, c, d, e, f, g

Diaz, R. J. and Rosenberg, R.: Spreading Dead Zones and Consequences for Marine Ecosystems, Science, 321, 926–929,, 2008. a

Doney, S. C., Lindsay, K., Fung, I., and John, J.: Natural Variability in a Stable, 1000-yr Global Coupled Climate–Carbon Cycle Simulation, J. Climate, 19, 3033–3054,, 2006. a

Eby, M., Zickfeld, K., Montenegro, A., Archer, D., Meissner, K. J., and Weaver, A. J.: Lifetime of Anthropogenic Climate Change: Millennial Time Scales of Potential CO2 and Surface Temperature Perturbations, J. Climate, 22, 2501–2511,, 2009. a

Frölicher, T. L., Joos, F., Plattner, G.-K., Steinacher, M., and Doney, S. C.: Natural variability and anthropogenic trends in oceanic oxygen in a coupled carbon cycle–climate model ensemble, Global Biogeochem. Cy., 23, gB1003,, 2009. a

Garcia, H. E., Locarnini, R. A., Boyer, T. P., Antonov, J. I., Baranova, O. K., Zweng, M. M., Reagan, J. R., and Johnson, D. R.: World Ocean Atlas 2013, in: Vol. 3: Dissolved Oxygen, Apparent Oxygen Utilization, and Oxygen Saturation, NOAA Atlas NESDIS 75, National oceanic and atmospheric administration (NOAA), (last access: June 2018), 2014. a, b

Gattuso, J.-P., Magnan, A., Billé, R., Cheung, W. W. L., Howes, E. L., Joos, F., Allemand, D., Bopp, L., Cooley, S. R., Eakin, C. M., Hoegh-Guldberg, O., Kelly, R. P., Pörtner, H.-O., Rogers, A. D., Baxter, J. M., Laffoley, D., Osborn, D., Rankovic, A., Rochette, J., Sumaila, U. R., Treyer, S., and Turley, C.: Contrasting futures for ocean and society from different anthropogenic CO2 emissions scenarios, Science, 349, 6243,, 2015. a, b, c, d

Gruber, N.: Warming up, turning sour, losing breath: ocean biogeochemistry under global change, Philos. T. Roy. Soc. A, 369, 1980–1996,, 2011. a

Heuzé, C., Heywood, K. J., Stevens, D. P., and Ridley, J. K.: Changes in Global Ocean Bottom Properties and Volume Transports in CMIP5 Models under Climate Change Scenarios, J. Climate, 28, 2917–2944,, 2015. a

Hofmann, M. and Schellnhuber, H.-J.: Oceanic acidification affects marine carbon pump and triggers extended marine oxygen holes, P. Natl. Acad. Sci. USA, 106, 3017–3022,, 2009. a

IPCC: Climate Change 2013: The Physical Science Basis, in: Contribution of Working Group I to the Fifth Assessment Report of the Intergovernmental Panel on Climate Change, Cambridge Univ. Press, Cambridge, UK and New York, NY, USA, p. 1535,, 2013. a, b

Ito, T., Follows, M. J., and Boyle, E. A.: Is AOU a good measure of respiration in the ocean?, Geophys. Res. Lett., 31, L17305,, 2004. a

Ito, T., Minobe, S., Long, M. C., and Deutsch, C.: Upper ocean O2 trends: 1958–2015, Geophys. Res. Lett., 44, 4214–4223, 2017. a

Joeri, R., d. Michel, E., Niklas, H., Taryn, F., Hanna, F., Harald, W., Roberto, S., Fu, S., Keywan, R., and Malte, M.: Paris Agreement climate proposals need a boost to keep warming well below 2 C, Nature, 534, 631–639, 2016. a

Kalnay, E., Kanamitsu, M., Kistler, R., Collins, W., Deaven, D., Gandin, L., Iredell, M., Saha, S., White, G., Woollen, J., Zhu, Y., Chelliah, M., Ebisuzaki, W., Higgins, W., Janowiak, J., Mo, K. C., Ropelewski, C., Wang, J., Leetmaa, A., Reynolds, R., Jenne, R., and Joseph, D.: The NCEP/NCAR 40-year reanalysis project, B. Am. Meteorol. Soc., 77, 437–471,<0437:TNYRP>2.0.CO;2, 1996. a

Keeling, R. F., Körtzinger, A., and Gruber, N.: Ocean Deoxygenation in a Warming World, Annu. Rev. Mar. Sci., 2, 199–229, 2010. a

Key, R. M., Kozyr, A., Sabine, C. L., Lee, K., Wanninkhof, R., Bullister, J. L., Feely, R. A., Millero, F. J., Mordy, C., and Peng, T.-H.: A global ocean carbon climatology: Results from Global Data Analysis Project (GLODAP), Global Biogeochem. Cy., 18, GB4031,, 2004. a

Laufkötter, C., John, J. G., Stock, C. A., and Dunne, J. P.: Temperature and oxygen dependence of the remineralization of organic matter, Global Biogeochem. Cy., 31, 1038–1050,, 2017. a

Li, C., von Storch, J.-S., and Marotzke, J.: Deep-ocean heat uptake and equilibrium climate response, Clim. Dynam., 40, 1071–1086,, 2013. a

Lord, N. S., Ridgwell, A., Thorne, M. C., and Lunt, D. J.: An impulse response function for the “long tail” of excess atmospheric CO2 in an Earth system model, Global Biogeochem. Cy., 30, 2–17,, 2016. a

Magnan, A. K., Colombier, M., Billé, R., Joos, F., Hoegh-Guldberg, O., Pörtner, H.-O., Waisman, H., Spencer, T., and Gattuso, J.-P.: Implications of the Paris agreement for the ocean, Nat. Clim. Change, 6, 732–735,, 2016. a

Martin, J. H., Knauer, G. A., Karl, D. M., and Broenkow, W.: VERTEX: Carbon cycling in the northeast Pacific, Deep-Sea Res., 34, 267–285,, 1987. a

Matear, R. J.: Climate change impacts on marine systems, Aust. Microbiol., 21, 17–20, 2000. a

Matear, R. J. and Hirst, A. C.: Long-term changes in dissolved oxygen concentrations in the ocean caused by protracted global warming, Global Biogeochem. Cy., 17, 1125,, 2003. a, b, c

Mathesius, S., Hofmann, M., Caldeira, K., and Schellnhuber Hans, J.: Long-term response of oceans to CO2 removal from the atmosphere, Nat. Clim. Change, 5, 1107–1113, 2015. a, b, c

McKay, M. D., Beckman, R. J., and Conover, W. J.: A comparison of three methods for selecting values of input variables in the analysis of output from a computer code, Technometrics, 21, 239–245, 1979. a

Meinshausen, M., Smith, S., Calvin, K., Daniel, J., Kainuma, M., Lamarque, J.-F., Matsumoto, K., Montzka, S., Raper, S., Riahi, K., Thomson, A., Velders, G., and van Vuuren, D. P.: The RCP greenhouse gas concentrations and their extensions from 1765 to 2300, Climatic Change, 109, 213–241,, 2011. a, b

Millar, R. J., Fuglestvedt, J. S., Friedlingstein, P., Rogelj, J., Grubb, M. J., Matthews, H. D., Skeie, R. B., Forster, P. M., Frame, D. J., and Allen, M. R.: Emission budgets and pathways consistent with limiting warming to 1.5 C, Nat. Geosci., 10, 741–748, 2017. a

Mislan, K. A. S., Deutsch, C. A., Brill, R. W., Dunne, J. P., and Sarmiento, J. L.: Projections of climate-driven changes in tuna vertical habitat based on species-specific differences in blood oxygen affinity, Global Change Biol., 23, 4019–4028,, 2017. a

Müller, S. A., Joos, F., Edwards, N. R., and Stocker, T. F.: Water mass distribution and ventilation time scales in a cost-efficient, three-dimensional ocean model, J. Climate, 19, 5479–5499,, 2006. a

Müller, S. A., Joos, F., Edwards, N. R., and Stocker, T. F.: Modeled natural and excess radiocarbon: Sensitivities to the gas exchange formulation and ocean transport strength, Global Biogeochem. Cy., 22, GB3011,, 2008. a

Najjar, R. G., Orr, J., Sabine, C. L., and Joos, F.: Biotic-HOWTO, Internal OCMIP Report, Tech. rep., LSCE/CEA Saclay, Gif-sur-Yvette, France, 1999. a

Niemeyer, D., Kemena, T. P., Meissner, K. J., and Oschlies, A.: A model study of warming-induced phosphorus–oxygen feedbacks in open-ocean oxygen minimum zones on millennial timescales, Earth Syst. Dynam., 8, 357–367,, 2017. a

Orr, J. C. and Epitalon, J.-M.: Improved routines to model the ocean carbonate system: mocsy 2.0, Geosci. Model Dev., 8, 485–499,, 2015. a

Orr, J. C. and Najjar, R. G.: Abiotic-HOWTO. Internal OCMIP Report, Tech. rep., LSCE/CEA Saclay, Gif-sur-Yvette, France, 1999. a

Oschlies, A., Duteil, O., Getzlaff, J., Koeve, W., Landolfi, A., and Schmidtko, S.: Patterns of deoxygenation: sensitivity to natural and anthropogenic drivers, Philos. T. Roy. Soc. Lond. A, 375, 2102,, 2017. a

Parekh, P., Joos, F., and Müller, S. A.: A modeling assessment of the interplay between aeolian iron fluxes and iron-binding ligands in controlling carbon dioxide fluctuations during Antarctic warm events, Paleoceanography, 23, PA4202,, 2008. a

Paris Agreement: United Nations Framework Convention on Climate Change, Adoption of The Paris Agreement,, last access: 10 April 2018. a, b, c

Pfister, P. L. and Stocker, T. F.: Earth system commitments due to delayed mitigation, Environ. Res. Lett., 11, 014010,, 2016. a

Plattner, G.-K., Joos, F., Stocker, T. F., and Marchal, O.: Feedback mechanisms and sensitivities of ocean carbon uptake under global warming, Tellus B, 53, 564–592,, 2001. a

Pörtner, H.-O.: Oxygen- and capacity-limitation of thermal tolerance: a matrix for integrating climate-related stressor effects in marine ecosystems, J. Exp. Biol., 213, 881–893,, 2010. a

Randerson, J. T., Lindsay, K., Munoz, E., Fu, W., Moore, J. K., Hoffman, F. M., Mahowald, N. M., and Doney, S. C.: Multicentury changes in ocean and land contributions to the climate-carbon feedback, Global Biogeochem. Cy., 29, 744–759,, 2015. a

Ridgwell, A. and Schmidt, D. N.: Past constraints on the vulnerability of marine calcifiers to massive carbon dioxide release, Nat. Geosci., 3, 196–200,, 2010. a, b, c

Ritz, S. P., Stocker, T. F., and Joos, F.: A coupled dynamical ocean-energy balance atmosphere model for paleoclimate studies, J. Climate, 24, 349–375,, 2011. a, b

Roth, R., Ritz, S. P., and Joos, F.: Burial-nutrient feedbacks amplify the sensitivity of atmospheric carbon dioxide to changes in organic matter remineralisation, Earth Syst. Dynam., 5, 321–343,, 2014. a, b

Rugenstein, M. A. A., Sedláček, J., and Knutti, R.: Nonlinearities in patterns of long-term ocean warming, Geophys. Res. Lett., 43, 3380–3388,, 2016. a

Sanderson, B. M., O'Neill, B. C., and Tebaldi, C.: What would it take to achieve the Paris temperature targets?, Geophys. Res. Lett., 43, 7133–7142,, 2016. a

Schmidtko, S., Stramma, L., and Visbeck, M.: Decline in global oceanic oxygen content during the past five decades, Nature, 542, 335–339,, 2017. a

Schmittner, A., Oschlies, A., Matthews, H. D., and Galbraith, E. D.: Future changes in climate, ocean circulation, ecosystems, and biogeochemical cycling simulated for a business-as-usual CO2 emission scenario until year 4000 AD, Global Biogeochem. Cy., 22, GB1013,, 2008. a, b, c, d

Shaffer, G., Olsen, S. M., and Pedersen, J. O. P.: Long-term ocean oxygen depletion in response to carbon dioxide emissions from fossil fuels, Nat. Geosci., 2, 105–109,, 2009. a, b

Steinacher, M. and Joos, F.: Transient Earth system responses to cumulative carbon dioxide emissions: linearities, uncertainties, and probabilities in an observation-constrained model ensemble, Biogeosciences, 13, 1071–1103,, 2016. a, b

Steinacher, M., Joos, F., Frölicher, T. L., Bopp, L., Cadule, P., Cocco, V., Doney, S. C., Gehlen, M., Lindsay, K., Moore, J. K., Schneider, B., and Segschneider, J.: Projected 21st century decrease in marine productivity: a multi-model analysis, Biogeosciences, 7, 979–1005,, 2010. a

Steinacher, M., Joos, F., and Stocker, T. F.: Allowable carbon emissions lowered by multiple climate targets, Nature, 499, 197–201,, 2013. a, b, c

Storch, D., Menzel, L., Frickenhaus, S., and Pörtner, H. O.: Climate sensitivity across marine domains of life: limits to evolutionary adaptation shape species interactions, Global Change Biol., 20, 3059–3067, 2014. a

Sweetman, A. K., Thurber, A., Smith, C., Levin, L., Mora, C., and Wei, C. L.: Major impacts of climate change on deep-sea benthic ecosystems, Elem. Sci. Anth., 5, 4,, 2017. a, b, c, d

Taucher, J., Bach, L. T., Riebesell, U., and Oschlies, A.: The viscosity effect on marine particle flux: A climate relevant feedback mechanism, Global Biogeochem. Cy., 28, 415–422, 2014. a

Tschumi, T., Joos, F., and Parekh, P.: How important are Southern Hemisphere wind changes for low glacial carbon dioxide? A model study, Paleoceanography, 23, PA4208,, 2008. a

Tschumi, T., Joos, F., Gehlen, M., and Heinze, C.: Deep ocean ventilation, carbon isotopes, marine sedimentation and the deglacial CO2 rise, Clim. Past, 7, 771–800,, 2011. a

UNFCCC: United Nations Framework Convention on Climate Change,, last access: 6 February 2018. a, b

Weaver, A. J., Sedláček, J., Eby, M., Alexander, K., Crespin, E., Fichefet, T., Philippon-Berthier, G., Joos, F., Kawamiya, M., Matsumoto, K., Steinacher, M., Tachiiri, K., Tokos, K., Yoshimori, M., and Zickfeld, K.: Stability of the Atlantic meridional overturning circulation: A model intercomparison, Geophys. Res. Lett., 39, L20709,, 2012. a

Weiss, R.: Carbon dioxide in water and seawater: The solubility of a non-ideal gas, Mar. Chem., 2, 203–215, 1974. a

Yamamoto, A., Abe-Ouchi, A., Shigemitsu, M., Oka, A., Takahashi, K., Ohgaito, R., and Yamanaka, Y.: Global deep ocean oxygenation by enhanced ventilation in the Southern Ocean under long-term global warming, Global Biogeochem. Cy., 29, 1801–1815,, 2015.  a, b, c, d, e

Zickfeld, K., Eby, M., Weaver, A. J., Alexander, K., Crespin, E., Edwards, N. R., Eliseev, A. V., Feulner, G., Fichefet, T., Forest, C. E., Friedlingstein, P., Goosse, H., Holden, P. B., Joos, F., Kawamiya, M., Kicklighter, D., Kienert, H., Matsumoto, K., Mokhov, I. I., Monier, E., Olsen, S. M., Pedersen, J. O. P., Perrette, M., Philippon-Berthier, G., Ridgwell, A., Schlosser, A., Deimling, T. S. V., Shaffer, G., Sokolov, A., Spahni, R., Steinacher, M., Tachiiri, K., Tokos, K. S., Yoshimori, M., Zeng, N., and Zhao, F.: Long-Term Climate Change Commitment and Reversibility: An EMIC Intercomparison, J. Climate, 26, 5782–5809,, 2013. a, b

Short summary
Human-caused, climate change hazards in the ocean continue to aggravate over a very long time. For business as usual, we project the ocean oxygen content to decrease by 40 % over the next thousand years. This would likely have severe consequences for marine life. Global warming and oxygen loss are linked, and meeting the warming target of the Paris Climate Agreement effectively limits related marine hazards. Developments over many thousands of years should be considered to assess marine risks.
Final-revised paper