The Relative Importance of Phytoplankton Light Absorption and Ecosystem Complexity in an Earth System Model

We investigate the relative importance of ecosystem complexity and phytoplankton light absorption for climate studies. While the complexity of Earth System models (ESMs) with respect to marine biota has increased over the past years, the relative importance of biological processes in driving climate‐relevant mechanisms such as the biological carbon pump and phytoplankton light absorption is still unknown. The climate effects of these mechanisms have been studied separately, but not together. To shed light on the role of biologically mediated feedbacks, we performed different model experiments with the EcoGENIE ESM. The model experiments have been conducted with and without phytoplankton light absorption and with two or 12 plankton functional types. For a robust comparison, all simulations are tuned to have the same primary production. Our model experiments show that phytoplankton light absorption changes ocean physics and biogeochemistry. Higher sea surface temperature decreases the solubility of CO2 which in turn increases the atmospheric CO2 concentration, and finally the atmospheric temperature rises by 0.45°C. An increase in ecosystem complexity increases the export production of particulate organic carbon but decreases the amount of dissolved organic matter. These changes in the marine carbon cycling, however, hardly reduces the atmospheric CO2 concentrations and slightly decreases the atmospheric temperature by 0.034°C. Overall we show that phytoplankton light absorption has a higher impact on the carbon cycle and on the climate system than a more detailed representation of the marine biota.

increases the sea surface temperature (SST) while a higher ecosystem complexity leads to a slightly deeper downward flux of organic matter • Our results suggest that phytoplankton light absorption has a higher impact on the climate system than increasing ecosystem complexity

Supporting Information:
Supporting Information may be found in the online version of this article.et al., 2005) with a focus on the carbon cycle.Indeed marine biota play an essential role in the carbon cycle through the biological pump, defined as the uptake of carbon dioxide at the surface of the ocean and the sinking of the organic carbon to the abyssal ocean.Over geological time the biological pump has shaped ocean chemistry, biogeochemical cycling and ecosystem structure (Meyer et al., 2016).Additionally the biological pump has contributed to past variations of atmospheric CO 2 , influencing the glacial/interglacial episodes during the ice ages of the Pleistocene period (Turner, 2015).Watson and Liss (1998) used a simple model of Sarmiento & Toggweiler, (1984) to determine the importance of the biological pump on atmospheric CO 2 .Their study suggests that if all marine life were to die, atmospheric CO 2 would increase by ca 450 ppm after a few hundred years.
Marine biota have not only an impact on the climate system through the carbon cycle but also affect, for example, the ocean's thermal, optical, and mechanical properties via biogeophysical mechanisms (Hense et al., 2013).Among these biogeophysical mechanisms, phytoplankton light absorption is particularly important.Previous observations between 1972 and 2006 reveal that during spring, when large blooms occur in coastal mid-latitude regions, 31%-42% of the light can be absorbed by the phytoplankton (Fleming- Lehtinen & Laamanen, 2012).The heat distribution in the upper ocean is then altered (Lewis et al., 1990;Sonntag & Hense, 2011), changing the sea surface temperature (SST) (Kahru et al., 1993;Patara et al., 2012), an important climate variable to understand interactions between the ocean and the atmosphere.Indeed, the SST affects atmospheric temperature, atmospheric humidity content, precipitation, as well as heat transfer between the ocean and the atmosphere (Jang et al., 2016;Lim et al., 2016).These changes in atmospheric physics and chemistry can even alter the Walker and Hadley circulation (Gnanadesikan & Anderson, 2009;Paulsen et al., 2018).A recent study (Patara et al., 2012) shows that phytoplankton light absorption increases the SST, enhancing evaporation, and atmospheric humidity, and therefore the greenhouse effect.Changes in SST, in turn impact the oceanic circulation (Manizza et al., 2008).These changes in ocean physics will cause feedbacks on the climate system.For instance, Manizza et al. (2005) used an Ocean General Circulation Model (OGCM) to show that phytoplankton light absorption intensifies the seasonal cycle of temperature, mixed layer depth, and ice cover by roughly 10%, leading to an increase in phytoplankton biomass and thus to an amplification of the initial physical perturbations.But these feedbacks depend on the spatial and seasonal scale (Oschlies, 2004).Finally, previous model studies focusing on phytoplankton light absorption report a strengthening (Marzeion et al., 2005;Paulsen et al., 2018) or weakening (Anderson et al., 2007;Jochum et al., 2010) of El Niño-Southern Oscillation as well as changes in its oscillation periods (Zhang et al., 2009).
None of these model studies have quantified the climate response to biologically-driven feedback mechanisms combined.In particular, the relative importance of ecosystem complexity and the effect on carbon export production compared to SST changes induced by phytoplankton light absorption has not been investigated.Additionally the climate response has been only evaluated for a specific climate.The potential of marine biota in altering the climate system, in which atmospheric CO 2 is allowed to evolve freely, has not been demonstrated.To investigate these aspects we use an ESM of intermediate complexity (EMIC) (Claussen et al., 2002), EcoGENIE (Ward et al., 2018).For this purpose, we have modified the model by implementing phytoplankton light absorption.

Model Description
We used the GENIE model (Lenton et al., 2007), consisting of several modules describing the dynamics of the individual Earth system components (Figure 1).GENIE is widely used for past climate and carbon-cycle studies (Gutjahr et al., 2017;Ödalen et al., 2018).The new ecosystem component (ECOGEM) and GENIE form the recent EcoGENIE model (Ward et al., 2018).We use the simplified atmosphere and carbon-centric version (cGENIE) that has been previously applied to analyze interactions between marine biological productivity, biogeochemistry, and climate (Gibbs et al., 2016;Meyer et al., 2016;Tagliabue et al., 2016).Eco-GENIE contains different Earth system components and related processes including ocean circulation and marine biogeochemistry, atmospheric circulation, the marine ecosystem component, and sea-ice dynamics.We modified the ecosystem (ECOGEM) and ocean component (GOLDSTEIN) to account for phytoplankton light absorption.

Ocean Physics Component
The GOLDSTEIN component describes the ocean physics.It is a 3-D frictional-geostrophic ocean model (Edwards & Marsh, 2005;Marsh et al., 2011) based on the reduced physics of the thermocline equations, described for a single-basin configuration in Edwards et al. (1998).This model is similar to general circulation models, except that it neglects momentum advection and acceleration but incorporates eddy-induced and isopycnal mixing.This oceanic component includes a surface mixed layer scheme based on the scheme by Kraus and Turner (1967).The parameters for GOLDSTEIN are calibrated against annual mean climatological observations of temperature, salinity, surface air temperature, and humidity using the ensemble Kalman filter (EnKF) methodology (Annan et al., 2005;Hargreaves et al., 2004).

Sea-Ice Component
The GOLDSTEINSEAICE component describes the sea-ice dynamics.It is a 2-D model and solves the equations for part of the ocean surface covered by ice (Edwards & Marsh, 2005).A diagnostic equation is solved for the ice surface temperature.The sea-ice growth or decay depends on the net heat flux into the ice.The sea-ice dynamic consists of advection by surface currents and diffusion.The sea-ice component acts as a coupling module between the ocean and the atmosphere.

Atmospheric Component
The Energy and Moisture Balance Model (EMBM) describes the atmospheric dynamics (Edwards & Marsh, 2005).EMBM is based on the UVic ESM (Weaver et al., 2001).It is a vertically integrated 2-D atmospheric model with air temperature and specific humidity as prognostic variables.Heat and moisture are advected by winds and mixing.Precipitation instantaneously removes all moisture corresponding to the excess above a relative humidity.Atmosphere, ocean, and sea-ice exchange moisture and heat fluxes.

Ocean Biogeochemistry Component
The BIOGEM module considers marine biogeochemical processes (Ridgwell et al., 2007).This module calculates the transformations and spatial redistribution of biogeochemical quantities, plus the sea-air gas exchange of CO 2 and O 2 .In this model, the state variables are the inorganic resources and the dead organic matter.The biological pump is parametrized by an implicit biological community: Biological uptake is limited by light, temperature, and nutrient availability.Any uptake returns instantly to exported particulate organic matter (POM) or dissolved organic matter (DOM) in the ocean interior.DOM is converted back to nutrients in the upper four oceanic layers (0-590 m) while POM turns into nutrients in the deeper layers (590 m-seafloor).Eventually all organic matter is remineralized back to nutrients.
The model includes iron (Fe) and phosphate (PO 4 ) as limiting nutrients which is sufficient to realistically describe the distribution pattern of phytoplankton and zooplankton, nitrate (NO 3 ) is not considered here.

Ecosystem Component
The ECOGEM component represents the marine plankton community and associated interactions within the ecosystem (Ward et al., 2018).The state variables are not subject to physical transport, there are only local sources and sinks.Biological uptake is limited by light, temperature, and nutrient availability replacing the BIOGEM formulation and ECOGEM also considers iron-light co-limitation.The ecosystem community is composed of different plankton populations, so-called PFTs.They are described by their size, their taxonomic position, and thus by a set of specific and size-dependent traits.We consider different degrees of ecosystem complexity where the community is divided in two classes of PFTs: phytoplankton and zooplankton.The number of phytoplankton or zooplankton is not fixed and can be subdivided into different size classes with size-dependent traits.The phytoplankton populations are characterized by nutrient uptake and photosynthesis traits whereas the zooplankton populations by predation traits.Moreover, all populations are subject to respiration, mortality, and internal trophic interactions.Plankton mortality and grazing are the two sources of organic matter, with a partitioning between non-sinking dissolved and sinking POM.According to Ward et al. (2018) the partitioning of organic matter is a size-based sigmoidal function following Equation 1: where, β is the parameter of partitioning between DOM and POM, ESD is the equivalent spherical diameter used to calculate the plankton cell volume, β a is the maximum fraction to DOM as ESD approaches zero, β b is the minimum fraction to DOM as ESD approaches infinity and β c is the size at which the partitioning is 50:50 between DOM and POM.Please note that for the simulations with different numbers of PFTs, the average ESD for the entire population is the same but β will still be different, because of the nonlinear equation.

Coupling Between BIOGEM and ECOGEM
The calculation in BIOGEM are performed 48 times per model year while the calculation in ECOGEM takes 20 time-steps for each BIOGEM time-step (Ward et al., 2018).At the beginning of each ECOGEM time-step the concentration of inorganic matter and important properties of the physical environment are imported from BIOGEM.The marine biota through photosynthesis transforms the inorganic compounds into living biomass.The rate of change in living biomass is used only to update the living biomass concentrations in ECOGEM (Ward et al., 2018).At the end of each ECOGEM time-step, the rates of change in inorganic and dead organic matter are passed back to BIOGEM and is used to update dissolved inorganic carbon (DIC), DOM, PO 4 , Fe, oxygen, and alkalinity.POM is instantly remineralized at depth using the standard export production functions.

Grid Resolution
The horizontal grid (36 × 36) is constructed to be uniform in longitude and uniform in sine latitude, giving ∼3.2° latitudinal increments at the equator increasing to 19.2° in the highest latitude.This horizontal mesh has been widely used for very large ensembles (Marsh et al., 2004) and biogeochemical simulations (Cameron et al., 2005) with focus on the carbon cycle (Colbourn, 2011).To better resolve the light absorption effect on ocean physics by biota, we use a higher vertical resolution than previous configurations with BIOGEM or EcoGENIE (e.g., Ward et al., 2018).We consider 32 vertical layers that increase logarithmically from 29.38 m for the surface layer to 456.56 m for the deepest layer.All physical and biological parameters in the model are unchanged from the tuning of Ward et al. (2018) with exceptions of β a in Equation 1 (see below for an explanation) and the Atlantic-Pacific moisture (APM) flux correction parameter.We increase the APM flux correction parameter from 0.8 to 1.53 Sv in order to simulate a realistic Atlantic meridional overturning circulation (AMOC) of 14.2 Sv.This flux correction is required because the atmospheric component is too simple to transport moisture across the American continent, and this parameter adjusts the salt balance between the different ocean basins.We need to validate the ocean circulation of our newly configured model against observations.We first focus on the AMOC because it is the main driver of the worldwide oceanic circulation.
The increase in vertical resolution improves the representation of the AMOC (Figure 2).The AMOC is closer to observations (based on World Ocean Circulation Experiment (WOCE), Lumpkin and Speer (2007)) and model results of Atmosphere-OGCMs (AOGCMs) with a higher spatial resolution (e.g., Boulton et al., 2014).In our new model setup the clockwise overturning cell in the upper ocean is more condensed and thus more realistic (Figures 2b and 2c).In addition, the model now reproduces a deeper counterclockwise overturning cell which was absent in the coarse resolution setup (Figure 2a).A change in AMOC and thus ocean circulation affects also the distribution and magnitude of the biogeochemical quantities.While in some regions, such as the subtropical ocean, the model representation is improved, there are other regions, for example, the Arctic Ocean with less agreement between the model and observations.Overall, the biogeochemical fields are different between the model runs with a coarse and a finer vertical resolution.Note that we are mainly interested in resolving the effects of light absorption and the relative differences between our selected experiments (see below); therefore the absolute values are less relevant.The comparison of the biogeochemical quantities between the coarse, the finer vertical resolution and observations are presented in the appendix (Supporting information S1).

Light Absorption in the Ocean
The incoming short-wave radiation is taken from the climate component of the model and varies seasonally (Edwards & Marsh, 2005;Marsh et al., 2011).The model takes into account the inhibition of light penetration due to the presence of organic and inorganic particles as well as dissolved molecules.The vertical light attenuation is described by the scheme Equation 2: where I(z) is the radiation at depth z, I 0 is the radiation at the surface of the ocean, k w is the light absorption by water (0.04 m −1 ), k Chl is the light absorption by chlorophyll (0.03 m −1 (mg Chl) −1 ) and Chl(z) is the chlorophyll concentration at depth z.The values for k w and k Chl are adopted from Ward et al. (2018).In our model, I 0 is always negative because it is a downward flux from the sun to the surface of the ocean.Solar radiation decreases exponentially with depth through attenuation; maximum absorption occurs in the upper ocean layer and the minimum absorption takes place in the lowest ocean layer.In the standard model (Ward et al., 2018), all of this solar energy is absorbed in the surface layer while we allow light to penetrate until the sixth oceanic layer (221.84 m deep).
Phytoplankton changes the optical properties of the ocean (Sonntag & Hense, 2011) through phytoplankton light absorption.The absorption process can cause a radiative heating and change in ocean temperature.
We implemented the warming of water by light absorption of phytoplankton into the model following, for example, Hense (2007) and Patara et al. ( 2012) Equation 3: ∂T/∂t denotes the temperature changes, c p is the specific heat capacity of water, ρ is the ocean density, I is the solar radiation incident at the ocean surface, and z is the depth.We assume that the whole light absorption leads to heating of the water.

Model Setup and Experiments
We performed four different model experiments (Table 1) to study the impact of varying ecosystem complexity and phytoplankton light absorption on the climate.We performed a BIOGEM spin-up for 10,000 years to allow for a realistic nutrient distribution in the ocean.The spin-up is run with a constant atmospheric CO 2 concentration of 278 ppm.The experiments restart after this period for 1,000 years with ECOGEM, meaning that all experiments consider zooplankton and phytoplankton.Furthermore, in the simulations the atmospheric CO 2 is not fixed and can evolve freely over time.After 700 years of simulation the model reaches a steady-state but we analyze the outputs after model runs of 1,000 years.The model setup we used is very similar to the existing model setup used by Ward et al. (2018) to describe ECOGEM, but with five differences between the model setups: First we use a different grid than Ward et al. (2018).They used the "worlg4″ topography, while we use the topography "ra32lv".In "worlg4″ the Torres strait between Australia and Papua New Guinea is open, permitting an oceanic connection between the Pacific and Indian ocean, while in "ra32lv" the Torres strait is closed.The main difference between these two configurations, however, is the vertical resolution of the oceanic component.Indeed, in this study the ocean has 32 layers with a surface layer of 29.38 m while Ward et al. (2018) have 16 oceanic layers with a surface layer of 80.84 m.For a comparison of the biogeochemical quantities between both configurations the reader is referred to the appendix (Supporting Information S1).
Second, Ward et al. (2018) used 17 oceanic tracers while we only use 14 oceanic tracers relevant for the climate; we remove SO 4 , H 2 S, Mg and the atmospheric tracer pH 2 S.  Third, we modify the ecosystem community to conduct our simulations.Our ecosystem community is based on the community used by Ward et al. (2018) to describe the ECOGEM model (Supporting Information S2).We run ECOGEM with different complexity of one and six phytoplankton and zooplankton species depending on the specific experiment.We additionally account for the phytoplankton light absorption for our selected simulations.
Fourth, we modify the vertical light attenuation scheme and apply Equation 2 and 3 for all the experiments.With this change, the absorption of the solar radiation can occur in all the layers of the ocean and not only in the uppermost layer as in the standard setup used by Ward et al. (2018).But for the experiments without phytoplankton light absorption k Chl = 0 m −1 (mg Chl) −1 , so the light is attenuated only by k w (see Equation 2).We assume a closed carbon cycle for all the model runs.Thus, there is no input of carbon through volcanic fluxes or anthropogenic activities and only the relative size of the reservoirs (atmospheric CO 2 , oceanic CO 2 , particulate organic carbon (POC), DIC…) can vary, not their total.
We therefore allow that different climates can develop.Depending on the strength of the ecosystem response in the respective experiment more or less CO 2 might be emitted into the atmosphere leading to altered air temperature.Thus, our setup allows changes in the atmospheric CO 2 concentrations and Earth's energy budget.
All experiments are forced with the same constant flux of dissolved iron into the ocean surface (Mahowald et al., 2006).The incoming shortwave radiation varies seasonally but no trend is considered.The longwave radiation emitted by the surface of the planet is absorbed by the atmosphere and re-emitted upward and downward (Weaver et al., 2001).The re-emission depends on the greenhouse gases concentration in the atmosphere.

Comparable State of the Experiments
In this study, we are interested in the relative importance of the processes regarding phytoplankton light absorption and ecosystem complexity.Thus, a common basis is needed to compare our suite of experiments.In Earth system modeling this is often achieved by tuning the export production or the nutrients to obtain comparable model simulations.But since we are more interested in the climate impact, we make sure that the primary production (PP) is comparable.
We are aware that PP is much less constrained by observations than for example nutrients.However, first we are more interested in the relative differences between the experiments and not in the absolute values.Second, adjusting the nutrients fields (e.g., DIC) would automatically adjust the carbon fluxes and will mask any change in carbon dynamics among the experiments.
For a robust comparison of the climate system between the different experiments, we adjusted the parameter βa (see Equation 1) in the partitioning function to obtain the same values for PP (Supporting Information S1).The parameter β a is used to tune the model because it is not constrained by observation and has already been changed between different studies (Ward & Follows, 2016;Ward et al., 2018).PP determines the two variables that drive the climate relevant feedback mechanisms we are interested in.First, PP determines the amount of phytoplankton biomass, and this biomass affects phytoplankton light absorption and thus the sea surface and air temperature.Second, PP is directly and indirectly the source of all different forms of organic biomass of which a part is sinking, leading to carbon export or carbon production and thus changes in atmospheric CO 2 concentrations and air temperature (Figure 3).After tuning the model, the values of our adjusted PP show only minor differences from an average value of 35.51 ± 0.61 Gt/yr (Table 2).The standard deviation is only 0.61 Gt/yr.We consider this value small enough to compare the climate systems between the experiments.
Please note that we evaluate the total climate impact encompassing changes in the heat budget and in the carbon cycle due to the biogeophysical or biogeochemical feedbacks.

Biogeochemical Properties
To compare the biogeochemical effects of phytoplankton light absorption and increasing ecosystem complexity on the climate system, we look at changes in biogeochemical properties.We first present the changes in surface phytoplankton biomass due to phytoplankton light absorption.Second, we compare the changes in the downward flux of organic matter.Third, we look at changes in the atmospheric CO 2 , that is allowed to evolve freely, and is strongly determined by the environmental conditions at the air-sea interface (partial pressure of CO 2 , sea surface temperature, the downward flux of organic matter, and the air-sea flux of CO 2 ).

Surface Phytoplankton Biomass
We first compare phytoplankton concentration in chlorophyll units of the different simulations (    comparable PP, the chlorophyll concentrations of the model runs with and without light absorption are different (1P1ZLA and 1P1Z: Figure 4c, 6P6ZLA and 6P6Z: Supporting Information S1).The general patterns between the simulations with low and high ecosystem complexity are similar except in the equatorial region (Supporting Information S1) where chlorophyll is lower; globally chlorophyll concentration differences are slightly lower between the simulations 6P6ZLA-6P6Z compared to 1P1Z-LA-1P1Z.Overall, the chlorophyll concentration is higher of 0.014 mg-Chl/m 3 in the model simulation 1P1ZLA compared to 1P1Z, similar to Manizza et al. (2005) but there are pronounced regional differences.The largest differences occur in the northern and southern polar regions where changes in sea-ice lead to a strong response in phytoplankton concentration.The presence or absence of sea-ice and thus light availability as well as the coarse resolution explain the rather sharp patterns of chlorophyll concentrations.In 1P1ZLA, the global averaged sea-ice thickness is reduced by ∼0.02 m and the sea-ice cover is diminished by 4.73% due to higher sea surface temperature (Figure 8c.As a consequence, light availability increases, stimulating the growth of phytoplankton.Moreover, the upwelling and mid-latitude regions show higher chlorophyll concentration in the simulation with phytoplankton light absorption in contrast to subtropical gyres where almost no differences occur.These regional patterns of higher chlorophyll concentration are controlled by the vertical velocity and the distribution of the nutrients.For instance, the upward vertical velocity in the upwelling regions along the western African coast and the south-western American coast are enhanced by 0.11 m/yr and by 0.54 m/yr, respectively (Supporting Information S1).These local increases of vertical velocity bring more nutrients to the surface and reduce the phosphate limitation in these regions (Supporting Information S1).The reduced phosphate limitation in 1P1ZLA permits an increase of phytoplankton growth and therefore a higher phytoplankton biomass.

Phytoplankton Light Absorption
To compare the downward flux of organic matter between our different simulations, we apply the approach by Toggweiler et al. (2003).In their approach, the two water masses North Atlantic Deep water (NADW) and Antarctic Bottom Water (AABW) are used as indicators for the nutrient "turnover" and thus the downward flux of organic matter.In models, this flux has upper and lower limits, constrained by the initial PO 4 concentrations in the deep water that is formed in the North Atlantic and Southern Ocean.The downward flux of organic matter can vary between these two limits via changes in deep water formation.The initial PO 4 concentration of the modeled AABW of 2.09 µmol/kg is slightly lower compared to average observed PO 4 concentrations in the deep water of 2.28 µmol/kg (World Ocean Atlas 2001, Conkright et al., 2002) and we use the AABW as the lower limit.The initial PO 4 concentration of the modeled NADW of 1.78 µmol/ kg is even lower and therefore the NADW represents the upper limit.Setting the "remineralization trajectories" for NADW and AABW (black diagonal lines) and plotting the DIC versus PO 4 concentrations of all grid points in the deep water (2,823 m) for the model runs 1P1ZLA and 1P1Z indicate the differences of the downward flux of organic matter at this specific depth.Figure 5a shows that in 1P1ZLA all points are significantly closer to the AABW remineralization trajectory, pointing toward a shallower downward flux of organic matter.On average, at 2,823 m depth, the PO 4 concentration decreases by 0.0318 µmol/kg and the DIC concentration decreases by 11.75 µmol/kg in 1P1ZLA compared to the reference simulation.
A closer look at the global DIC concentrations in the water column confirms a more shallower downward flux of organic matter in the simulation with phytoplankton light absorption.The concentration of DIC is higher in the surface ocean and the DIC gradient between surface and deep ocean ocean is smaller in 1P1Z-LA compared to 1P1Z (Table 3).While our parameter changes (β a ) in Equation 1 calculating the partitioning between dissolved and POM do not explain the differences in DIC (see Supporting Information S1), changes in the dynamics associated with phytoplankton light absorption are responsible for the higher DIC concentration at the surface.As a result of a positive feedback, the atmospheric and oceanic CO 2 concentrations increase in 1P1ZLA (see below).Additional effects contribute to the higher DIC concentration in 1P1ZLA.

Table 3 Concentration of DIC (mmol/kg) at the Surface and Benthic ocean
The upward vertical velocity, specifically in the upwelling regions is enhanced which reduces the penetration depth of sinking material and "traps" organic matter closer to the surface.This process strengthen the "short" near-surface loop of carbon cycling ("microbial loop").The deeper carbon cycling is weaker and less carbon is stored in the nonliving carbon pools and more in the living organic carbon pools.Not only phytoplankton but also zooplankton biomass is higher in the surface ocean; the latter is increased by ∼4%.Interestingly, a higher zooplankton biomass does not result in higher (deeper) export production of fecal pellets which would lead to higher DIC concentrations in the deep ocean.Instead, zooplankton respiration is higher and together with the overall increase in atmospheric CO 2 concentrations, this leads to higher DIC concentrations in the surface ocean.In our steady state system with a closed carbon cycle, carbon that is higher in one pool is lower in another pool.Accounting for light absorption clearly accelerates carbon turnover.

Increasing Ecosystem Complexity
Differences also occur in the distribution of the points between the model runs 1P1Z and 6P6Z (see Figure 5b).For 6P6Z, several points are closer to the NADW remineralization trajectory while others are closer to the AABW remineralization trajectory.We therefore calculate the average distance between all the points of 6P6Z and the remineralization trajectories.The average distance between the AABW limit and the points is 0.1807 µmol/kg while the average distance between the NADW limit and the points is 0.1486 µmol/kg.Our results show that the points of 6P6Z are slightly closer to the NADW remineralization trajectory, pointing toward a deeper downward flux of organic matter compared to 1P1Z.Along with the higher export production of POC (see Table 4) and lower dissolved organic phosphorus (Supporting Information S1), the concentration of DIC increases by 0.0123 µmol/kg because all organic matter is eventually remineralized close to the seafloor (see also [Ward et al., 2018], [Ridgwell et al., 2007]).The higher export production of POC is due to the different number of plankton groups and the different surface to volume ratio between 1P1Z and 6P6Z.Although the average equivalent spherical diameter for the entire population is the same between both simulations, the fraction that goes into DOM is higher for smaller organisms and lower for larger organisms.Due to the nonlinearity of the system (see Equation 1), a higher proportion of organic matter ends up in POC in the model run 6P6Z (Table 4).Again, our changes in the parameter (β a ) in Equation 1that calculates the partitioning between DOM and POM do not explain the differences (see Supporting Information S1).
Since the changes in ocean circulation are negligible in the simulations with different ecosystem complexity, we additionally compute the Apparent Oxygen Utilization (AOU) to further study the downward flux of organic matter on a global scale.AOU is a measure for biological activities of a water parcel since the last time it has been in contact with the atmosphere.It is computed as the oxygen saturation minus the oxygen concentration under the same temperature and salinity (Weiss, 1970).The higher the AOU, the greater the amount of oxygen removed by heterotrophic biological processes (respiration and remineralization).Since there is no sediment layer in the model, all organic matter is eventually Note.The third column is the partition of organic matter going into the dissolved phase.Abbreviations: DOC, dissolved organic carbon; DOM, dissolved organic matter; POC, particulate organic carbon.remineralized back to nutrients albeit at different depth levels depending on the organic matter form; DOM is remineralized in the upper ocean while POM is remineralized in the deeper ocean.We thus calculate the AOU in such a way that we take the difference between the oxygen saturation minus the oxygen concentration of the respective oceanic layer.At the surface, the AOU is similar between the simulations 1P1Z and 6P6Z but from 1,000 m depth downwards the AOU distribution diverge (Figure 6).At depth, a higher ecosystem complexity increases the AOU, indicating a deeper downward flux and a greater amount of oxygen being removed by biological processes.

Atmospheric CO 2 Concentration
Atmospheric CO 2 concentration is significantly higher in 1P1ZLA compared to 1P1Z (Table 5).A closer look at the spin-up phase reveals the mechanisms behind.Phytoplankton light absorption immediately affects the stratification, increases SST, reduces the solubility of gases and therefore increases the air-sea CO 2 flux.
To estimate the effect of the reduced CO 2 solubility on the air-sea CO 2 flux we used the values of 1P1Z for our calculations.The CO 2 fluxes are a function of several factors including solubility, atmospheric CO 2 , DIC, and the proportion of sea ice cover.By adopting the values for these factors from 1P1Z except for one, for which we take the value from 1P1ZLA, we are able to separate the individual effects during the spinup time.We find that by far the solubility has the largest effect on CO 2 fluxes and that the air-sea CO 2 flux increases from the beginning onwards with roughly 10% increase already after 500 years.As a consequence, the atmospheric CO 2 concentration rises and thereby generates a positive feedback, leading to a global temperature rise (see below).The changes in solubility and hence CO 2 fluxes are the immediate response of light absorption that increases the stratification and SST; other second order effects with changes in ocean physics arise during the course of the simulation.Yet, these changes such as a reduced sea-ice cover and enhanced upwelling only slightly modify the global CO 2 fluxes during the spin-up phase (<1% after 500 years).Finally, in steady state the difference in the atmospheric CO 2 concentrations between 1P1ZLA and 1P1Z is 18 ppm (Table 5), corresponding to 38 GtC.The largest contribution is due to the changes in solubility but an altered biogeochemistry additionally contribute to changes in atmospheric CO 2 concentration.As shown above, the downward flux of organic matter is shallower with phytoplankton light absorption, affecting the carbon cycle.Together with a shallower remineralization and zooplankton respiration, the surface concentration of DIC increases (Table 3).Taking the difference of the vertically-integrated DIC between 1P1ZLA and 1P1Z, the changes are equivalent to 35 Gt of CO 2 .This decrease in DIC inventory is mainly responsible for the higher atmospheric concentration in 1P1ZLA.Compared to solubility changes, the changes in the downward flux of organic matter, however, may explain a smaller part of the climate system's response.
In contrast, in the model runs with a higher ecosystem complexity changes in atmospheric concentrations are caused exclusively by biogeochemistry.Increasing the ecosystem complexity leads to a deeper downward flux of organic matter.Since all organic matter is remineralized at depth, the concentration of DIC increases at depth as well (see Figure 5b; Table 3).However, the changes are rather small and so the atmospheric CO 2 concentration only slightly decreases with a higher ecosystem complexity, as also suggested by the preformed PO 4 nutrients analysis (Supporting Information S1).The small decrease in atmospheric CO 2 concentration drives the small heat loss in the heat budget when ecosystem complexity increases (Table 6).

Phytoplankton Light Absorption
The changes in atmospheric CO 2 concentration significantly affect the system's heat budget (Table 6).The global surface atmospheric temperature (SAT) is higher by 0.45°C in the model runs with phytoplankton light absorption (Figure 7d).Along with a higher SAT, the ocean temperature increases globally in the model runs with phytoplankton light absorption as well (Figure 8d).Thus, the well-known effects of phytoplankton light absorption on ocean stratification (e.g., Manizza et al., 2008;Sonntag & Hense, 2011) are not visible as such.Yet, within the first 100 years of our simulations, we indeed find in the simulation with In steady state, we see the greenhouse effect with more longwave radiation being trapped in the atmosphere (Supporting Information S1).The higher greenhouse gas concentration changes therefore the overall heat budget of the climate system leading to global warming of the atmosphere and ocean.
The maximum difference of SAT between 1P1ZLA-1P1Z is 1.08°C in the Southern Ocean (Figure 7d).This value is slightly higher compared to previous model studies that show maximum values between 0.5-1K in model runs with phytoplankton light absorption compared to those without (Patara et al., 2012;Shell et al., 2003).However, in contrast to our model setup, Shell et al. (2003) use an ocean general circulation model with an uncoupled atmospheric model, thereby neglecting any interaction between the atmosphere and the ocean.Patara et al. (2012) use a coupled ocean-atmosphere general circulation model with a constant and prescribed atmospheric CO 2 concentration for their simulations.We argue that this strong response in SAT in our results is due to our experimental setup.We use a fully twoway coupled ocean-atmosphere model and assume a closed carbon cycle in which CO 2 is allowed to evolve freely.
The general pattern of differences in SST between 1P1ZLA-1P1Z is similar to SAT (Figures 7d and 8c) and we also find the same features in 6P6ZLA-6P6Z although the magnitude is lower in the latter (Supporting Information S1).Here the SST is strictly speaking not SST but the mean temperature in the upper surface layer, from the surface to 29.38 m depth, which for practical purposes is called SST.The global average heating of the ocean surface is 0.33°C between 1P1ZLA and 1P1Z, which is in accordance with previous modeling studies (Lengaigne et al., 2009;Manizza et al., 2005;Patara et al., 2012;Wetzel et al., 2006).We find a higher ASSELOT ET AL.SST even in the regions where differences in chlorophyll are small (Figure 4).This missing spatial coincidence between the chlorophyll and SST patterns can be explained by the model setup.Chlorophyll biomass is not subject to physical transport, but heat is transported by ocean currents, explaining why the patterns of the physical quantities are more smooth.Regional differences in SST are hardly visible.Both polar regions are characterized by minor changes in SST because the sea-ice dynamics limit the redistribution of heat.
The maximum local increase in SST is 0.47°C and is close to reported values in previous modeling studies (Lengaigne et al., 2007;Patara et al., 2012;Wetzel et al., 2006).On the other hand, observations show a local heating effect of 1.5-4°C caused by the absorption of light by phytoplankton surface blooms (Kahru et al., 1993;Sathyendranath et al., 1991); thus our model underestimates the local heating effect due to phytoplankton light absorption.

Increasing Ecosystem Complexity
An increase in ecosystem complexity (6P6Z-1P1Z) results in a global SAT decrease by 0.034°C (Figure 7e).This slight change is driven by the small decrease in atmospheric CO 2 concentration.The changes in SAT are less pronounced than the changes in SAT in the simulation with phytoplankton light absorption.The Southern Ocean is characterized by the largest SAT fluctuations with a local cooling of 1.64°C when ecosystem complexity increases (Figure 7e).The other world regions show a less pronounced cooling effect with a higher ecosystem complexity.Only a few grid cells in the Southern Ocean and in the North Atlantic show an increase in SAT.The global and regional changes in SST between 6P6Z-1P1Z follow the SAT patterns but are almost negligible (Supporting Information S1).

Summary and Conclusions
To study the relative importance of biogeophysical and biogeochemical climate relevant mechanisms, we implemented phytoplankton light absorption in the EcoGENIE model (Ward et al., 2018), and varied the complexity of the ecosystem by increasing the number of phytoplankton and zooplankton groups.In our simulations, the atmospheric CO 2 can evolve freely, affecting therefore the global heat budget.To obtain comparable primary productivity in all model runs, it was necessary to adjust the partitioning between the labile and refractory organic matter.We show that the climate system responds differently to our modifications in adding phytoplankton light absorption or increasing ecosystem complexity (Table 6).
A higher ecosystem complexity affects the plankton community, influencing the partitioning of organic matter going into the dissolved and particulate phase.Changes in the surface-to-volume ratio reduce the part of organic matter going into the dissolved phase and increase the part going into the particulate phase.As a result, the export production of POC increases while DOC decreases in 6P6Z.These changes in the downward flux of organic matter slightly affect the carbon cycle and with it the air-sea CO 2 flux.Hence, the atmospheric CO 2 concentration slightly decreases and the atmosphere cools down by 0.034°C when the ecosystem complexity increases.These small changes in atmospheric temperature and carbon cycle slightly reduces the sea surface temperature and slightly decreases the chlorophyll biomass.
Phytoplankton light absorption affects the climate system in various ways.Most notably, we find enhanced stratification, a higher SST and a reduced solubility of CO 2 that increases the air-sea CO 2 flux during the spin-up phase leading to higher atmospheric CO 2 concentration in steady state.In addition, the downward carbon flux is shallower with a stronger "microbial loop", contributing to higher atmospheric CO 2 concentration.Reduced sea-ice cover and enhanced upwelling only slightly affect the climate system.The sensitivity analysis indicate that by far the changes in CO 2 solubility have the largest effect on the climate system.All these changes lead to an increase by 0.45°C of the SAT with phytoplankton light absorption.
This study shows clearly that phytoplankton light absorption has a higher impact on the climate system than a higher ecosystem complexity.Therefore, we conclude that ESMs should include phytoplankton light absorption by default for climate change scenarios.
importance of ecosystem complexity and phytoplankton light absorption with the Earth System model (ESM) Eco Grid-ENabled Integrated ESM (EcoGENIE) • Phytoplankton light absorption

Figure 1 .
Figure 1.Schematic representation of EcoGENIE Earth system model modules used for this study.Black arrows represent the links between the different components.The addition of ECOGEM to the previous GENIE components forms the EcoGENIE Earth system model.

Figure 2 .
Figure 2. Modeled and observed Atlantic meridional overturning circulation (AMOC) stream function (Sv).Positive values represent clockwise overturning and negative values represent counter-clockwise overturning.Simulated AMOC for 1P1Z with the model configuration with: (a) 16 oceanic layers and (b) 32 oceanic layers.(c) Mean AMOC estimated from hydrographic sections from the World Ocean Circulation Experiment (Lumpkin & Speer, 2007).The gray line indicates the crest of the Mid-Atlantic Ridge and the white line represents the climatological mixed layer depth (Figure from Buckley and Marshall (2016)).
phyto-and 1 zooplankton species 1P1ZLA Simulation with additional phytoplankton light absorption 6P6Z Simulation with 6 phyto-and 6 zooplankton species 6P6ZLA Simulation with additional phytoplankton light absorption Fifth and last, Ward et al. (2018) allow primary production (PP) only in their surface layer, from the surface to 80.84 m deep while we allow PP until the sixth oceanic layer, from the surface to 221.84 m deep.

Figure 3 .
Figure 3. Schematic representation of the different links between the climate variables we analyze.The global atmospheric temperature is impacted by a biogeophysical mechanism (left, red color) and a biogeochemical mechanism (right, green color).Phytoplankton light absorption impacts directly sea surface temperature, biogeochemical properties and atmospheric CO 2 concentration, leading indirectly to changes in chlorophyll biomass and global atmospheric temperature.Increasing ecosystem complexity affects directly the export production of particulate organic carbon and thus the biogeochemical pump.As a consequence, atmospheric CO 2 concentrations and therefore global atmospheric temperature can be altered.
Figures 4a and 4b) as this variable affects phytoplankton light absorption and thus heat distribution.Despite ASSELOT ET AL.

Figure 4 .
Figure 4. Chlorophyll biomass (mgChl/m 3 ) for the model runs (a) 1P1ZLA and (b) 1P1Z.(c) Chlorophyll biomass difference between the two simulations; blue colors indicate lower while red colors indicate higher chlorophyll concentrations in the model run with light absorption.

Figure 5 .
Figure 5. Dissolved inorganic carbon (DIC) versus PO 4 (µmol/kg) composition for all grid cells at the ocean depth level 26 (2,823 m depth) showing the influence of North Atlantic deep water (NADW) and Antartic bottom water (AABW) on the composition of the deep water.The gray horizontal line gives the surface equilibrium of DIC for the simulations.The black diagonal lines show the phosphorus trajectories of NADW and AABW and as such, they define the upper and lower limits of the downward flux of organic matter.The results shown here are exclusively for the oceanic depth level 26 (2,823 m depth) but same results can be shown from depth level 24 to 29.(a) Red points indicate the simulation 1P1ZLA and blue points indicate the simulation 1P1Z.(b) Red points indicate the simulation 6P6Z and blue points indicate the simulation 1P1Z.

Figure 6 .
Figure 6.Vertical profile of Apparent Oxygen Utilization (µmol/kg) for the different simulations.The red curve represents the vertical profile of the simulation 1P1Z.The gray curve represents the vertical profile of the simulation 6P6Z.
Atmospheric CO 2 Concentration (ppm) After 1,000 years-Simulation for the Four Experiments With and Without Phytoplankton Light Absorption phytoplankton light absorption that the upper surface ocean is warmer and the deeper part cooler (Supporting Information S1), because phytoplankton absorbs light and shades the water column below.During the course of the model simulation into steady state, however, the decrease in the temperature-dependent solubility of CO 2 and changes in biogeochemistry (see above) lead to a rise of atmospheric CO 2 concentration in 1P1ZLA.
The values represent the average differences between the biogeophysical and biogeochemical scenarios.

Figure 7 .
Figure 7. Surface air temperature (°C) for the model runs (a) 1P1ZLA, (b) 1P1Z and (c) 6P6Z.Surface air temperature difference for the biogeophysical scenario: (d) Difference between 1P1ZLA and 1P1Z.Dark/brown colors indicate large while white/yellow colors represent small differences.Note that the differences between these two simulations are always positive.Surface air temperature difference for the biogeochemical scenario: (e) Difference between 6P6Z and 1P1Z.Blue colors indicate a higher surface air temperature for the model run 1P1Z while red/yellow colors indicate a higher surface air temperature for the model run 6P6Z.

Figure 8 .
Figure 8. Sea surface temperature (°C) for the simulations (a) 1P1ZLA and (b) 1P1Z.Blue colors indicate low while red colors represent high ocean temperatures.(c) Differences of sea surface temperature between the two simulations.White color indicates small and orange/brown colors represent large differences.Note that the differences between the model simulations are always positives.(d) Variation of the oceanic temperature with depth.The red line corresponds to the model run with phytoplankton light absorption.The blue line represents the model run without.The green line represents the difference between the model simulation with and without phytoplankton light absorption.

Table 1
Name and Description of the Four Simulations Conducted

Table 2
Global Values of Primary Production (Gt/Yr) for the Four Simulations The third column represents surface-to-deep gradient of DIC.

Table 6
Summary of the Impact of Phytoplankton Light Absorption and Increasing Ecosystem Complexity on the Different Climate Variables Described previously