# Inconsistent strategies to spin up models in CMIP5: implications for ocean biogeochemical model performance assessment

Roland Séférian, Marion Gehlen, Laurent Bopp, Laure Resplandy, James C. Orr, Olivier Marti, John P. Dunne, James R. Christian, Scott C. Doney, Tatiana Ilyina, et al.

<figure><img src="image_1.png" /></figure>

## To cite this version:

Roland Séférian, Marion Gehlen, Laurent Bopp, Laure Resplandy, James C. Orr, et al. Inconsistent strategies to spin up models in CMIP5: implications for ocean biogeochemical model performance assessment. Geoscientific Model Development, 2016, 9 (5), pp.1827-1851. 10.5194/gmd-9-1827-2016. hal-01339121

HAL Id: hal-01339121

https://hal.sorbonne-universite.fr/hal-01339121

Submitted on 29 Jun 2016

HAL is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published or not. The documents may come from teaching and research institutions in France or abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est destinée au dépôt et à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

<figure><img src="image_2.png" /></figure>
<figure><img src="image_3.png" /></figure>

# Inconsistent strategies to spin up models in CMIP5: implications for ocean biogeochemical model performance assessment

Roland Séférian<fnref n="1" />, Marion Gehlen<fnref n="2" />, Laurent Bopp<fnref n="2" />, Laure Resplandy<fnref n="3" />,<fnref n="2" />, James C. Orr<fnref n="2" />, Olivier Marti<fnref n="2" />, John P. Dunne<fnref n="4" />, James R. Christian<fnref n="5" />, Scott C. Doney<fnref n="6" />, Tatiana Ilyina<fnref n="7" />, Keith Lindsay<fnref n="8" />, Paul R. Halloran<fnref n="9" />, Christoph Heinze<fnref n="10" />,<fnref n="11" />, Joachim Segschneider<fnref n="12" />, Jerry Tjiputra<fnref n="11" />, Olivier Aumont<fnref n="13" />, and Anastasia Romanou<fnref n="14" />,<fnref n="15" />

Correspondence to: Roland Séférian (rseferian.cnrm@gmail.com)

Received: 31 August 2015 – Published in Geosci. Model Dev. Discuss.: 13 October 2015

Revised: 19 April 2016 – Accepted: 22 April 2016 – Published: 12 May 2016

**Abstract.** During the fifth phase of the Coupled Model Intercomparison Project (CMIP5) substantial efforts were made to systematically assess the skill of Earth system models. One goal was to check how realistically representative marine biogeochemical tracer distributions could be reproduced by models. In routine assessments model historical hindcasts were compared with available modern biogeochemical observations. However, these assessments considered neither how close modeled biogeochemical reservoirs were to equilibrium nor the sensitivity of model performance to initial conditions or to the spin-up protocols. Here, we explore how the large diversity in spin-up protocols used for marine biogeochemistry in CMIP5 Earth system models (ESMs) contributes to model-to-model differences in the simulated

fields. We take advantage of a 500-year spin-up simulation of IPSL-CM5A-LR to quantify the influence of the spin-up protocol on model ability to reproduce relevant data fields. Amplification of biases in selected biogeochemical fields (O₂, NO₃, Alk-DIC) is assessed as a function of spin-up duration. We demonstrate that a relationship between spin-up duration and assessment metrics emerges from our model results and holds when confronted with a larger ensemble of CMIP5 models. This shows that drift has implications for performance assessment in addition to possibly aliasing estimates of climate change impact. Our study suggests that differences in spin-up protocols could explain a substantial part of model disparities, constituting a source of model-to-model uncertainty. This requires more attention in future

<footnote n="1">¹CNRM, Centre National de Recherches Météorologiques, Météo-France/CNRS, 42 Avenue Gaspard Coriolis, 31057 Toulouse, France</footnote>

<footnote n="2">²LSCE/IPSL, Laboratoire des Sciences du Climat et de l'Environnement, Orme des Merisiers, CEA/Saclay 91198 Gif-sur-Yvette CEDEX, France</footnote>

<footnote n="3">³Scripps Institution of Oceanography, UCSD, La Jolla, CA, USA</footnote>

<footnote n="4">⁴Geophysical Fluid Dynamics Laboratory, NOAA, Princeton, NJ, USA</footnote>

<footnote n="5">⁵Fisheries and Oceans Canada and Canadian Centre for Climate Modelling and Analysis, Victoria, B.C., Canada</footnote>

<footnote n="6">⁶Marine Chemistry and Geochemistry Department, Woods Hole Oceanographic Institution, Woods Hole, MA, USA</footnote>

<footnote n="7">⁷Max Planck Institute for Meteorology, Bundesstraße 53, 20146 Hamburg, Germany</footnote>

<footnote n="8">⁸Climate and Global Dynamics Division, National Center for Atmospheric Research, Boulder, CO, USA</footnote>

<footnote n="9">⁹College of Life and Environmental Sciences, University of Exeter, Exeter, EX4 4RJ, UK</footnote>

<footnote n="10">¹⁰Geophysical Institute, University of Bergen, Bergen, Norway</footnote>

<footnote n="11">¹¹Uni Research Climate, Bjerknes Centre for Climate Research, Bergen, Norway</footnote>

<footnote n="12">¹²Department of Geosciences, University of Kiel, Kiel, Germany</footnote>

<footnote n="13">¹³Sorbonne Universités (UPMC, Univ Paris 06)-CNRS-IRD-MNHN, LOCEAN-IPSL Laboratory, 4 Place Jussieu, 75005 Paris, France</footnote>

<footnote n="14">¹⁴Dept. of Applied Math. and Phys., Columbia University, 2880 Broadway, New York, NY 10025, USA</footnote>

<footnote n="15">¹⁵NASA-Goddard Institute for Space Studies at Columbia University, New York, NY, USA</footnote>
model intercomparison exercises in order to provide quan-
titatively more correct ESM results on marine biogeochem-
istry and carbon cycle feedbacks.

# 1 Introduction

## 1.1 Context

Earth system models (ESMs) are recognized as the current state-of-the-art global coupled models used for climate research (e.g., Hajima et al., 2014; IPCC, 2013). They expand the numerical representation of the climate system used during the 4th IPCC assessment report (AR4) that was limited to coupled physical general circulation models, to the inclusion of biogeochemical and biophysical interactions between the physical climate system and the biosphere. The ESMs that contributed to CMIP5 substantially differed from each other in terms of their simulations of physical and biogeochemical components of the Earth system. These differences in design translate into a significant variability between the skill with which the different models reproduce the observed biogeochemistry and carbon cycle, which in turn may impact projected climate change responses (IPCC, 2013).

In the typical objective evaluation and intercomparison of these models, a suite of standardized statistical metrics (e.g., correlation, root-mean-squared errors) are applied to quantify differences between modeled and observed variables (e.g., Doney et al., 2009; Rose et al., 2009; Stow et al., 2009; Romanou et al., 2013, 2014). With the goal of constraining future projections, statistical metrics are often used for model ranking (e.g., Anav et al., 2013), weighting of model projections (e.g., Steinacher et al., 2010) or selection of the most skillful models across a wider ensemble (e.g., Cox et al., 2013; Massonnet et al., 2012; Wenzel et al., 2014). Most of these approaches can be considered as “blind” given that they are routinely applied without considering models’ specific characteristics and treat models a priori as equivalently independent of observations. However, since these models are typically initialized from observations, the spin-up procedure (e.g. the length of time for which the model has been run since initialization, and therefore the degree to which it has approached its own equilibrium) has the potential to exert a significant control over the statistical metrics calculated for each model, using those observations.

## 1.2 Initialization of biogeochemical fields and spin-up protocols in CMIP5

Ocean initialization protocols aim at obtaining stable and equilibrated distributions of model state variables, such as temperature or concentrations of dissolved tracers. Most commonly used initialization protocols consist of initializing both physical and biogeochemical variables from either climatologies (derived from the observed fields or previous

model simulations) or spatially constant values before running the model to equilibrium. In theory, equilibrium corresponds to steady-state and, hence, temporal derivatives of tracer fields close to zero. The time needed to equilibrate tracer distributions or, in other words, the integration time needed by the model to converge towards its own attractor (which is different from the true state of the climate system) varies greatly between components of the climate system. It spans from several weeks for the atmosphere (e.g., Phillips et al., 2004) to several centuries for ocean and sea ice components (e.g., Stouffer et al., 2004). The equilibration of ocean biogeochemical tracers across the entire water column amounts to several thousands of years (e.g., Heinze et al., 1999; Wunsch and Heimbach, 2008) and depends on the state of background ocean circulation as well as the turbulent mixing and eddy stirring parameterizations (e.g., Auumont et al., 1998; Bryan, 1984; Gnanadesikan, 2004; Marinov et al., 2008). The equilibration time can be different in a coupled model configuration (i.e., ocean-atmosphere general circulation models or ESMs) compared to stand-alone climate components due to leaks in the energy budget (Hobbs et al., 2016). In practice, these simulations, called “spin-ups”, often span only several hundreds of years, at the end of which a quasi-equilibrium state is assumed for the interior ocean tracers.

The present degree of complexity and spatial as well as temporal resolution of marine biogeochemical ESM components (as well as other physical and chemical components), however, often precludes a spin-up to reach adequate equilibration of biogeochemical tracers. This is a consequence of the large number of state variables present in most of the current generation of biogeochemical models (e.g., for each tracer a separate advection equation has to be solved via a numerical CPU time demanding algorithm), more complex process descriptions (e.g., including more plankton functional types than before), and spatial as well as temporal resolution. This number of state variables has continuously increased from simple biogeochemical models (e.g., HAMOCC3, Maier-Reimer and Hasselmann, 1987) to marine biodiversity models (e.g., Follows et al., 2007). Current generation biogeochemical models embedded in CMIP5 ESMs contain roughly 2 to 4 times more state variables than physical models (e.g., atmosphere, ocean, sea-ice), which makes their equilibration computationally costly and difficult. The initialization of biogeochemical state variables is further complicated by the scarcity of biogeochemical observations as compared to observations of physical variables (e.g., temperature, salinity). So far, three-dimensional observation-based climatologies exist for macro-nutrients, oxygen, dissolved carbon and alkalinity. For other tracers such as dissolved iron, dissolved organic carbon and biomass of the various plankton functional types data are still sparse in space and time in spite of considerable efforts such as the GEOTRACES program for trace elements, or MARE-DAT for biomasses of plankton functional types. The lat-
ter set of variables is initialized either with constant values (e.g., global average estimates) or with output from a previous model run. An additional difficulty stems from the use of modern climatologies to initialize the ocean state, implicitly assuming a long-term steady state, which does not necessarily represent the preindustrial state of the ocean. These climatologies incorporate the ongoing anthropogenic perturbation of marine biogeochemical fields, be it the uptake of anthropogenic CO₂ or the excess of nutrients inputs and pollutants (e.g., Doney, 2010). Although methods exist to remove the anthropogenic perturbation from some observed ocean carbon tracer fields, their use is still debated since they lead to non-unique results (e.g., Tanhua et al., 2007; Yool et al., 2010).

The equilibration of marine biogeochemical tracer distributions is driven not only by the ocean circulation, but also by numerous internal biogeochemical processes acting at various timescales. For example, while the transport and degradation of sinking organic matter spans days to perhaps several months, the associated impact on deep water chemistry accumulates over several decades to centuries as zones of differential remineralization are mixed across water masses and follows the ocean circulation (Wunsch and Heimbach, 2008). For models including interactive sediment modules, the sediment equilibration takes even longer ($O(10^4)$ years; e.g., Archer et al., 2009, and Heinze et al., 1999). As a consequence of the interplay between ocean circulation and biogeochemical processes, biogeochemical models require long spin-up times to equilibrate (e.g., Khatiwala et al., 2005; Wunsch and Heimbach, 2008). Modeling studies of paleo-oceanographic passive tracers such as $δ^{18}O$ or $Δ^{14}C$ (Duplessy et al., 1991), or global ocean passive tracers (Wunsch and Heimbach, 2008), as well as more recently available modern global-scale data compilations (e.g., Key et al., 2004; Sarmiento and Gruber, 2006) and GEOTRACES Intermediate Data product 2014 (Mawji et al., 2015) provide an estimate of the time required for the ocean biogeochemical reservoir to equilibrate with the climate systems (excluding continental weathering and reaction with marine sediments). For the deep water masses, this time is about 1500 years in the Atlantic Ocean and reaches up to 10 000 years in the North Pacific Ocean (Wunsch and Heimbach, 2008).

In a context of model-to-model intercomparison, this time range contributes to the model uncertainty. Lessons from the previous Ocean Carbon Model Intercomparison Project phase 2 (OCMIP-2) exercise have demonstrated that some models required ~10 000 years to reach a state where the global sea-air carbon flux is about 0.01 Pg C y⁻¹.

While it is recognized that long timescale processes influence the length of spin-up to equilibrium, the spin-up duration is usually defined ad hoc based on external constraints or internal biogeochemical criteria. The computational cost is commonly invoked as external constraint to shorten and limit the spin-up duration. It is directly related to model com-

plexity (e.g., Tjiputra et al., 2013; Vichi et al., 2011; Yool et al., 2013) and spatial resolution (Ito et al., 2010). The internal biogeochemical criteria applied to derive the duration of the spin-up simulations are generally defined by (i) reaching a steady-state, quasi equilibrium of the long-term global CO₂ fluxes between the ocean and the atmosphere (e.g., Dunne et al., 2013; Ilyina et al., 2013; Lindsay et al., 2014; Romanou et al., 2013; Séférian et al., 2013), (ii) determining the amount of carbon stored in the ocean at preindustrial state (e.g., Dunne et al., 2013; Vichi et al., 2011) or (iii) representing relevant biogeochemical tracer patterns (e.g., oxygen minimum zone in Ito and Deutsch, 2013).

Despite its importance, only limited information on spin-up procedures is available through the CMIP5 metadata portal (http://metaforclimate.eu/trac). Information on spin-up protocols and model initialization is usually not taken into account in model intercomparison studies (e.g., Andrews et al., 2013; Bopp et al., 2013; Cocco et al., 2013; Frölicher et al., 2014; Gehlen et al., 2014; Keller et al., 2014; Resplandy et al., 2013, 2015; Rodgers et al., 2015; Séférian et al., 2014). This information, if available, can only be found separately in the reference papers of individual models (e.g., Adachi et al., 2013; Arora et al., 2011; Collins et al., 2011; Dunne et al., 2013; Ilyina et al., 2013; Lindsay et al., 2014; Romanou et al., 2013; Séférian et al., 2013, 2016; Tjiputra et al., 2013; Vichi et al., 2011; Volodin et al., 2010; Watanabe et al., 2011; Wu et al., 2013). The duration of spin-up simulations of CMIP5 ocean biogeochemical components spans from one hundred years (e.g., CMCC-CESM) to several thousand years (e.g., MPI-ESM-LR, MPI-ESM-MR) (Fig. 1 and Table 1). Model initialization and spin-up procedures are equally variable across the model ensemble (Fig. 1 and Table 1). Four different sources of initialization and four different procedures of model equilibration emerge from the 24 ESMs reviewed for this study.

Biogeochemical state variables were mostly initialized from observations, although from various releases of the same World Ocean Atlas global climatology (WOA1994, WOA2001, WOA2006, WOA2010). A small subset of ESMs relied either on a mix between previous model output and observations or solely on model output from a previous simulation for initialization. Similarly, spin-up procedures fall into two categories. The first one may be called “sequential”: it consists in decomposing the spin-up integration into one long offline simulation (~200–10 000 years) and one shorter online simulation (~100–1000 years). During the offline simulation, the biogeochemical model is forced by dynamical fields from the climate model or from reanalysis (CanESM2, MRI-ESM, Fig. 1 and Table 1). Some modeling groups have adopted a “direct” strategy, which consists in running solely one online or coupled spin-up simulation (e.g., CNRM-ESM1, GFDL-ESM2M, GFDL-ESM2G, GISS-E2-H-CC, GISS-E2-R-CC, NorESM1-ME). Finally, a spin-up “acceleration” procedure is used by CMCC-CESM. This technique consists of enhancing the ocean carbon out-
<figure><img src="image_4.png" /><figcaption>**Figure 1.** Spin-up protocols of CMIP5 Earth system models. Color shading represents strategies of the various modeling groups. *Online* and *Offline* steps refer to runs performed with coupled climate model and with stand-alone ocean biogeochemistry model, respectively. Sources of initial conditions for biogeochemical components of CMIP5 Earth system models are indicated as hatching below the barplot.</figcaption></figure>

gassing to remove anthropogenic carbon from the ocean, a legacy from initialization with modern data (Global Data Analysis Project or GLODAP following Key et al., 2004). None of these spin-up procedures, durations and sources of initialization can be considered as “standard”; each of them is unique and subjectively employed by one modeling group.

Objective arguments and hypotheses justifying the choice of one method of spin-up rather than the others have been the focus of previous studies (e.g., Dunne et al., 2013; Heinze and Ilyina, 2015; Tjiputra et al., 2013). Similarly, individual modeling groups have discussed the impacts of their particular spin-up procedure on model performance individually (e.g., Dunne et al., 2013; Lindsay et al., 2014; Séférian et al., 2013; Vichi et al., 2011). However, no study has addressed the potential for the large diversity of spin-up procedures found across the CMIP5 ensemble to translate into model-to-model differences in terms of comparative model performance assessments or model evaluations in terms of future projections.

## 1.3 Objectives of this study

This study assesses the role of the spin-up protocol in controlling the “final” representation of biogeochemical fields, and subsequent model skill assessment, providing a complementary analysis from the studies of Sen Gupta et al. (2012, 2013). It relies on a 500-year long spin-up simulation from a state-of-the-art Earth system model, IPSL-CM5A-LR to investigate the impacts of spin-up strategy on selected biogeochemical tracers and residual model drift across the various ESMs of the CMIP5 ensemble. We demonstrate that the duration of the spin-up has implications for the determination of robust and meaningful skill-score metrics that should improve future intercomparison studies such as CMIP6 (Meehl et al., 2014).

Section 2 describes the model, the observations, the model experiments, as well as the methods used for assessing the impacts of spin-up protocols on the representation of biogeochemical fields in IPSL-CM5A-LR, as well as across the ensemble of CMIP5 ESMs. Section 3 presents the analysis
<table><caption>**Table 1.** Summary of spin-up strategy, sources of initial conditions, offline/online durations and references used to equilibrate ocean biogeochemistry in CMIP5 ESMs. The so-called direct and sequential strategies inform whether the spin-up of the ocean biogeochemical model is run directly in online/coupled mode or first in offline (ocean biogeochemistry only) and then in online/coupled mode. DIC* refers to the observation-derived estimates of preindustrial dissolved inorganic carbon concentration using the ΔC* method. w/acc. and forced w/obs. indicate the strategy using “acceleration” and observed atmospheric forcings during the spin-up, respectively.</caption><thead><tr><th>Models</th><th>Spin-up procedure</th><th>Initial conditions</th><th>Offline time</th><th>Online time</th><th>Total spin-up duration</th><th>References</th></tr></thead><tbody><tr><td>BCC-CSM1-1</td><td>Sequential</td><td>WOA2001, GLODAP</td><td>200</td><td>100</td><td>300</td><td>Wu et al. (2013)</td></tr><tr><td>BCC-CSM1-1-m</td><td>Sequential</td><td>WOA2001, GLODAP</td><td>200</td><td>100</td><td>300</td><td>Wu et al. (2013)</td></tr><tr><td>CanESM2</td><td>Sequential (forced w/obs.)</td><td>OCMIP profiles, CanESM1</td><td>6000</td><td>600</td><td>6600</td><td>Arora et al. (2011)</td></tr><tr><td>CESM1-BGC</td><td>Direct</td><td>CCSM4</td><td>0</td><td>1000</td><td>1000</td><td>Lindsay et al. (2014)</td></tr><tr><td>CMCC-CESM</td><td>Sequential (w/acc.)</td><td>WOA2001, GLODAP</td><td>100</td><td>100</td><td>200</td><td>Vichi et al. (2011)</td></tr><tr><td>CNRM-CM5</td><td>Sequential</td><td>WOA1994, GLODAP, IPSL</td><td>3000</td><td>100</td><td>3100</td><td>Séférían et al. (2013)</td></tr><tr><td>CNRM-CM5-2</td><td>Sequential</td><td>WOA1994, GLODAP, CNRM</td><td>3000</td><td>100</td><td>3100</td><td>Schwinger et al. (2014)</td></tr><tr><td>CNRM-ESM1</td><td>Sequential</td><td>CNRM-CM5</td><td>0</td><td>1300</td><td>1300</td><td>Séférían et al. (2016)</td></tr><tr><td>GFDL-ESM2G</td><td>Direct</td><td>WOA2005, GLODAP</td><td>0</td><td>1000</td><td>1000</td><td>Dunne et al. (2013)</td></tr><tr><td>GFDL-ESM2M</td><td>Direct</td><td>WOA2005, GLODAP</td><td>0</td><td>1000</td><td>1000</td><td>Dunne et al. (2013)</td></tr><tr><td>GISS-E2-H-CC</td><td>Direct</td><td>WOA2005, GLODAP DIC*</td><td>0</td><td>3300</td><td>3300</td><td>Romanou et al. (2013)</td></tr><tr><td>GISS-E2-R-CC</td><td>Direct</td><td>WOA2005, GLODAP DIC*</td><td>0</td><td>3300</td><td>3300</td><td>Romanou et al. (2013)</td></tr><tr><td>HadGEM2-CC</td><td>Sequential</td><td>HadCM3LC, WOA2011</td><td>400</td><td>100</td><td>500</td><td>Collins et al. (2011), Wassmann et al. (2010)</td></tr><tr><td>HadGEM2-ES</td><td>Sequential</td><td>HadCM3LC, WOA2010</td><td>400</td><td>100</td><td>500</td><td>Collins et al. (2011)</td></tr><tr><td>INMCM4</td><td>Sequential</td><td>Uniform DIC</td><td>3000</td><td>200</td><td>3200</td><td>Volodin et al. (2010)</td></tr><tr><td>IPSL-CM5A-LR</td><td>Sequential</td><td>WOA1994, GLODAP, IPSL</td><td>3000</td><td>600</td><td>3600</td><td>Séférían et al. (2013)</td></tr><tr><td>IPSL-CM5A-MR</td><td>Sequential</td><td>WOA1994, GLODAP, IPSL</td><td>3000</td><td>300</td><td>3300</td><td>Dufresne et al. (2013)</td></tr><tr><td>IPSL-CM5B-LR</td><td>Sequential</td><td>IPSL-CM5A-LR</td><td>0</td><td>300</td><td>300</td><td>Dufresne et al. (2013)</td></tr><tr><td>MIROC-ESM</td><td>Sequential</td><td>GLODAP/constant values</td><td>1245</td><td>480</td><td>1725</td><td>Watanabe et al. (2011)</td></tr><tr><td>MIROC-ESM-CHEM</td><td>Sequential</td><td>GLODAP/constant values</td><td>1245</td><td>484</td><td>1729</td><td>Watanabe et al. (2011)</td></tr><tr><td>MPI-ESM-LR</td><td>Sequential</td><td>HAMOCC/constant values</td><td>10 000</td><td>1900</td><td>11 900</td><td>Ilyina et al. (2013)</td></tr><tr><td>MPI-ESM-MR</td><td>Sequential</td><td>HAMOCC/constant values</td><td>10 000</td><td>1500</td><td>11 500</td><td>Ilyina et al. (2013)</td></tr><tr><td>MRI-ESM1</td><td>Sequential (forced w/obs.)</td><td>GLODAP</td><td>550</td><td>395</td><td>945</td><td>Adachi et al. (2013)</td></tr><tr><td>NorESM</td><td>Direct</td><td>WOA2010, GLODAP</td><td>0</td><td>900</td><td>900</td><td>Tjiputra et al. (2013)</td></tr></tbody></table>

developed for the assessment of the impact of spin-up duration on the representation of biogeochemical structures. Implications and recommendations are discussed in Sects. 4 and 5, respectively.

## 2 Methods

### 2.1 Model simulations

This study exploits in particular results from one simulation performed with IPSL-CM5A-LR (Dufresne et al., 2013), considered here to be representative of the likely behavior of other CMIP5 Earth system models. Like other current generation of ESMs, IPSL-CM5A-LR combines the major components of the climate system (Chap. 9, Table 9.1 of IPCC, 2013). The atmosphere is represented by the atmospheric general circulation model LMDZ (Hourdin et al., 2006) with a horizontal resolution of 3.75° × 1.87° and 39 levels. The land surface is simulated with ORCHIDEE (Krinner et al.,

2005). The oceanic component is NEMOv3.2 in its ORCA2 global configuration (Madec, 2008). It has a horizontal resolution of about 2° with enhanced resolution at the equator (0.5°) and 31 levels. NEMOv3.2 includes the sea-ice model LIM2 (Fichefet and Maqueda, 1997), and the marine biogeochemistry model PISCES (Aumont and Bopp, 2006). PISCES simulates the biogeochemical cycles of oxygen, carbon and the main nutrients with 24 state variables. The model simulates dissolved inorganic carbon and total alkalinity (carbonate alkalinity + borate + water) and the distributions of macronutrients (nitrate and ammonium, phosphate, and silicate) and the micronutrient iron. PISCES represents two sizes of phytoplankton (i.e., nanophytoplankton and diatoms) and two zooplankton size classes: microzooplankton and mesozooplankton. PISCES simulates semi-labile dissolved organic matter, and small and large sinking particles with different sinking speeds (3 and 50 to 200 m day⁻¹, respectively). While fixed elemental stoichiometric C:N:P:ΔO₂ ratios after Takahashi et al. (1985) are imposed for
these three compartments the internal concentrations of iron, silica and calcite are simulated prognostically. The carbon system is represented by dissolved inorganic carbon, alkalinity and calcite. Calcite is prognostically simulated following Maier-Reimer (1993) and Moore et al. (2002). Alkalinity in the model system includes the contribution of carbonate, bicarbonate, borate, protons, and hydroxide ions. Oxygen is prognostically simulated. The model distinguishes between oxic and suboxic remineralization pathways, the former relying on oxygen as electron acceptor, the latter on nitrate. For carbon and oxygen pools, air-sea exchange follows the Wanninkhof (1992) formulation.

The model's boundary conditions account for nutrient supplies from three different sources: atmospheric dust deposition for iron, phosphorus and silica (Jickells and Spokes, 2001; Moore et al., 2004; Tegen and Fung, 1995), rivers for nutrients, alkalinity and carbon (Ludwig et al., 1996) and sediment mobilization for sedimentary iron (de Baar and de Jong, 2001; Johnson et al., 1999). To ensure conservation of nitrogen in the ocean, annual total nitrogen fixation is adjusted to balance losses from denitrification. For the other macronutrients, alkalinity and organic carbon, the conservation is ensured by tuning the sedimental burial to the total external input from rivers and dust. In PISCES, an adequate treatment of external boundary conditions has been demonstrated to be essential for the accurate simulation of nutrient distributions (Aumont and Bopp, 2006; Aumont et al., 2003). Riverine carbon inputs induce a natural outgassing of carbon of $0.6 \text{ Pg C y}^{-1}$ that has been shown to be essential to model the inter-hemispheric gradient of atmospheric $CO_2$ under a preindustrial state (Aumont et al., 2001).

The core simulation used within this study is a 500-year long coupled preindustrial run. It uses the same atmospheric, land surface and ocean configurations as IPSL-CM5A-LR (Dufresne et al., 2013) for which the marine biogeochemistry has been extensively evaluated (see, e.g., Séférían et al., 2013, for modern-state evaluation). The only difference between the “standard” preindustrial simulation contributed to CMIP5 and the present one is the initial conditions. While the CMIP5 preindustrial simulation starts from an ocean circulation after several thousand years of online physical adjustment, the present simulation starts from an ocean at rest using the January temperature and salinity fields from the World Ocean Atlas (Levitus and Boyer, 1994). Biogeochemical state variables were initialized from data compilations or climatologies as explained in the following section. Atmospheric $CO_2$ and other greenhouse gases, as well as natural aerosols, were set to their 1850 preindustrial values. The simulation is extensively described in terms of ocean physics by Mignot et al. (2013). Mignot and coworkers show that the strength of the Atlantic meridional overturning circulation and the Antarctic circumpolar current as well as the upper 300 m ocean heat content stabilize after 250 years of simulation.

Although the spin-up protocol used to conduct this 500-year long simulation is not readily comparable to the one used to produce the initial conditions for the CMIP5 preindustrial simulation, its duration is greater than the median length of online adjustment computed from the multiple spin-up protocols applied during CMIP5 (~395 years, Fig. 1 and Table 1). Besides, the methodology of initializing biogeochemical state variables from data fields is not broadly employed by the various modeling groups that have contributed to CMIP5. Despite the above-mentioned methodological shortcuts, we take this 500-year long preindustrial simulation as a representative example of a spin-up protocol given the diversity of approaches used by CMIP5 models.

## 2.2 Observations for initialization and evaluation

Two streams of data sets were used in this study. The first stream combines data from the World Ocean Atlas 1994 (WOA94, Levitus and Boyer, 1994, and Levitus et al., 1993) for the initialization of three-dimensional fields of temperature and salinity, dissolved nitrate, silicate, phosphate and oxygen, and data from GLODAP (Key et al., 2004) for preindustrial dissolved inorganic carbon and total alkalinity. This stream of data was chosen purposely in our experimental setup to be slightly different than the second stream of data, World Ocean Atlas 2013 (WOA2013, Levitus et al., 2013), the evaluation data set.

A second stream of data was used to compare modeled biogeochemical fields. It includes up-to-date observed climatologies of nitrate and oxygen from the WOA2013. This database is based on samples collected since 1965, and including data more recently collected than that made use of in WOA94. For the concentrations of preindustrial dissolved inorganic carbon and total alkalinity, we still use GLODAP. The second stream of data was selected to be as close as possible to the “standard” evaluation procedure of skill-assessment protocols found in CMIP5 model reference papers (Adachi et al., 2013; Arora et al., 2011; Collins et al., 2011; Dunne et al., 2013; Ilyina et al., 2013; Lindsay et al., 2014; Romanou et al., 2013; Séférían et al., 2013, 2016; Tjiputra et al., 2013; Vichi et al., 2011; Volodin et al., 2010; Watanabe et al., 2011; Wu et al., 2013). Differences between these two streams of data are minor and are further detailed below.

## 2.3 Approach and statistical analysis

To quantify the impacts of a large diversity of spin-up procedures on the representation of biogeochemical fields in CMIP5, we employ a three-fold approach.

1. The 500-year long spin-up simulation described in Sect. 2.1 is used to determine the influence of the spin-up procedure on the representation of biogeochemical fields in IPSL-CM5A-LR.
2. In the next step, relationships between biases in modeled fields, model-data mismatches and the duration of the spin-up simulation are identified across the CMIP5 ensemble. For this step, drifts in biogeochemical fields are determined from the first century of the preindustrial simulation (referred to as *piControl*) of each CMIP5 ESM.

3. Finally, the ensembles of industrial-revolution to present-day simulation (referred to as *historical*) from each available CMIP5 ESM are used to estimate the impact of these drifts in biogeochemical fields on the ability of models to replicate modern observations. For a given model, we use the ensemble average of the available “historical” members if several realizations are available.

For this purpose, several statistical skill score metrics are computed following Rose et al. (2009) and Stow et al. (2009) from model fields interpolated on a regular 1° grid and to fixed depth levels. The skill score metrics are (1) the globally averaged concentrations for overall drift; (2) the error or bias between modeled and observed fields at each grid cell; (3) spatial correlation between model and observations to assess mismatches between modeled and observed large-scale structures; (4) the root-mean squared error (RMSE) to assess the total cumulative errors between modeled and observed fields. These statistical metrics are computed at different depth levels, but for clarity we focus on surface, 150 m (thermocline) and 2000 m (deep) levels. These statistical metrics were chosen among those described in the literature, because they proved to yield the most indicative scores for tracking model errors or improvement along the various intercomparison exercises (IPCC, 2013).

The drift is determined for either concentrations in simulated biogeochemical fields or for skill score metrics (e.g., RMSE) using a linear regression fit over a time window of 100 years. This time window of 100 years was chosen as a trade-off between a longer time window (> 200 years) that smoothes the drift signal and a shorter time window (< 100 years) that introduces fluctuations due to internal variability and hence impacting the quality of the fit (see the assessment performed with the millennial-long CMIP5 *piControl* simulation of IPSL-CM5A-LR in Fig. S1 in the Supplement).

The drift is assumed to decrease exponentially during the spin-up simulation and is described by a simple drift model:

$$drift(t) = drift(t = 0) \times \exp(-\frac{1}{\tau}t) \quad (1)$$

where $τ$ is the relaxation time of the respective field at a given depth level. It corresponds to the time required to nullify the drift.

Our analyses focus on the global distribution of nitrate ($NO_3$), dissolved oxygen ($O_2$) and the difference between total alkalinity and dissolved inorganic carbon (Alk-DIC).

The latter serves as an approximation of carbonate ion concentration following Zeebe and Wolf-Gladrow (2001). We use this approximation of the carbonate ion concentration rather than its concentration, $[CO_3^{2-}]$, since the latter was poorly assessed in CMIP5 reference papers and was not provided by a majority of ESMs. These three biogeochemical tracers were chosen because (1) most current biogeochemical models simulate Alk, DIC, $NO_3$ and $O_2$ prognostically and (2) they are frequently used in state-of-the-art model performance assessment (e.g., Anav et al., 2013; Bopp et al., 2013; Doney et al., 2009; Friedrichs et al., 2009, 2007; Stow et al., 2009), and (3) DIC and Alk are both used as “master tracers” for the carbonate system in the ocean biogeochemistry models (while $[CO_3^{2-}]$, e.g., is not explicitly transported as a tracer with the velocity fields but diagnosed from temperature, salinity, DIC, Alk, $[H^+]$, and $pCO_2$ when needed). Modeled distributions of $NO_3$, $O_2$ and Alk-DIC reflect the representation of biogeochemical processes related to the biological pump ($CO_2$, $NO_3$, $O_2$), the air-sea gas exchange and ocean ventilation ($CO_2$ and $O_2$), as well as carbonate chemistry (Alk-DIC). These biogeochemical processes are of particular relevance for investigating the impact of climate change on marine productivity (e.g., Henson et al., 2010), ocean deoxygenation (e.g., Gruber, 2011; Keeling et al., 2009) and the ocean carbon sink, processes for which future projections with the current generation of ESMs yield large inter-model spreads (e.g., Friedlingstein et al., 2013; Resplandy et al., 2015; Séférían et al., 2014; Tjiputra et al., 2014).

# 3 Results

## 3.1 Comparison of observational data sets

Our review of spin-up protocols for CMIP5 ESM shows that several modeling groups have employed different streams of data sets to initialize their biogeochemical models (e.g., WOA1994, WOA2001), while model evaluation relies on the most up-to-date stream of data. To investigate the differences between the two data streams used for initializing and assessing, respectively, $NO_3$ and $O_2$ concentrations are analyzed. Table 2 summarizes the RMSE and correlation between WOA1994 and WOA2013 for these two biogeochemical fields.

Table 2 indicates that differences between the two streams of data are fairly small. The total difference (RMSE) represents a departure between 5 to 10% from the global average concentrations of WOA2013 across depth levels. It is generally lower in regions where the sampling density has not increased markedly between the two releases. These values can be used as a baseline for model-to-model comparison assuming that errors attributed to the various sources of initialization cannot be larger than 10%. Considering that some models have used outputs from previous model simulations
<table><caption>**Table 2.** Differences between the oxygen ($O_2$, µmol L$^{-1}$) and nitrate ($NO_3$, µmol L$^{-1}$) data sets used for initializing IPSL-CM5A-LR (WOA1994) and the data sets used for assessing its performance (WOA2013).</caption><thead><tr><th rowspan="2">Depth</th><th colspan="3">O<sub>2</sub></th><th colspan="3">NO<sub>3</sub></th></tr><tr><th>Surface</th><th>150 m</th><th>2000 m</th><th>Surface</th><th>150 m</th><th>2000 m</th></tr></thead><tbody><tr><td>RMSE</td><td>7.19</td><td>8.75</td><td>5.50</td><td>2.07</td><td>2.90</td><td>2.08</td></tr><tr><td>R<sup>2</sup></td><td>0.98</td><td>0.98</td><td>0.99</td><td>0.96</td><td>0.92</td><td>0.94</td></tr></tbody></table>

or globally averaged concentrations as initial conditions, we acknowledge that this baseline is not a perfect criterion for benchmarking model performance. There is, however, no ideal solution to address this issue since there is no standardized set of initial conditions in CMIP5 except some recommendations for the decadal prediction exercise in which specific attention was paid to initialization (e.g., Keenlyside et al., 2008; Kim et al., 2012; Matei et al., 2012; Meehl et al., 2013, 2009; Servonnat et al., 2014; Smith et al., 2007; Swingedouw et al., 2013).

## 3.2 Equilibration state metrics in IPSL-CM5A-LR

The global mean sea surface temperature (SST) is a common metric to quantify the energetic equilibrium of the model. This metric has been widely used in various papers referenced in this study to determine the equilibration of ESM physical components. Figure 2a shows the evolution of this metric during the 500-year long spin-up simulation. The global average SST sharply decreases during the first 250 years of the simulation. In the last 250 years of the simulation, the global averaged SST displays a small residual drift of ~$-10^{-4} \text{°C y}^{-1}$, which falls into the range of the drifts reported for CMIP5 ESMs (Sen Gupta et al., 2013). The evolution over the last 250 years is comparable to those of other physical equilibration metrics, such as the ocean heat content or the meridional overturning circulation (Mignot et al., 2013).

To assess if ocean carbon cycle reservoirs are equilibrated, we track the temporal evolution of sea-to-air $CO_2$ fluxes during the spin-up simulation. This metric was used in phase 2 of the Ocean Carbon Model Intercomparison Project (OCMIP-2; Orr, 2002) and has still widely been used during CMIP5 as an equilibration metric for the marine biogeochemistry. Figure 2b presents its evolution in the 500-year long spin-up simulation. The global ocean sea-to-air $CO_2$ flux is ~$-0.7 \text{ Pg C y}^{-1}$ over the last decades of the spin-up simulation (negative values indicate ocean $CO_2$ uptake).

We use the range of values estimated from preindustrial natural ocean carbon flux inversions (e.g. Gerber and Joos, 2010, or Mikaloff Fletcher et al., 2007) to evaluate the global sea-to-air carbon flux simulated by IPSL-CM5A-LR. Since, these estimates do not account for the preindustrial carbon outgassing induced by the river input, while our model does, we have added a constant outgassing of $0.45 \text{ Pg C y}^{-1}$

to the range of $0.03 \pm 0.08 \text{ Pg C y}^{-1}$ (Mikaloff Fletcher et al., 2007). This value of $0.45 \text{ Pg C y}^{-1}$ corresponds to the global open-ocean river-induced carbon outgassing according to IPCC (2013) or Le Quéré et al. (2015). Consequently, in our modeling framework, the target value of the global sea-to-air carbon flux ranges between 0.4 and $0.56 \text{ Pg C y}^{-1}$.

Figure 2b shows that the global sea-to-air carbon flux is still lower than the range of values estimated from preindustrial natural ocean carbon flux inversions (~0.4–0.56 Pg C y$^{-1}$). Besides, Fig. 2b shows that the drift in the global sea-to-air carbon flux becomes smaller more slowly after a strong decline during the first 50 years of the spin-up simulation. From year 250–500 this drift is about 0.001 Pg C y$^{-2}$ and still weaker over the last century of the simulation ($7 \times 10^{-4} \text{ Pg C y}^{-2}$). A one-sided t-test indicates that the two drifts differ from each other with a *p* value <2 × 10$^{-16}$. When fitted with drifts computed from overlapping time segments of 100 years, our simple drift model (Eq. 1) gives a relaxation time of around 160 years. We use this relaxation time and the drift of $7 \times 10^{-4} \text{ Pg C y}^{-2}$ to estimate the additional spin-up time required for the model to reach an outgassing of 0.4–0.56 Pg C y$^{-1}$, as 1100 to 1300 years. However, even after this integration time, the drift in global sea-to-air carbon flux estimated with our simple drift model still ranges from $2 \times 10^{-7}$ to $7 \times 10^{-7} \text{ Pg C y}^{-2}$.

These estimates do not account for the non-linearity of the ocean carbon cycle and the associated process uncertainties (Schwinger et al., 2014), and hence potentially underestimate the time required to equilibrate the ocean carbon cycle and sea-to-air carbon fluxes in the range of inversion estimates. The drift of 0.001 Pg C y$^{-2}$ is, however, much smaller than the oceanic sink for anthropogenic carbon. Even if not fully equilibrated in terms of carbon balance, it is likely that this run would have given consistent estimates of anthropogenic carbon uptake in transient historical hindcasts.

## 3.3 Temporal evolution of model errors in IPSL-CM5A-LR

Figure 3 shows the temporal evolution of globally averaged concentrations for $O_2$, $NO_3$ and Alk-DIC at the surface (Fig. 3a, b and c), 150 m (Fig. 3d, e and f) and 2000 m (Fig. 3g, h, and i). Globally averaged concentrations of $O_2$, $NO_3$ and Alk-DIC (solid lines) reach steady state after 100 to 250 years of spin-up at the surface. While modeled nominal
<figure><img src="image_5.png" /></figure>

**Figure 2.** Time series of two climate indices over the 500-year spin-up simulation of IPSL-CM5A-LR. They represent the globally averaged sea surface temperature (a) and the global mean sea-air carbon flux (b). For sea-air carbon flux, negative value indicates uptake of carbon. Steady state equilibrium of physical components as described in Mignot et al. (2013) is reached at ~250 years and is indicated with a vertical dashed line. Drifts in sea surface temperature and global carbon flux are indicated with dashed blue lines. They are computed using a linear regression fit over years 250 to 500. Hatching on panel (b) represents the range of inverse modeling estimates for preindustrial global carbon flux as described in Mikaloff Fletcher et al. (2007), i.e., $0.03 \pm 0.08$ plus $0.45 \text{ Pg C y}^{-1}$ corresponding to the riverine-induced natural $\text{CO}_2$ outgassing outside of near-shore regions consistently with Le Quéré et al. (2015).

values for $\text{O}_2$ concentration converge toward the observed concentration (i.e., $172.3 \text{ µmol L}^{-1}$), that of $\text{NO}_3$ presents persistent deviations from WOA2013. At the surface, the convergence of the simulated oxygen to observed values is expected since the dominant governing process of thermodynamic saturation (through the air-sea gas exchange) is well understood and modeled. The deviation in surface $\text{NO}_3$ highlights uncertainty related to near surface biological processes and upper ocean physics. Below the surface, concentrations of biogeochemical tracers drift away from the globally averaged concentrations computed from WOA2013 or GLODAP (Fig. 3d-i). At 150 and 2000 m, the drift in global averaged concentrations for these fields, computed over the last 250 years, is still significant with $p < 10^{-4}$ (Table 3). Except for the surface fields, Fig. 3 shows that RMSE, indicated with dashed lines in Fig. 3, globally increases with time for all biogeochemical fields. The linear drift in RMSE over the last 250 years of the spin-up simulation falls within the $2-3 \text{ % ky}^{-1}$ range at the surface. It is much larger at 2000 m ($144-280 \text{ % ky}^{-1}$; Table 3). This is also the case regionally, because the latitudinal maximum in RMSE ($\text{RMSE}_{\text{max}}$) is similar to the global RMSE. Table 3 also shows that the magnitude of drift in RMSE for $\text{O}_2$, $\text{NO}_3$ and Alk-DIC differs at

a given depth as different processes affect the interior distribution of these biogeochemical fields.

## 3.4 Evolution of geographical mismatches in IPSL-CM5A-LR

To further explore the evolution of mismatch in biogeochemical distributions, we analyze differences ($\epsilon$) between simulated and observed fields of $\text{O}_2$ and $\text{NO}_3$ from WOA2013 and Alk-DIC from GLODAP after the initialization and at the end of the spin-up, i.e., the first year and the last year of the core spin-up simulation performed with the IPSL-CM5A-LR model (Figs. 4, 5 and 6).

Figure 4a, c, e shows that surface concentrations of biogeochemical fields are associated with small biases at initialization. This error represents less than 5% of the observed surface concentrations for $\text{O}_2$, $\text{NO}_3$ and Alk-DIC and reflects the weak difference between the data stream employed for initialization and validation. After 500 years of spin-up, deviations between the modeled and observed fields at the surface have increased locally by up to ~40% (Fig. 4b, d, and f). The largest deviations are found in high-latitude oceans for $\text{O}_2$ and $\text{NO}_3$ and also to some extent in the tropics for $\text{NO}_3$ and Alk-DIC.
<figure><img src="image_6.png" /><figcaption>**Figure 3.** Time series of globally averaged concentration (in solid lines) and globally averaged root-mean squared error (RMSE in dashed lines) for dissolved oxygen (O₂), nitrate (NO₃) and difference between alkalinity and dissolved inorganic carbon (Alk-DIC) as simulated by IPSL-CM5A-LR. Globally averaged concentration and RMSE are given at surface (a, b, c), 150 m (d, e, f), and 2000 m (g, h, i) for these three biogeochemical fields. Their values are indicated on the left-side and right-side y axis, respectively. Hatching represents the ±σ observational uncertainty due to optimal interpolation of in situ concentrations around the observed globally averaged concentration.</figcaption></figure>
<table><caption>**Table 3.** Drift in % ky⁻¹ for oxygen (O₂), nitrate (NO₃) and total alkalinity minus DIC (Alk-DIC) at surface, 150 and 2000 m as simulated by the IPSL-CM5A-LR model. The drift has been computed over the last 250 years of the spin-up simulation using a linear regression fit of the globally averaged concentrations, root-mean squared error (RMSE) and latitudinal maximum root-mean squared error (RMSEmax) with respect to the values at year 250.</caption><thead><tr><th rowspan="2">Metrics</th><th colspan="3">O₂</th><th colspan="3">NO₃</th><th colspan="3">Alk-DIC</th></tr><tr><th>Mean</th><th>RMSE</th><th>RMSEmax</th><th>Mean</th><th>RMSE</th><th>RMSEmax</th><th>Mean</th><th>RMSE</th><th>RMSEmax</th></tr></thead><tbody><tr><td>Surf</td><td>-0.2</td><td>2.6</td><td>55.8</td><td>-0.1</td><td>-0.1</td><td>34.2</td><td>1.6</td><td>-0.1</td><td>-0.1</td></tr><tr><td>150 m</td><td>3.4</td><td>39.0</td><td>31.5</td><td>-15.9</td><td>33.4</td><td>55.2</td><td>6.1</td><td>27.9</td><td>24.7</td></tr><tr><td>2000 m</td><td>-30.4</td><td>144.3</td><td>-40.1</td><td>2</td><td>51.8</td><td>-34.8</td><td>-69.6</td><td>281.8</td><td>47.5</td></tr></tbody></table>

<figure><img src="image_7.png" /><figcaption>**Figure 4.** Snapshots of spatial biases, ε, in surface concentrations (µmol L⁻¹) in biogeochemical fields during the 500-year spin-up simulation of IPSL-CM5A-LR. ε in dissolved oxygen (O₂), nitrate (NO₃) and difference between alkalinity and dissolved inorganic carbon (Alk-DIC) is given for the first year (a, c, e, respectively) and for the last year of spin-up simulation (b, d, f, respectively).</figcaption></figure>

Below the surface, distributions of modeled biogeochemical fields compare well to the observations at 150 m at initialization with averaged errors close to zero (Fig. 5a, c, and e). This result was expected since WOA2013 and WOA1994 differ little at these depth levels. Subsurface distributions at initialization strongly contrast with the concentrations that resulted from 500 years of spin-up (Fig. 5b, d, and f). After 500 years of spin-up, substantial mismatches characterize the distribution of O₂, NO₃ and Alk-DIC fields in the high-latitude oceans and in the tropics. Figure 5 illustrates that patterns of errors for O₂, NO₃ and Alk-DIC fields are well correlated with each other (R > 0.6). This reflects that in PISCES carbon, nitrogen and oxygen concentrations are linked by the elemental C:N: -ΔO₂ stoichiometry fixed in space and time. Figure 6 shows that model-data deviations at 2000 m have substantially increased at a regional level after 500 years of simulation, showing large errors in the Southern Hemisphere oceans. This appears clearly in Fig. 6d and f for NO₃ and Alk-DIC fields, respectively.

The temporal evolution of the RMSE between modeled and observed fields of O₂, NO₃ and Alk-DIC over the whole water column is presented in Fig. 7 in terms of RMSE (Fig. 7a-c). As expected, Fig. 7 illustrates that there is a good match during the first years of simulation for all biogeochemical fields at all depth levels with low RMSE. After a few centuries, patterns of error evolve differently across depth for O₂, NO₃ and Alk-DIC.

The temporal evolution of RMSE shows that patterns of error have reached a steady state a few decades after 250 years of spin-up within the upper hundred meters of the ocean but continue to evolve at greater depths, even after 500 years. Patterns of errors within the thermocline and upper 1000 m water masses evolve relatively fast (within a few centuries) due to the relatively short mixing time in the upper ocean. Mid-depth (~1500–2500 m) RMSE evolves much slower because of the slow ocean circulation at these depth levels. Characteristic timescales here are thousands of years as evidenced by the observed radiocarbon age of seawater (e.g.,
<figure><img src="image_8.png" /><figcaption>**Figure 5.** As Fig. 4 but for concentrations at 150 m. Note that color shading does not represent the same amplitude in spatial biases as in Figs. 4 and 6.</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>**Figure 6.** As Fig. 4 but for concentrations at 2000 m. Note that color shading does not represent the same amplitude in spatial biases as in Figs. 4 and 5.</figcaption></figure>

Wunsch and Heimbach, 2007, 2008). This explains why, at the end of the spin-up simulation, two maxima of comparable amplitude are found for RMSE at 150 and 3750 m for O₂ and at 50 and 3800 m for Alk-DIC (Fig. 7).

## 3.5 Drifts in IPSL-CM5A-LR spin-up simulation

With the evolution of the RMSE established, we can use the simple drift model (Eq. 1) to determine the relaxation time, $\tau$, which characterizes the e-folding timescale of the

RMSE. To use this simple drift model, we compute the drift in RMSE determined from time segments of 100 years distributed evenly every 5 years from year 250 to 500 for O₂, NO₃ and Alk-DIC tracers. The drift model (magenta lines in Fig. 8) is fitted to the 80 drift values for each field and each depth level (colored crosses in Fig. 8).

The simple drift model fits well the evolution of the drift in RMSE for the biogeochemical variables along the spin-up simulation of IPSL-CM5A-LR (Fig. 8). Correlation coefficients are mostly significant at 90% confidence level ($r^* = $
<figure><img src="image_10.png" /><figcaption>**Figure 7.** Temporal-vertical evolution in root-mean squared error (RMSE) for biogeochemical tracers during the 500-year long spin-up simulation of IPSL-CM5A-LR. The RMSE is given for (a) dissolved oxygen O₂, (b) nitrate NO₃ and the (c) difference between alkalinity and dissolved inorganic carbon Alk-DIC.</figcaption></figure>

0.3 determined with a student distribution with significance level of 90% and ~15 effective degrees of freedom estimated with the formulation of Bretherton et al., 1999), except for NO₃ at surface and Alk-DIC at 150 and 2000 m. Another exception is found for NO₃ at 150 m where the drift does not correspond to an exponential decay of the drift as function of time. The large confidence interval of the fit indicates that the fit would have been considered as non-significant given a longer spin-up simulation or a higher confidence threshold.

When significant, estimates of τ for O₂ RMSE are ≈ 90, 564 and 1149 y at the surface, 150 and 2000 m, respectively. These values match reasonably well τ estimated for NO₃ RMSE at 2000 m (1130 y) and those for Alk-DIC RMSE at surface and 2000 m (137 and 1163 y). However, these estimates are sensitive to the time windows used to compute the drift. For a subset of time windows between 100 and 250 years by step of 50 years, τ estimates for O₂ RMSE are ≈ 114 ± 67, 375 ± 140 and 1116 ± 527 y at the surface, 150 and 2000 m depth. These large uncertainties associated with τ estimates are essentially due to the length of the spin-up simulation. A longer spin-up simulation would improve the quality of the fit (see Fig. S1).

## 3.6 Drifts in CMIP5 ESM preindustrial simulations

In this subsection, the analysis is extended to the CMIP5 archive. We focus on oxygen fields in the long preindustrial simulation, piControl, for the 15 available CMIP5 ESMs. From these simulations that span from 250 to 1000 years, we compute the drift in O₂ RMSE across depth from several time segments of 100 years distributed evenly every 5 years from the beginning until the end of the piControl simulation. These drifts are used as a surrogate for drift computed from the spin-up of each model since such simulations are not available through the data portal.

Figure 9 represents the drift in O₂ RMSE vs. the spin-up duration for each CMIP5 ESM. The analysis shows that the

drift in O₂ RMSE differs substantially between models. For a given model, drifts in other biogeochemical tracers (NO₃ and Alk-DIC) display similar features (not shown). The between-model differences in drift are not surprising since there are no reasons for different models to exhibit similar drift for a given field. Yet, Fig. 9 shows that a global relationship emerges from this ensemble when using the simple drift model to fit the drift in O₂ RMSE as a function of the spin-up duration (solid green lines in Fig. 9). With a 90% confidence level, this relationship suggests a general decrease in the drift as a function of spin-up duration for all depth levels. At the surface and at 2000 m depth, the quality of fits is low, with correlation coefficients of about 0.4. These are however significant at the 90% confidence level ($r^* = 0.34$ determined with a student distribution with a significance level of 90% and 15 models as a degree of freedom). The weakest correlation coefficient is found for the fit at 150 m depth and hence indicates that there is no link between the drift in O₂ RMSE and the duration of the spin-up simulation. This low significance level must be put into perspective given the large diversity of spin-up protocols and initial conditions (Fig. 1 and Table 1) that can deteriorate the drift-spin-up duration relationship in this ensemble of models.

The drift vs. spin-up duration relationship established from the 15 CMIP5 ESMs is nonetheless consistent with the results obtained with IPSL-CM5A-LR (the results in Fig. 8 have been reported in Fig. 9 with magenta crosses). Indeed, the drifts in RMSE decrease in the course of time at the various depth levels for the IPSL-CM5A-LR model, although their magnitudes differ. This difference in magnitude is not surprising if one considers that drift is highly model and protocol dependent and that the length of the IPSL-CM5A-LR spin-up simulation is potentially too short to determine accurate estimates of the long-term drift in O₂ RMSE. Despite these differences, our analyses show that a relationship between the drift in O₂ RMSE vs. the spin-up duration emerges
<figure><img src="image_11.png" /></figure>

**Figure 8.** Temporal evolution of drift in root-mean squared error (RMSE) for dissolved oxygen (O₂, blue crosses), nitrate (NO₃, green crosses) and difference between alkalinity and dissolved inorganic carbon (Alk-DIC, orange crosses) during the 500-year long spin-up simulation of IPSL-CM5A-LR. Drift in RMSE is given at surface (a, b, c), 150 m (d, e, f), and 2000 m (g, h, i) for these three biogeochemical fields. Drift in RMSE is computed from time segments of 100 years beginning every 5 years from the beginning until year 400 of the spin-up simulation for O₂, NO₃ and Alk-DIC tracers. The best-fit regressions between drifts in RMSE and spin-up duration over year 250 to 500 are indicated in solid magenta lines; their 90% confidence intervals are given by thin dashed envelopes. Least square correlation coefficients are tested against a one-tailed t-test with significance level of 90% and ~15 effective degrees of freedom estimated with the formulation of Bretherton et al. (1999); * indicates if a given least square correlation coefficient passes the test.
<figure><img src="image_12.png" /><figcaption>**Figure 9.** Scatterplot of drifts in root-mean squared error (RMSE) in O₂ concentration vs. the duration of the spin-up simulation for the available CMIP5 Earth system models. Drifts in O₂ RMSE are respectively given for surface (a), 150 m (b) and 2000 m (c) for oxygen concentrations. Drift in O₂ RMSE is computed from several time segments of 100 years beginning every 5 years from the beginning until the end of the piControl simulation for the available CMIP5 models. Colored symbols indicate the mean drift in O₂ RMSE, while vertical lines represent the associated 90% confidence interval. The best-fit regressions between models’ mean drifts in RMSE and spin-up duration are indicated as solid green lines; their 90% confidence intervals are given by thin dashed envelopes. Fits are assumed robust if correlation coefficients are significant at 90% (i.e., r* > 0.34). For comparison, drift in O₂ RMSE from our spin-up simulation with IPSL-CM5A-LR (Fig. 8) are represented by magenta crosses.</figcaption></figure>
from an ensemble of models and is broadly consistent with our theoretical framework of a drift model established from the results of the IPSL-CM5A-LR model (Fig. 8).

## 3.7 Impact of the drift on model skill score assessment metrics across CMIP5 ESMs

In the following, we investigate the influence of model drift on skill score assessment metrics that are routinely used to benchmark model performance. For this purpose, we use the ensemble-mean $O_2$ RMSE as a metric to assess the distance between the biogeochemical observations and model results. For this purpose, we compute $O_2$ RMSE from each ensemble member of the CMIP5 models averaged from 1986 to 2005 with respect to WOA2013 observations. The model-data distance is then determined for each CMIP5 model using the mean across the available ensemble members.

The left-hand side panels of Fig. 10 present the performance of available CMIP5 models in terms of distance to oxygen observations at the surface, 150 and 2000 m, respectively. In these panels, the various CMIP5 models are ordered as function of their distance to the oxygen observations. Following Knutti et al. (2013), either the ensemble mean or the ensemble median is used to identify groups of models with similar skill within the CMIP5 ensemble. The left-hand side panels of Fig. 10 show that the ability of models to reproduce oxygen observations varies across depth levels. The RMSE in the simulated $O_2$ fields in CESM1-BGC, HadGEM2-ES, HadGEM2-CC, GFDL-ESM2M, MPI-ESM-LR and MPI-ESM-MR is generally smaller than the ensemble mean or ensemble median RMSE across the various depth levels (Fig. 10a, b and c). On the other side of the ranking, CMCC-CESM, CNRM-CM5, CNRM-CM5-2, IPSL-CM5B-LR and NorESM1-ME exhibit RMSE generally higher than the ensemble mean and median RMSE across the various depth levels. The other models, i.e., CNRM-ESM1, GFDL-ESM2G, IPSL-CM5A-LR and IPSL-CM5A-MR display $O_2$ RMSE that is generally close to the ensemble mean or the ensemble median.

To assess the impact of model's drift inherited from the diversity of spin-up strategies (Fig. 1 and Table 1) on the performance metrics, we use a simple additive assumption to incorporate an incremental error due to the drift, $\Delta$RMSE, to the above-mentioned RMSE. This incremental error due to the drift is computed using the relaxation time $\tau$ determined from the piControl simulations of each CMIP5 model at each depth level (Eq. 1 and Fig. 9) and a common duration of $T = 3000$ years for all models ($m$):

$$ \Delta \text{RMSE}_m(z) = \int_0^T \text{drift}_m(z, t=0) \times \exp\left(-\frac{1}{\tau(z)}t\right) dt \quad (2) $$

where $\Delta$RMSE has the same unit as RMSE.

The common duration $T$ is used to bring model drift close to zero and hence to make models comparable to each other.

We employ $\Delta$RMSE to penalize the distance from the observations assuming that this drift-induced deviation in tracer fields can be added to RMSE. This means that the effect of the penalty is to increase the distance giving a consistent measure of the equilibration error.

The right-hand side panels of Fig. 10 show the influence of this penalization approach on the model ranking at the various depth levels. They show that several models have been upgraded in the ranking while others have not. For example, both MPI-ESM-LR, MPI-ESM-MR have been upgraded at the surface and 2000 m. On the other hand, the rank of HadGEM2-ES and HadGEM2-CC has been downgraded to the 5th and 3rd position due to the large drift in surface oxygen concentrations in comparison to that of the other models. The surface drift might be attributed to drivers in oxygen fluxes (e.g., SST, SSS). The ranking of GFDL-ESM2G and GFDL-ESM2M slightly changes with penalization but both models stay close to the ensemble mean or the ensemble median. At the bottom of the ranking, models with large deviation from the oxygen observations (i.e., CMCC-CESM, IPSL-CM5B-LR, NorESM1-ME, CNRM-CM5) are found. For these models, the computed $\Delta$RMSE and RMSE result in similar ranking, because even a small drift and hence relatively low $\Delta$RMSE cannot compensate for their large RMSE.

# 4 Discussion

## 4.1 Implications for biogeochemical processes

Our results show that errors in ocean biogeochemical fields amplify during the spin-up simulation but not at the same rate at all depths. These differences in error evolution are consistent with an increasing contribution of biogeochemical processes in setting the distribution of tracers at depth. Indeed, Mignot et al. (2013) with the same model simulation showed that the main physical climate fields as well as the large-scale ocean circulation reach quasi-equilibrium after 250 years of spin-up, but our analyses indicate that biogeochemical tracers do not (Fig. 3).

Besides, our analysis demonstrates that drifts in biogeochemical fields are highly model dependent. For example, despite having the same initialization strategy and comparable spin-up duration, the GFDL-ESM2G, GFDL-ESM2M, and NorESM1-ME models display considerable difference in drift (Figs. 9 and 10) that mirror large differences in model performance and properties (e.g., resolution, simulated processes).

The identification of the dynamical or biogeochemical processes responsible for these errors is not within the scope of this study and would require (or “would have required”) additional long simulations with additional tracers targeted for attribution of the various biogeochemical processes and the underlying ocean physics (e.g., Doney et al., 2004) involved (e.g. using abiotic, passive tracers as suggested in
<figure><img src="image_13.png" /><figcaption>**Figure 10.** Rankings of CMIP5 Earth system models based on the standard and penalized version of the distance from oxygen observations. The standard distance metric is calculated as the ensemble-mean root-mean squared error (RMSE) for $O_2$ concentrations at surface (a), 150 m (b) and 2000 m (c). The penalized distance metric incorporates drift-induced changes in $O_2$ RMSE ($\Delta$RMSE) into $O_2$ RMSE at surface (d), 150 m (e) and 2000 m (f). Ensemble-mean RMSEs are calculated using available ensemble members of Earth system model oxygen concentrations averaged over the 1986–2005 historical period relative to WOA2013 observations. $\Delta$RMSE is determined using Eq. (2) and fits derived from the first century of the CMIP5 piControl simulations. Solid red and magenta lines indicate the multi-model mean standard and penalized distance from $O_2$ observations, respectively. With the same color pattern, dashed lines are indicative of the multi-model median for the standard and penalized distance from $O_2$ observations.</figcaption></figure>

Walin et al., 2014). Some mechanisms can be nonetheless invoked to explain differences or similarities in behavior between biogeochemical fields. For example, the evolution of surface concentrations for $O_2$ and Alk-DIC is controlled in part by the solubility of $O_2$ and $CO_2$ in seawater and the concentration of these gases in the atmosphere (set to the observed values and kept constant in all experiments performed with IPSL-CM5A-LR discussed here) and the biological soft-tissue and calcium carbonate counter pumps (in relation with the vertical transport of nutrients and alkalinity). Therefore, the equilibration of the $O_2$ and Alk-DIC surface fields once the physical equilibrium is to a large degree reached (~250 years of spin-up) is expected (Figs. 3a, c and 7). Nevertheless, spatial errors could increase depending

on the physical state of the model (Fig. 4b and f). By contrast, the evolution of $NO_3$ concentration is predominantly determined by ocean circulation, biological processes, and to a lesser extent by external supplies from rivers and atmosphere. Below the surface, concentrations of $O_2$, $NO_3$, and Alk-DIC evolve in response to the combined effect of ocean circulation and biogeochemical processes. The combination of dynamical and biogeochemical processes on the one hand, and the spin-up strategy on the other hand both shape the modeled distributions of large-scale biogeochemical tracers.

Consequences of the difficulty in achieving the correct equilibration procedure have important implications for biogeochemical features that are defined by regional characteristics in tracer concentrations, such as high nutrient/low
chlorophyll regions, oxygen minimum zones and nutrient-to-light colimitation patterns. This point is illustrated by recent studies focusing on future changes in phytoplankton productivity (e.g., Vancoppenolle et al., 2013, and Laufkötter et al., 2015). Vancoppenolle and co-workers report a wide spread of surface mean NO₃ concentrations (1980–1999) in the Arctic with a range from 1.7 to 8.9 µmol L⁻¹ across a subset of 11 CMIP5 models. The spread in present-day NO₃ concentrations translates into a large model-to-model uncertainty in future net primary production. Laufkötter et al. (2015) determined limitation terms of phytoplankton production for a subset of CMIP5 and MAREMIP (Marine Ecosystem Model Intercomparison Project) models. The authors demonstrate that nutrient-to-light colimitation patterns differ in strength, location and type between models and arise from large differences in the simulated nutrient concentrations. Although Vancoppenolle et al. (2013) and Laufkötter et al. (2015) explain a part of the difference in simulated nutrient concentration by the differences in the spatial resolution and the complexity of the models, the authors of both studies qualitatively invoked differences in spin-up duration to explain the remaining differences in simulated concentrations. Besides, a recent assessment of interannual to decadal variability of ocean–atmosphere CO₂ and O₂ fluxes in CMIP5 models, suggests that decadal variability can range regionally from 10 to 50% of the total natural variability among a subset of 6 ESMs (Resplandy et al., 2015). In that study, the authors demonstrate that, despite the robustness of driving mechanisms (mostly related to vertical transport of water masses) across the model ensemble, model-to-model spread can be related to differences in modeled carbon and oxygen concentrations. In light of present results, it appears likely that differences in spin-up strategy and sources of initialization could also contribute to the amplitude of the natural variability of the ocean CO₂ and O₂ fluxes.

## 4.2 Implications for future projections

The inconsistent strategies used to spin up models in CMIP5 have resulted in a significant source of uncertainty to the multi-model spread. It needs to be better constrained in order to draw robust conclusions on the impact of climate change on the carbon cycle as well as on climate feedback (e.g., Arora et al., 2013; Friedlingstein et al., 2013; Roy et al., 2011; Schwinger et al., 2014; Séférian et al., 2012) and on marine ecosystems (e.g., Bopp et al., 2013; Boyd et al., 2015; Cheung et al., 2012; Doney et al., 2012; Gattuso et al., 2015; Lehodey et al., 2006). So far, the most frequently used approach relies on long preindustrial control simulations running parallel to a transient simulations, allowing the “removal” of the drift in the simulated fields over the historical period or future projections (e.g., Bopp et al., 2013; Cocco et al., 2013; Friedlingstein et al., 2013, 2006; Frölicher et al., 2014; Gehlen et al., 2014; Keller et al., 2014; Steinacher et al., 2010; Tjiputra et al., 2014). Although this approach

allows one to determine relative changes, it does not allow one to investigate the underlying reasons for the spread between models in terms of processes, variability and response to climate change. The “drift-correction” approach, much as the one used for this study, assumes that drift-induced errors in the simulated fields can be isolated from the signal of interest. Verification of this fundamental hypothesis would require a specific experimental set-up consisting of the perturbation of model fields (e.g., nutrients or carbon-related fields) to assess by how much the model projections would be modified. So far, several modeling groups have generated ensemble simulations in CMIP5 using a perturbation approach. However, the perturbations were applied either to physical fields only or to both the physical and marine biogeochemical fields. To assess impacts of different spin-up strategies and/or initial conditions on future projections of marine biogeochemical tracer distributions, ensemble simulations in which only biogeochemical fields are perturbed would be needed.

## 4.3 Implications for multi-model skill-score assessments

While the importance of spin-up protocols is well accepted in the modeling community, the link between spin-up strategy and the ability of a model to reproduce modern observations remains to be addressed.

Most of the recent CMIP5 skill assessment approaches were based on *historical* hindcasts that were started from preindustrial runs of varying duration and from various spin-up strategies. Therefore, in typical intercomparison exercises, Earth system models with a short spin-up, and hence modeled distributions still close to initial fields, are confronted with Earth system models with a longer spin-up duration and modeled distributions that have drifted further away from their initial states. Our study highlights that such inconsistencies in spin-up protocols and initial conditions across CMIP5 Earth system models (Fig. 1 and Table 1) could significantly contribute to model-to-model spread in performance metrics. The analysis of the first century of CMIP5 *piControl* simulations demonstrated a significant spread of drift between CMIP5 models (Fig. 9). An approximate exponential relationship between the amplitude of the drift and the spin-up duration emerges from the ensemble of CMIP5 models, which is consistent with results from IPSL-CM5A-LR. For example, while the global average root-mean square error increased up to 70% during a 500-year spin-up simulation with IPSL-CM5A-LR, its rate of increase (or drift) decreased with time to a very small rate (0.001 Pg C y⁻²). Combining a simple drift model and this relationship, we propose a penalization approach in an effort to assess more objectively the influence of documented model differences on model-data biases. Figure 10 compares the standard approach to assessing model performance (left-hand side panels) to the drift-penalized approach (right-hand side panels).
This novel approach penalizes models with larger drift without affecting the models with smaller drift. Taking into account drift in modeled fields results in subtle adjustments in ranking, which reflect differences in spin-up and initialization strategies.

## 4.4 Limitations of the framework

In this work, the analyses focus on the globally averaged O₂ RMSE across a diverse ensemble of CMIP5 models, which differ in terms of represented processes, spatial resolution and performance in addition to differences in spin-up protocols. Major limitations of the framework are presented below.

Due to their specificities in terms of processes and resolution (e.g., Cabré et al., 2015; Laufkötter et al., 2015), regional drift in CMIP5 models may differ from the drift computed from globally averaged skill-score metrics (see Figs. S2 and S3). These differences may lead to different estimates of the relaxation time τ at a regional scale. Moreover, the combination of regional ocean physics and biogeochemical processes in each individual model may drive an evolution of a regional drift in RMSE that does not fit the hypothesis of an exponential decay of the drift during the course of the spin-up simulation.

Besides, a difference in the simulated processes and resolution in the different models can explain the relatively low confidence level of the fit to drift across the multi-model CMIP5 ensemble (Fig. 9). The relatively low significance level of the fit reflects not only the large diversity of spin-up protocols and initial conditions (Fig. 1 and Table 1), but also the large diversity of processes and resolution of the CMIP5 models. Indeed, as shown in Kriest and Oschlies (2015), various parameterizations of the particle sinking speed in a common physical framework may lead to a similar evolution of the globally averaged RMSE in the first century of the spin-up simulation but display very different behavior within a timescale of O(10³) years. As such, drift and τ estimates need to be used with caution when computed from short spin-up simulations because they can be subject to large uncertainties. An improved derivation of the penalization would require access to output from spin-up simulations for each individual model or, at least, a better quantification of model–model differences in terms of initial conditions.

# 5 Conclusions and recommendation for future intercomparison exercises

Skill-score metrics are expected to be widely used in the framework of the upcoming CMIP6 (Meehl et al., 2014) with the development of international community benchmarking tools like the ESMValTool (http://www.pa.op.dlr.de/ESMValTool; see also Eyring et al., 2015). The assessment of model skill to reproduce observations will focus on the mod-

ern period. Complementary to this approach, our results call for the consideration of spin-up and initialization strategies in the determination of skill assessment metrics (e.g., Friedrichs et al., 2009; Stow et al., 2009) and, by extension, in model weighting (e.g., Steinacher et al., 2010) and model ranking (e.g., Anav et al., 2013). Indeed, the use of equilibrium-state metrics of the model like the three-dimensional drift of relevant skill-score metrics (e.g., RMSE) could be employed to increase the reliability of these traditional metrics and, as such, should be included in the set of standard assessment tools for CMIP6.

In an effort to better represent interactions between marine biogeochemistry and climate (Smith et al., 2014), future generations of Earth system models are likely to include more complex ocean biogeochemical models, be it in terms of processes (e.g., Tagliabue and Völker, 2011; Tagliabue et al., 2011) or interactions with other biogeochemical cycles (e.g., Gruber and Galloway, 2008) or increased spatial resolution (e.g., Dufour et al., 2013; Lévy et al., 2012) in order to better represent mesoscale biogeochemical dynamics. These developments will go along with an increase in the diversity and complexity of spin-up protocols applied to Earth system models, especially those including an interactive atmospheric CO₂ or interactive nitrogen cycle (e.g., Dunne et al., 2013; Lindsay et al., 2014). The additional challenge of spinning-up emission-driven simulations with interactive carbon cycle will also require us to extend the assessment of the impact of spin-up protocols to the terrestrial carbon cycle. Processes such as soil carbon accumulation, peat formation as well as shift in biomes such as tropical and boreal ecosystems for dynamic vegetation models require several long timescales to equilibrate (Brovkin et al., 2010; Koven et al., 2015). In addition, the terrestrial carbon cycle has large uncertainties in terms of carbon sink/source behavior (Anav et al., 2013; Dalmonech et al., 2014; Friedlingstein et al., 2013), which might affect ocean CO₂ uptake (Brovkin et al., 2010). A novel numerical algorithm to accelerate the spin-up integration time for computationally expensive ocean biogeochemical models has emerged (Khatiwala, 2008), which could help to disentangle physical from biogeochemical contributions to the inter-model spread, but at the same time, could also potentially complicate the determination of inter-model spread by increasing the diversity of spin-up protocols.

To evaluate the contribution of variable spin-up and initialization strategies to model performance, these should be documented extensively and the corresponding model output should be archived. Ideally, for future coupled model intercomparison exercises (i.e., CMIP6, CMIP7, Meehl et al., 2014), the community should agree on a set of simple recommendations for spin-up protocols, following past projects such as OCMIP-2. In parallel, any trade-off between model equilibration and computationally efficient spin-up procedures has to be linked with efforts to reduce model errors due to the physical and biogeochemical parameterizations.
The Supplement related to this article is available online
at doi:10.5194/gmd-9-1827-2016-supplement.

Acknowledgements. We sincerely thank I. Kriest, F. Joos, the anonymous reviewer and A. Yool for their useful comments on this paper. This work was supported by H2020 project CRESCENDO “Coordinated Research in Earth Systems and Climate: Experiments, kKnowedge, Dissemination and Outreach”, which received funding from the European Union’s Horizon 2020 research and innovation programme under grant agreement no. 641816 and by the EU FP7 project CARBOCHANGE “Changes in carbon uptake and emissions by oceans in a changing climate” which received funding from the European community’s Seventh Framework Programme under grant agreement no. 264879. Supercomputing time was provided by GENCI (Grand Equipement National de Calcul Intensif) at CCRT (Centre de Calcul Recherche et Technologie), allocation 016178. Finally, we are grateful to the ESGF project which makes data available for all the community. Roland Séférian is grateful to Aurélien Ribes for his kind advices on statistics. Jerry Tjiputra acknowledges ORGANIC project (239965/F20) funded by the Research Council of Norway. Christoph Heinze and Jerry Tjiputra are grateful for support through project EVA – Earth system modelling of climate variations in the Anthropocene (229771/E10) funded by the Research Council of Norway, as well as CPU-time and mass storage provided through NOTUR project NN2345K as well as NorStore project NS2345K. Keith Lindsay and Scott C. Doney acknowledge support from the National Science Foundation.

Edited by: A. Yool

References

Adachi, Y., Yukimoto, S., Deushi, M., Obata, A., Nakano, H., Tanaka, T. Y., Hosaka, M., Sakami, T., Yoshimura, H., Hirabara, M., Shindo, E., Tsujino, H., Mizuta, R., Yabu, S., Koshiro, T., Ose, T., and Kitoh, A.: Basic performance of a new earth system model of the Meteorological Research Institute (MRI-ESM1), Pap. Meteorol. Geophys., 64, 1–18, doi:10.2467/mripapers.64.1, 2013.

Anav, A., Friedlingstein, P., Kidston, M., Bopp, L., Ciais, P., Cox, P., Jones, C., Jung, M., Myneni, R., and Zhu, Z.: Evaluating the Land and Ocean Components of the Global Carbon Cycle in the CMIP5 Earth System Models, J. Climate, 26, 6801–6843, doi:10.1175/JCLI-D-12-00417.1, 2013.

Andrews, O. D., Bindoff, N. L., Halloran, P. R., Ilyina, T., and Le Quéré, C.: Detecting an external influence on recent changes in oceanic oxygen using an optimal fingerprinting method, Biogeosciences, 10, 1799–1813, doi:10.5194/bg-10-1799-2013, 2013.

Archer, D., Buffett, B., and Brovkin, V.: Ocean methane hydrates as a slow tipping point in the global carbon cycle, P. Natl. Acad. Sci., 106, 20596–20601, 2009.

Arora, V. K., Scinocca, J. F., Boer, G. J., Christian, J. R., Denman, K. L., Flato, G. M., Kharin, V. V., Lee, W. G., and Merryfield, W.

J.: Carbon emission limits required to satisfy future representative concentration pathways of greenhouse gases, Geophys. Res. Lett., 38, L05805, doi:10.1029/2010GL046270, 2011.

Arora, V. K., Boer, G. J., Friedlingstein, P., Eby, M., Jones, C. D., Christian, J. R., Bonan, G., Bopp, L., Brovkin, V., Cadule, P., Hajima, T., Ilyina, T., Lindsay, K., Tjiputra, J. F., and Wu, T.: Carbon-Concentration and Carbon-Climate Feedbacks in CMIP5 Earth System Models, J. Climate, 26, 5289–5314, doi:10.1175/JCLI-D-12-00494.1, 2013.

Aumont, O. and Bopp, L.: Globalizing results from ocean in situ iron fertilization studies, Global Biogeochem. Cy., 20, GB2017, doi:10.1029/2005GB002591, 2006.

Aumont, O., Orr, J., Jamous, D., Monfray, P., Marti, O., and Madec, G.: A degradation approach to accelerate simulations to steady-state in a 3-D tracer transport model of the global ocean, Clim. Dynam., 14, 101–116, 1998.

Aumont, O., Orr, J. C., Monfray, P., Ludwig, W., Amiotte-Suchet, P., and Probst, J.-L.: Riverine-driven interhemispheric transport of carbon, Global Biogeochem. Cy., 15, 393–405, doi:10.1029/1999GB001238, 2001.

Aumont, O., Maier-Reimer, E., Blain, S., and Monfray, P.: An ecosystem model of the global ocean including Fe, Si, P colimitations, Global Biogeochem. Cy., 17, 1060, doi:10.1029/2001GB001745, 2003.

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, doi:10.5194/bg-10-6225-2013, 2013.

Boyd, P. W., Lennartz, S. T., Glover, D. M., and Doney, S. C.: Biological ramifications of climate-change-mediated oceanic multi-stressors, Nature Clim. Change, 5, 71–79, 2015.

Bretherton, C. S., Widmann, M., Dymnikov, V. P., Wallace, J. M., and Bladé, I.: The Effective Number of Spatial Degrees of Freedom of a Time-Varying Field, J. Climate, 12, 1990–2009, doi:10.1175/1520-0442(1999)012<1990:TENOSD>2.0.CO;2, 1999.

Brovkin, V., Lorenz, S. J., Jungclaus, J., Raddatz, T., Timmreck, C., Reick, C. H., Segschneider, J., and Six, K.: Sensitivity of a coupled climate-carbon cycle model to large volcanic eruptions during the last millennium, Tellus B, 62, 674–681, doi:10.1111/j.1600-0889.2010.00471.x, 2010.

Bryan, K.: Accelerating the Convergence to Equilibrium of Ocean-Climates Models, J. Phys. Oceanogr., 14, 666–673, doi:10.1175/1520-0485(1984)014<0666:ATCTEO>2.0.CO;2, 1984.

Cabré, 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, doi:10.5194/bg-12-5429-2015, 2015.

Cheung, W. W. L., Sarmiento, J. L., Dunne, J. P., Frölicher, T. L., Lam, V. W. Y., Palomares, M. L. D., Watson, R., and Pauly, D.: Shrinking of fishes exacerbates impacts of global ocean changes on marine ecosystems, Nature Climate Change, 2, 1–5, doi:10.1038/nclimate1691, 2012.

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, doi:10.5194/bg-10-1849-2013, 2013.

Collins, W. J., Bellouin, N., Doutriaux-Boucher, M., Gedney, N., Halloran, P., Hinton, T., Hughes, J., Jones, C. D., Joshi, M., Lidicoat, S., Martin, G., O’Connor, F., Rae, J., Senior, C., Sitch, S., Totterdell, I., Wiltshire, A., and Woodward, S.: Development and evaluation of an Earth-System model – HadGEM2-Geosci. Model Dev., 4, 1051–1075, doi:10.5194/gmd-4-1051-2011, 2011.

Cox, P. M., Pearson, D., Booth, B. B., Friedlingstein, P., Huntingford, C., Jones, C. D., and Luke, C. M.: Sensitivity of tropical carbon to climate change constrained by carbon dioxide variability, Nature, 494, 341–344, doi:10.1038/nature11882, 2013.

Dalmonech, D., Foley, A. M., Anav, A., Friedlingstein, P., Friend, A. D., Kidston, M., Willeit, M., and Zaehle, S.: Challenges and opportunities to reduce uncertainty in projections of future atmospheric CO₂: a combined marine and terrestrial biosphere perspective, Biogeosciences Discuss., 11, 2083–2153, doi:10.5194/bgd-11-2083-2014, 2014.

de Baar, H. J. W. and de Jong, J. T. M.: The biogeochemistry of iron in seawater, edited by: Turner, D. R. and Hunter, K. A., John Wiley, Hoboken, NJ, USA, 2001.

Doney, S. C.: The Growing Human Footprint on Coastal and Open-Ocean Biogeochemistry, Science, 328, 1512–1516, doi:10.1126/science.1185198, 2010.

Doney, S. C., Lindsay, K., Caldeira, K., Campin, J.-M., Drange, H., Dutay, J.-C., Follows, M., Gao, Y., Gnanadesikan, A., Gruber, N., Ishida, A., Joos, F., Madec, G., Maier-Reimer, E., Marshall, J. C., Matear, R. J., Monfray, P., Mouchet, A., Najjar, R., Orr, J. C., Plattner, G.-K., Sarmiento, J., Schlitzer, R., Slater, R., Totterdell, I. J., Weirig, M.-F., Yamanaka, Y., and Yool, A.: Evaluating global ocean carbon models: The importance of realistic physics, Global Biogeochem. Cy., 18, GB3017, doi:10.1029/2003GB002150, 2004.

Doney, S. C., Lima, I., Moore, J. K., Lindsay, K., Behrenfeld, M. J., Westberry, T. K., Mahowald, N., Glover, D. M., and Takahashi, T.: Skill metrics for confronting global upper ocean ecosystem-biogeochemistry models against field and remote sensing data, J. Marine Syst., 76, 95–112, doi:10.1016/j.jmarsys.2008.05.015, 2009.

Doney, S. C., Ruckelshaus, M., Emmett Duffy, J., Barry, J. P., Chan, F., English, C. A., Galindo, H. M., Grebmeier, J. M., Hollowed, A. B., Knowlton, N., Polovina, J., Rabalais, N. N., Sydeman, W. J., and Talley, L. D.: Climate Change Impacts on Marine Ecosystems, Annu. Rev. Marine. Sci., 4, 11–37, doi:10.1146/annurev-marine-041911-111611, 2012.

Dufour, C. O., Sommer, J. L., Gehlen, M., Orr, J. C., Molines, J.-M., Simeon, J., and Barnier, B.: Eddy compensation and controls of the enhanced sea-to-air CO₂ flux during positive phases of the Southern Annular Mode, Global Biogeochem. Cy., 27, 950–961, doi:10.1002/gbc.20090, 2013.

Dufresne, J.-L., Foujols, M. A., Denvil, S., Caubel, A., Marti, O., Aumont, O., Balkanski, Y., Bekki, S., Bellenger, H., Benshila, R., Bony, S., Bopp, L., Braconnot, P., Brockmann, P., Cadule, P., Cheruy, F., Codron, F., Cozic, A., Cugnet, D., Noblet, N., Duvel, J. P., Ethe, C., Fairhead, L., Fichefet, T., Flavoni, S., Friedlingstein, P., Grandpeix, J. Y., Guez, L., Guilyardi, E., Hauglustaine, D., Hourdin, F., Idelkadi, A., Ghattas, J., Joussaume, S.,

Kageyama, M., Krinner, G., Labetoulle, S., Lahellec, A., Lefebvre, M.-P., Lefèvre, F., Lévy, C., Li, Z. X., Lloyd, J., Lott, F., Madec, G., Mancip, M., Marchand, M., Masson, S., Meurdesoif, Y., Mignot, J., Musat, I., Parouty, S., Polcher, J., Rio, C., Schulz, M., Swingedouw, D., Szopa, S., Talandier, C., Terray, P., Viovy, N., and Vuichard, N.: Climate change projections using the IPSL-CM5 Earth System Model: from CMIP3 to CMIP5, Clim. Dynam., 40, 2123–2165, doi:10.1007/s00382-012-1636-1, 2013.

Dunne, J. P., John, J. G., Adcroft, A. J., Griffies, S. M., Hallberg, R. W., Shevliakova, E., Stouffer, R. J., Cooke, W., Dunne, K. A., Harrison, M. J., Krasting, J. P., Malyshev, S. L., Milly, P. C. D., Phillipps, P. J., Sentman, L. A., Samuels, B. L., Spelman, M. J., Winton, M., Wittenberg, A. T., and Zadeh, N.: GFDL’s ESM2 Global Coupled Climate–Carbon Earth System Models. Part I: Physical Formulation and Baseline Simulation Characteristics, J. Climate, 25, 6646–6665, doi:10.1175/JCLI-D-11-00560.1, 2013.

Duplessy, J. C., Bard, E., Arnold, M., Shackleton, N. J., Duprat, J., and Labeyrie, L.: How fast did the ocean—atmosphere system run during the last deglaciation?, Earth Planet. Sc. Lett., 103, 27–40, doi:10.1016/0012-821X(91)90147-A, 1991.

Eyring, V., Righi, M., Evaldsson, M., Lauer, A., Wenzel, S., Jones, C., Anav, A., Andrews, O., Cionni, I., Davin, E. L., Deser, C., Ehbrecht, C., Friedlingstein, P., Gleckler, P., Gottschaldt, K.-D., Hagemann, S., Juckes, M., Kindermann, S., Krasting, J., Kunert, D., Levine, R., Loew, A., Mäkelä, J., Martin, G., Mason, E., Phillips, A., Read, S., Rio, C., Roehrig, R., Senftleben, D., Sterl, A., van Ulft, L. H., Walton, J., Wang, S., and Williams, K. D.: ESMValTool (v1.0) – a community diagnostic and performance metrics tool for routine evaluation of Earth System Models in CMIP, Geosci. Model Dev. Discuss., 8, 7541–7661, doi:10.5194/gmdd-8-7541-2015, 2015.

Fichefet, T. and Maqueda, M. A. M.: Sensitivity of a global sea ice model to the treatment of ice thermodynamics and dynamics, J. Geophys. Res., 102, 12609–12646, 1997.

Follows, M. J., Dutkiewicz, S., Grant, S., and Chisholm, S. W.: Emergent Biogeography of Microbial Communities in a Model Ocean, Science, 315, 1843–1846, doi:10.1126/science.1138544, 2007.

Friedlingstein, P., Cox, P., Betts, R., Bopp, L., Bloh, Von, W., Brovkin, V., Cadule, P., Doney, S., Eby, M., Fung, I., Bala, G., John, J., Jones, C., Joos, F., Kato, T., Kawamiya, M., Knorr, W., Lindsay, K., Matthews, H. D., Raddatz, T., Rayner, P., Reick, C., Roeckner, E., Schnitzler, K. G., Schnur, R., Strassmann, K., Weaver, A. J., Yoshikawa, C., and Zeng, N.: Climate–Carbon Cycle Feedback Analysis: Results from the C⁴MIP Model Intercomparison, J. Climate, 10, 3337–3353, doi:10.1175/JCLI3800.1, 2006.

Friedlingstein, P., Meinshausen, M., Arora, V. K., Jones, C. D., Anav, A., Liddicoat, S. K., and Knutti, R.: Uncertainties in CMIP5 climate projections due to carbon cycle feedbacks, J. Climate, 27, 511–526, doi:10.1175/JCLI-D-12-00579.1, 2013.

Friedrichs, M. A. M., Dusenberry, J. A., Anderson, L. A., Armstrong, R. A., Chai, F., Christian, J. R., Doney, S. C., Dunne, J. P., Fujii, M., Hood, R., McGillicuddy Jr., D. J., Moore, J. K., Schartau, M., Spitz, Y. H., and Wiggert, J. D.: Assessment of skill and portability in regional marine biogeochemical models: Role of multiple planktonic groups, J. Geophys. Res., 112, C08001, doi:10.1029/2006JC003852, 2007.
Friedrichs, M. A. M., Carr, M.-E., Barber, R. T., Scardi, M., An-
toine, D., Armstrong, R. A., Asanuma, I., Behrenfeld, M. J.,
Buitenhuis, E. T., Chai, F., Christian, J. R., Ciotti, A. M., Doney,
S. C., Dowell, M., Dunne, J. P., Gentili, B., Gregg, W., Hoepffner,
N., Ishizaka, J., Kameda, T., Lima, I., Marra, J., Mélin, F., Moore,
J. K., Morel, A., O'Malley, R. T., O'Reilly, J., Saba, V. S.,
Schmeltz, M., Smyth, T. J., Tjiputra, J., Waters, K., Westberry,
T. K., and Winguth, A.: Assessing the uncertainties of model es-
timates of primary productivity in the tropical Pacific Ocean, J.
Marine Syst., 76, 113–133, doi:10.1016/j.jmarsys.2008.05.010,
2009.

Frölicher, T. L., Sarmiento, J. L., Paynter, D. J., Dunne, J. P., Krast-
ing, J. P., and Winton, M.: Dominance of the Southern Ocean
in anthropogenic carbon and heat uptake in CMIP5 models, J.
Climate, 28, 862–886, doi:10.1175/JCLI-D-14-00117.1, 2014.

Gattuso, J. P., Magnan, A., Bille, 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., Portner, 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 an-
thropogenic CO₂ emissions scenarios, Science, 349, 6243,
doi:10.1126/science.aac4722, 2015.

Gehlen, M., Séférian, R., Jones, D. O. B., Roy, T., Roth, R., Barry, J., Bopp, L., Doney, S. C., Dunne, J. P., Heinze, C., Joos, F., Orr, J. C., Resplandy, L., Segschneider, J., and Tjiputra, J.: Projected pH reductions by 2100 might put deep North Atlantic biodiversity at risk, Biogeosciences, 11, 6955–6967, doi:10.5194/bg-11-6955-2014, 2014.

Gerber, M. and Joos, F.: Carbon sources and sinks from an Ensem-
ble Kalman Filter ocean data assimilation, Global Biogeochem.
Cy., 24, GB3004, doi:10.1029/2009GB003531, 2010.

Gnanadesikan, A.: Oceanic ventilation and biogeochemical cycling: Understanding the physical mechanisms that produce realistic distributions of tracers and productivity, Global Biogeochem. Cy., 18, GB4010, doi:10.1029/2003GB002097, 2004.

Gruber, N.: Warming up, turning sour, losing breath: ocean biogeo-
chemistry under global change, Philos. T. R. Soc. A, 369, 1980–
1996, doi:10.1098/rsta.2011.0003, 2011.

Gruber, N. and Galloway, J. N.: An Earth-system perspec-
tive of the global nitrogen cycle, Nature, 451, 293–296,
doi:10.1038/nature06592, 2008.

Hajima, T., Kawamiya, M., Watanabe, M., Kato, E., Tachiiri, K.,
Sugiyama, M., Watanabe, S., Okajima, H., and Ito, A.: Modeling
in Earth system science up to and beyond IPCC AR5, Progress in
Earth and Planetary Science, 1, 29–25, doi:10.1186/s40645-014-
0029-y, 2014.

Heinze, C., Maier-Reimer, E., Winguth, A., and Archer, D.: A
global oceanic sediment model for long-term climate studies,
Global Biogeochem. Cy., 13, 221–250, 1999.

Heinze, M. and Ilyaina, T.: Ocean biogeochemistry in the
warm climate of the late Paleocene, Clim. Past, 11, 63–79,
doi:10.5194/cp-11-63-2015, 2015.

Henson, S. A., Sarmiento, J. L., Dunne, J. P., Bopp, L., Lima, I.,
Doney, S. C., John, J., and Beaulieu, C.: Detection of anthro-
pogenic climate change in satellite records of ocean chlorophyll
and productivity, Biogeosciences, 7, 621–640, doi:10.5194/bg-7-
621-2010, 2010.

Hobbs, W., Palmer, M. D., and Monselesan, D.: An Energy Con-
servation Analysis of Ocean Drift in the CMIP5 Global Cou-
pled Models, J. Climate, 29, 1639–1653, doi:10.1175/JCLI-D-
15-0477.1, 2015.

Hourdin, F., Musat, I., Bony, S., Braconnot, P., Codron, F.,
Dufresne, J.-L., Fairhead, L., Filiberti, M.-A., Friedlingstein, P.,
Grandpeix, J.-Y., Krinner, G., LeVan, P., Li, Z.-X., and Lott, F.:
The LMDZ4 general circulation model: climate performance and
sensitivity to parametrized physics with emphasis on tropical
convection, Clim. Dynam., 27, 787–813, doi:10.1007/s00382-
006-0158-0, 2006.

Ilyina, T., Six, K. D., Segschneider, J., Maier-Reimer, E., Li,
H., and Núñez-Riboni, I.: Global ocean biogeochemistry model
HAMOCC: Model architecture and performance as component
of the MPI-Earth system model in different CMIP5 experi-
mental realizations, J. Adv. Model. Earth Syst., 5, 287–315,
doi:10.1029/2012MS000178, 2013.

IPCC: Climate Change 2013: The Physical Science Basis, edited
by: Stoker, T. F., Qin, D., Plattner, G., Tignor, M., Allen, S. K.,
Boschung, J., Nauels, A., Xia, Y., Bex, V., and Midgley, P. M.,
Cambridge Univ. Press, Cambridge, UK, and New York, NY,
USA, 2013.

Ito, T. and Deutsch, C.: Variability of the Oxygen Minimum Zone in the Tropical North Pacific during the Late 20th Century, Global Biogeochem. Cy., 27, 1119–1128, doi:10.1002/2013GB004567, 2013.

Ito, T., Woloszyn, M., and Mazloff, M.: Anthropogenic carbon dioxide transport in the Southern Ocean driven by Ekman flow, Nature, 463, 80–83, doi:10.1038/nature08687, 2010.

Jickells, T. and Spokes, L.: The biogeochemistry of iron in seawater, edited by: Turner, D. R. and Hunter, K. A., John Wiley, Hoboken, NJ, USA, 2001.

Johnson, K., Chavez, F., and Friederich, G.: Continental-shelf sedi-
ment as a primary source of iron for coastal phytoplankton, Na-
ture, 398, 697–700, 1999.

Keeling, R. F., Körtzinger, A., and Gruber, N.: Ocean Deoxygenation in a Warming World, Annu. Rev. Marine. Sci., 2, 199–229, doi:10.1146/annurev.marine.010908.163855, 2009.

Keenlyside, N. S., Latif, M., Jungclaus, J., Kornblueh, L.,
and Roeckner, E.: Advancing decadal-scale climate pre-
diction in the North Atlantic sector, Nature, 453, 84–88,
doi:10.1038/nature06921, 2008.

Keller, K. M., Joos, F., and Raible, C. C.: Time of emergence of trends in ocean biogeochemistry, Biogeosciences, 11, 3647–3659, doi:10.5194/bg-11-3647-2014, 2014.

Key, R., Kozyr, A., Sabine, C., Lee, K., Wanninkhof, R., Bullis-
ter, J., Feely, R., Millero, F., Mordy, C., and Peng, T.: A global
ocean carbon climatology: Results from Global Data Analy-
sis Project (GLODAP), Global Biogeochem. Cy., 18, GB4031,
doi:10.1029/2004GB002247, 2004.

Khatiwala, S.: Fast spin up of ocean biogeochemical models using
matrix-free Newton-Krylov, Ocean Model., 23, 121–129, 2008.

Khatiwala, S., Visbeck, M., and Cane, M. A.: Accelerated sim-
ulation of passive tracers in ocean circulation models, Ocean
Model., 9, 51–69, doi:10.1016/j.ocemod.2004.04.002, 2005.

Kim, H.-M., Webster, P. J., and Curry, J. A.: Evaluation
of short-term climate change prediction in multi-model
CMIP5 decadal hindcasts, Geophys. Res. Lett., 39, L10701,
doi:10.1029/2012GL051644, 2012.
Knutti, R., Masson, D., and Gettelman, A.: Climate model genealog: Generation CMIP5 and how we got there, Geophys. Res. Lett., 40, 1194–1199, doi:10.1002/grl.50256, 2013.

Koven, C. D., Chambers, J. Q., Georgiou, K., Knox, R., Negron-Juarez, R., Riley, W. J., Arora, V. K., Brovkin, V., Friedlingstein, P., and Jones, C. D.: Controls on terrestrial carbon feedbacks by productivity versus turnover in the CMIP5 Earth System Models, Biogeosciences, 12, 5211–5228, doi:10.5194/bg-12-5211-2015, 2015.

Kriest, I. and Oschlies, A.: MOPS-1.0: towards a model for the regulation of the global oceanic nitrogen budget by marine biogeochemical processes, Geosci. Model Dev., 8, 2929–2957, doi:10.5194/gmd-8-2929-2015, 2015.

Krinner, G., Viovy, N., de Noblet-Ducoudré, N., Ogée, J., Polcher, J., Friedlingstein, P., Ciais, P., Sitch, S., and Prentice, I. C.: A dynamic global vegetation model for studies of the coupled atmosphere-biosphere system, Global Biogeochem. Cy., 19, 1–33, 2005.

Laufkötter, C., Vogt, M., Gruber, N., Aita-Noguchi, M., Aumont, O., Bopp, L., Buitenhuis, E., Doney, S. C., Dunne, J., Hashioka, T., Hauck, J., Hirata, T., John, J., Le Quéré, C., Lima, I. D., Nakano, H., Seferian, R., Totterdell, I., Vichi, M., and Völker, C.: Drivers and uncertainties of future global marine primary production in marine ecosystem models, Biogeosciences, 12, 6955–6984, doi:10.5194/bg-12-6955-2015, 2015.

Lehodey, P., Alheit, J., Barange, M., Baumgartner, T., Beaugrand, G., Drinkwater, K., Fromentin, J.-M., Hare, S. R., Ottersen, G., Perry, R. I., Roy, C., van der Lingen, C. D., and Werner, F.: Climate variability, fish, and fisheries, J. Climate, 19, 5009–5030, doi:10.1175/JCLI3898.1, 2006.

Le Quéré, C., Moriarty, R., Andrew, R. M., Peters, G. P., Ciais, P., Friedlingstein, P., Jones, S. D., Sitch, S., Tans, P., Arneth, A., Boden, T. A., Bopp, L., Bozec, Y., Canadell, J. G., Chini, L. P., Chevallier, F., Cosca, C. E., Harris, I., Hoppema, M., Houghton, R. A., House, J. I., Jain, A. K., Johannessen, T., Kato, E., Keeling, R. F., Kitidis, V., Klein Goldewijk, K., Koven, C., Landa, C. S., Landschützer, P., Lenton, A., Lima, I. D., Marland, G., Mathis, J. T., Metzl, N., Nojiri, Y., Olsen, A., Ono, T., Peng, S., Peters, W., Pfeil, B., Poulter, B., Raupach, M. R., Regnier, P., Rödenbeck, C., Saito, S., Salisbury, J. E., Schuster, U., Schwinger, J., Séférian, R., Segschneider, J., Steinhoff, T., Stocker, B. D., Sutton, A. J., Takahashi, T., Tilbrook, B., van der Werf, G. R., Viovy, N., Wang, Y.-P., Wanninkhof, R., Wiltshire, A., and Zeng, N.: Global carbon budget 2014, Earth Syst. Sci. Data, 7, 47–85, doi:10.5194/essd-7-47-2015, 2015.

Levitus, S. and Boyer, T.: World ocean atlas 1994, volume 4: Temperature, PB-95-270112/XAB, National Environmental Satellite, Data, and Information Service, Washington, DC, USA, 1994.

Levitus, S., Conkright, M. E., Reid, J. L., Najjar, R. G., and Mantyla, A.: Distribution of nitrate, phosphate and silicate in the world oceans, Prog. Oceanogr., 31, 245–273, 1993.

Levitus, S., Antonov, J. I., Baranova, O. K., Boyer, T. P., Coleman, C. L., Garcia, H. E., Grodsky, A. I., Johnson, D. R.,LOCarnini, R. A., Mishonov, A. V., Reagan, J. R., Sazama, C. L., Seidov, D., Smolyar, I., Yarosh, E. S., and Zweng, M. M.: The World Ocean Database TI, Data Science Journal, 12, WDS229–WDS234, 2013.

Lévy, M., Lengaigne, M., Bopp, L., Vincent, E. M., Madec, G., Ethe, C., Kumar, D., and Sarma, V. V. S. S.: Contribution of tropical cyclones to the air-sea CO₂ flux: A global view, Global Biogeochem. Cy., 26, GB2001, doi:10.1029/2011GB004145, 2012.

Lindsay, K., Bonan, G. B., Doney, S. C., Hoffman, F. M., Lawrence, D. M., Long, M. C., Mahowald, N. M., Moore, J. K., Randerson, J. T., and Thornton, P. E.: Preindustrial Control and 20th Century Carbon Cycle Experiments with the Earth System Model CESM1(BGC), J. Climate, 27, 8981–9005, doi:10.1175/JCLI-D-12-00565.1, 2014.

Ludwig, W., Probst, J., and Kempe, S.: Predicting the oceanic input of organic carbon by continental erosion, Global Biogeochem. Cy., 10, 23–41, 1996.

Madec, G.: NEMO ocean engine, Institut Pierre-Simon Laplace (IPSL), France, available at: http://www.nemo-ocean.eu/About-NEMO/Reference-manuals (last access: November 2013), 2008.

Maier-Reimer, E.: Geochemical cycles in an ocean general circulation model. Preindustrial tracer distributions, Global Biogeochem. Cy., 7, 645–677, doi:10.1029/93GB01355, 1993.

Maier-Reimer, E. and Hasselmann, K.: Transport and storage of CO₂ in the ocean – an inorganic ocean-circulation carbon cycle model, Clim. Dynam., 2, 63–90, doi:10.1007/BF01054491, 1987.

Marinov, I., Gnanadesikan, A., Sarmiento, J. L., Toggweiler, J. R., Follows, M., and Mignone, B. K.: Impact of oceanic circulation on biological carbon storage in the ocean and atmospheric pCO₂, Global Biogeochem. Cy., 22, GB3007, doi:10.1029/2007GB002958, 2008.

Massonnet, F., Fichefet, T., Goosse, H., Bitz, C. M., Philippon-Berthier, G., Holland, M. M., and Barriat, P.-Y.: Constraining projections of summer Arctic sea ice, The Cryosphere, 6, 1383–1394, doi:10.5194/tc-6-1383-2012, 2012.

Matei, D., Baehr, J., Jungclaus, J. H., Haak, H., Muller, W. A., and Marotzke, J.: Multiyear Prediction of Monthly Mean Atlantic Meridional Overturning Circulation at 26.5°N, Science, 335, 76–79, doi:10.1126/science.1210299, 2012.

Mawji, E., Schlitzer, R., Dodas, E. M., et al.: The GEOTRACES Intermediate Data Product 2014, Marine Chemistry, 177, Part 1, 1–8, doi:10.1016/j.marchem.2015.04.005.

Meehl, G. A., Goddard, L., Murphy, J., Stouffer, R. J., Boer, G., Danabasoglu, G., Dixon, K., Giorgetta, M. A., Greene, A. M., Hawkins, E., Hegerl, G., Karoly, D., Keenlyside, N., Kimoto, M., Kirtman, B., Navarra, A., Pulwarty, R., Smith, D., Stammer, D., and Stockdale, T.: Decadal Prediction, B. Am. Meteorol. Soc., 90, 1467–1485, doi:10.1175/2009BAMS2778.1, 2009.

Meehl, G. A., Goddard, L., Boer, G., Burgman, R., Branstator, G., Cassou, C., Corti, S., Danabasoglu, G., Doblas-Reyes, F., Hawkins, E., Karspeck, A., Kimoto, M., Kumar, A., Matei, D., Mignot, J., Msadek, R., Pohlmann, H., Rienecker, M., Rosati, T., Schneider, E., Smith, D., Sutton, R., Teng, H., van Oldenborgh, G. J., Vecchi, G., and Yeager, S.: Decadal Climate Prediction: An Update from the Trenches, B. Am. Meteorol. Soc., 95, 243–267, doi:10.1175/BAMS-D-12-00241.1, 2013.

Meehl, G. A., Moss, R., Taylor, K. E., Eyring, V., Stouffer, R. J., Bony, S., and Stevens, B.: Climate Model Intercomparisons: Preparing for the Next Phase, Eos Trans. AGU, 95, 77–78, doi:10.1002/2014EO090001, 2014.
Mignot, J., Swingedouw, D., Deshayes, J., Marti, O., Talandier, C., Séférían, R., Lengaigne, M., and Madec, G.: On the evolution of the oceanic component of the IPSL climate models from CMIP3 to CMIP5: A mean state comparison, Ocean Model., 72, 167–184, 2013.

Mikaloff Fletcher, S. E., Gruber, N., Jacobson, A. R., Gloor, M., Doney, S. C., Dutkiewicz, S., Gerber, M., Follows, M., Joos, F., Lindsay, K., Menemenlis, D., Mouchet, A., Müller, S. A., and Sarmiento, J. L.: Inverse estimates of the oceanic sources and sinks of natural CO₂ and the implied oceanic carbon transport, Global Biogeochem. Cy., 21, GB1010, doi:10.1029/2006GB002751, 2007.

Moore, J., Doney, S., Kleypas, J., Glover, D., and Fung, I.: An intermediate complexity marine ecosystem model for the global domain, Deep-Sea Res. Pt. II, 49, 403–462, 2002.

Moore, J., Doney, S., and Lindsay, K.: Upper ocean ecosystem dynamics and iron cycling in a global three-dimensional model, Global Biogeochem. Cy., 18, GB4028, doi:10.1029/2004GB002220, 2004.

Orr, J. C.: Global Ocean Storage of Anthropogenic Carbon (GOSAC), final report, EC Environ. and Clim. Programme, Inst. Pierre Simon Laplace, Paris, France, 117 pp., 2002.

Phillips, T. J., Potter, G. L., Williamson, D. L., Cederwall, R. T., Boyle, J. S., Fiorino, M., Hnilo, J. J., Olson, J. G., Xie, S., and Yio, J. J.: Evaluating Parameterizations in General Circulation Models: Climate Simulation Meets Weather Prediction, B. Am. Meteorol. Soc., 85, 1903–1915, doi:10.1175/BAMS-85-12-1903, 2004.

Resplandy, L., Bopp, L., Orr, J. C., and Dunne, J. P.: Role of mode and intermediate waters in future ocean acidification: Analysis of CMIP5 models, Geophys. Res. Lett., 40, 3091–3095, 2013.

Resplandy, L., Séférían, R., and Bopp, L.: Natural variability of CO₂ and O₂ fluxes: What can we learn from centuries-long climate models simulations? J. Geophys. Res.-Oceans, 120, 384–404, doi:10.1002/2014JC010463, 2015.

Rodgers, K. B., Lin, J., and Frölicher, T. L.: Emergence of multiple ocean ecosystem drivers in a large ensemble suite with an Earth system model, Biogeosciences, 12, 3301–3320, doi:10.5194/bg-12-3301-2015, 2015.

Romanou, A., Gregg, W. W., Romanski, J., Kelley, M., Bleck, R., Healy, R., Nazarenko, L., Russell, G., Schmidt, G. A., Sun, S., and Tausnev, N.: Natural air-sea flux of CO₂ in simulations of the NASA-GISS climate model: Sensitivity to the physical ocean model formulation, Ocean Model., 66, 26–44, doi:10.1016/j.ocemod.2013.01.008, 2013.

Romanou, A., Romanski, J., and Gregg, W. W.: Natural ocean carbon cycle sensitivity to parameterizations of the recycling in a climate model, Biogeosciences, 11, 1137–1154, doi:10.5194/bg-11-1137-2014, 2014.

Rose, K. A., Roth, B. M., and Smith, E. P.: Skill assessment of spatial maps for oceanographic modeling, J. Marine Syst., 76, 34–48, doi:10.1016/j.jmarsys.2008.05.013, 2009.

Roy, T., Bopp, L., Gehlen, M., Schneider, B., Cadule, P., Frölicher, T. L., Segschneider, J., Tjiputra, J., Heinze, C., and Joos, F.: Regional Impacts of Climate Change and Atmospheric CO₂ on Future Ocean Carbon Uptake: A Multi-model Linear Feedback Analysis, J. Climate, 24, 2300–2318, doi:10.1175/2010JCLI3787.1, 2011.

Sarmiento, J. L. and Gruber, N.: Ocean Biogeochemical Dynamics, Princeton University Press, Princeton, New Jersey, USA, 526 pp., 2006.

Schwinger, J., Tjiputra, J. F., Heinze, C., Bopp, L., Christian, J. R., Gehlen, M., Ilyina, T., Jones, C. D., Salas-Mélia, D., Segschneider, J., Séférían, R., and Totterdell, I.: Nonlinearity of Ocean Carbon Cycle Feedbacks in CMIP5 Earth System Models, J. Climate, 27, 3869–3888, doi:10.1175/JCLI-D-13-00452.1, 2014.

Séférían, R., Iudicone, D., Bopp, L., Roy, T., and Madec, G.: Water Mass Analysis of Effect of Climate Change on Air–Sea CO₂ Fluxes: The Southern Ocean, J. Climate, 25, 3894–3908, doi:10.1175/JCLI-D-11-00291.1, 2012.

Séférían, R., Bopp, L., Gehlen, M., Orr, J., Ethé, C., Cadule, P., Aumont, O., Salas y Mélia, D., Voldoire, A., and Madec, G.: Skill assessment of three earth system models with common marine biogeochemistry, Clim. Dynam., 40, 2549–2573, doi:10.1007/s00382-012-1362-8, 2013.

Séférían, R., Ribes, A., and Bopp, L.: Detecting the anthropogenic influences on recent changes in ocean carbon uptake, Geophys. Res. Lett., 41, 5968–5977, doi:10.1002/2014GL061223, 2014.

Séférían, R., Delire, C., Decharme, B., Voldoire, A., Salas y Melia, D., Chevallier, M., Saint-Martin, D., Aumont, O., Calvet, J.-C., Carrer, D., Douville, H., Franchistéguy, L., Joetzjer, E., and Sénési, S.: Development and evaluation of CNRM Earth system model – CNRM-ESM1, Geosci. Model Dev., 9, 1423–1453, doi:10.5194/gmd-9-1423-2016, 2016.

Sen Gupta, A. S., Muir, L. C., Brown, J. N., Phipps, S. J., Durack, P. J., Monselesan, D., and Wijffels, S. E.: Climate Drift in the CMIP3 Models, J. Climate, 25, 4621–4640, doi:10.1175/JCLI-D-11-00312.1, 2012.

Sen Gupta, A. S., Jourdain, N. C., Brown, J. N., and Monselesan, D.: Climate Drift in the CMIP5 models, J. Climate, 26, 8597–8615, doi:10.1175/JCLI-D-12-00521.1, 2013.

Servonnat, J., Mignot, J., Guilyardi, E., Swingedouw, D., Séférían, R., and Labetoulle, S.: Reconstructing the subsurface ocean decadal variability using surface nudging in a perfect model framework, Clim. Dynam., 44, 1–24, doi:10.1007/s00382-014-2184-7, 2014.

Smith, D. M., Cusack, S., Colman, A. W., Folland, C. K., Harris, G. R., and Murphy, J. M.: Improved Surface Temperature Prediction for the Coming Decade from a Global Climate Model, Science, 317, 796–799, doi:10.1126/science.1139540, 2007.

Smith, M. J., Palmer, P. I., Purves, D. W., Vanderwel, M. C., Lyutsarev, V., Calderhead, B., Joppa, L. N., Bishop, C. M., and Emmott, S.: Changing how Earth System Modelling is done to provide more useful information for decision making, science and society, B. Am. Meteorol. Soc., 95, 1453–1464, doi:10.1175/BAMS-D-13-00080.1, 2014.

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, doi:10.5194/bg-7-979-2010, 2010.

Stouffer, R. J., Weaver, A. J., and Eby, M.: A method for obtaining pre-twentieth century initial conditions for use in climate change studies, Clim. Dynam., 23, 327–339, doi:10.1007/s00382-004-0446-5, 2004.

Stow, C. A., Jolliff, J., McGillicuddy, D. J. J., Doney, S. C., Allen, J. I., Friedrichs, M. A. M., Rose, K. A., and
Wallhead, P.: Skill assessment for coupled biological/physical models of marine systems, J. Marine Syst., 76, 4–15, doi:10.1016/j.jmarsys.2008.03.011, 2009.

Swingedouw, D., Mignot, J., Labetoulle, S., Guilyardi, E., and Madec, G.: Initialisation and predictability of the AMOC over the last 50 years in a climate model, Clim. Dynam., 40, 2381–2399, doi:10.1007/s00382-012-1516-8, 2013.

Tagliabue, A. and Völker, C.: Towards accounting for dissolved iron speciation in global ocean models, Biogeosciences, 8, 3025–3039, doi:10.5194/bg-8-3025-2011, 2011.

Tagliabue, A., Bopp, L., and Gehlen, M.: The response of marine carbon and nutrient cycles to ocean acidification: Large uncertainties related to phytoplankton physiological assumptions, Global Biogeochem. Cy., 25, GB3017, doi:10.1029/2010GB003929, 2011.

Takahashi, T., Broecker, W., and Langer, S.: Redfield Ratio Based on Chemical-Data From Isopycnal Surfaces, J. Geophys. Res.-Oceans, 90, 6907–6924, 1985.

Tanhua, T., Koertzinger, A., Friis, K., Waugh, D. W., and Wallace, D. W. R.: An estimate of anthropogenic CO₂ inventory from decadal changes in oceanic carbon content, P. Natl. Acad. Sci. USA, 104, 3037–3042, doi:10.1073/pnas.0606574104, 2007.

Tegen, I. and Fung, I.: Contribution to the Atmospheric Mineral Aerosol Load From Land-Surface Modification, J. Geophys. Res.-Atmos., 100, 18707–18726, 1995.

Tjiputra, J. F., Roelandt, C., Bentsen, M., Lawrence, D. M., Lorentzen, T., Schwinger, J., Seland, Ø., and Heinze, C.: Evaluation of the carbon cycle components in the Norwegian Earth System Model (NorESM), Geosci. Model Dev., 6, 301–325, doi:10.5194/gmd-6-301-2013, 2013.

Tjiputra, J. F., Olsen, A., Bopp, L., Lenton, A., Pfeil, B., Roy, T., Segschneider, J., Totterdell, I., and Heinze, C.: Long-term surface pCO₂ trends from observations and models, Tellus B, 66, 151–168, doi:10.1007/s00382-007-0342-x, 2014.

Vancoppenolle, M., Bopp, L., Madec, G., Dunne, J. P., Ilyina, T., Halloran, P. R., and Steiner, N.: Future Arctic Ocean primary productivity from CMIP5 simulations: Uncertain outcome, but consistent mechanisms, Global Biogeochem. Cy., 27, 605–619, 2013.

Vichi, M., Manzini, E., Fogli, P. G., Alessandri, A., Patara, L., Scoccimarro, E., Masina, S., and Navarra, A.: Global and regional ocean carbon uptake and climate change: sensitivity to a substantial mitigation scenario, Clim. Dynam., 37, 1929–1947, doi:10.1007/s00382-011-1079-0, 2011.

Volodin, E. M., Dianskii, N. A., and Gusev, A. V.: Simulating present-day climate with the INMCM4.0 coupled model of the atmospheric and oceanic general circulations, Izv. Atmos. Ocean. Phy. +, 46, 414–431, doi:10.1134/S000143381004002X, 2010.

Walin, G., Hieronymus, J., and Nycander, J.: Source-related variables for the description of the oceanic carbon system, Geochem. Geophys. Geosyst., 15, 3675–3687, doi:10.1002/2014GC005383, 2014.

Wanninkhof, R.: A relationship between wind speed and gas exchange over the ocean, J. Geophys. Res., 97, 7373–7382, 1992.

Wassmann, P., Duarte, C. M., Agustí, S., and Sejr, M. K.: Footprints of climate change in the Arctic marine ecosystem, Glob. Change Biol., 17, 1235–1249, doi:10.1111/j.1365-2486.2010.02311.x, 2010.

Watanabe, S., Hajima, T., Sudo, K., Nagashima, T., Takemura, T., Okajima, H., Nozawa, T., Kawase, H., Abe, M., Yokohata, T., Ise, T., Sato, H., Kato, E., Takata, K., Emori, S., and Kawamiya, M.: MIROC-ESM 2010: model description and basic results of CMIP5-20c3m experiments, Geosci. Model Dev., 4, 845–872, doi:10.5194/gmd-4-845-2011, 2011.

Wenzel, S., Cox, P. M., Eyring, V., and Friedlingstein, P.: Emergent constraints on climate-carbon cycle feedbacks in the CMIP5 Earth system models, J. Geophys. Res.-Biogeosci., 119, 794–807, doi:10.1002/2013JG002591, 2014.

Wu, T., Li, W., Ji, J., Xin, X., Li, L., Wang, Z., Zhang, Y., Li, J., Zhang, F., Wei, M., Shi, X., Wu, F., Zhang, L., Chu, M., Jie, W., Liu, Y., Wang, F., Liu, X., Li, Q., Dong, M., Liang, X., Gao, Y., and Zhang, J.: Global carbon budgets simulated by the Beijing Climate Center Climate System Model for the last century, J. Geophys. Res.-Atmos., 118, 4326–4347, doi:10.1002/jgrd.50320, 2013.

Wunsch, C. and Heimbach, P.: Practical global oceanic state estimation, Physica D, 230, 197–208, doi:10.1016/j.physd.2006.09.040, 2007.

Wunsch, C. and Heimbach, P.: How long to oceanic tracer and proxy equilibrium?, Quaternary Sci. Rev., 27, 637–651, doi:10.1016/j.quascirev.2008.01.006, 2008.

Yool, A., Oschlies, A., Nurser, A. J. G., and Gruber, N.: A model-based assessment of the TrOCA approach for estimating anthropogenic carbon in the ocean, Biogeosciences, 7, 723–751, doi:10.5194/bg-7-723-2010, 2010.

Yool, A., Popova, E. E., and Anderson, T. R.: MEDUSA-2.0: an intermediate complexity biogeochemical model of the marine carbon cycle for climate change and ocean acidification studies, Geosci. Model Dev., 6, 1767–1811, doi:10.5194/gmd-6-1767-2013, 2013.

Zeebe, R. E. and Wolf-Gladrow, D. A.: CO₂ in seawater: equilibrium, kinetics, isotopes, Oceanography Book Series, vol. 65, Elsevier, Amsterdam, the Netherlands, 346 pp., 2001.